The steady Navier-Stokes equations governing the flow of an incompressible viscous fluid in the rotation form in -shaped domain is considered. The weighted finite element method based on the definition of an generalized solution is constructed. The advantage of the proposed approach over classical approximations is numerically established. The modern elements of computational technologies to find the optimal parameters
of the proposed method are used.
∈ ( , 2 ] belongs to the space 21+ − (Ω), < 1. Moreover, the approximate finite
In [1], it was proposed to define the solution of elliptic boundary value problems with a singularity as an -generalized one. The approach allows us to introduce a weight space or a set, depending on the geometry of the domain and input data (right-hand sides, equation coeficients, boundary and initial data) to which an -generalized solution belongs. In [2, 3, 4], the existence, uniqueness and diferential properties of the elliptic problems solution are proved. In [5], a weighted analogue of the Ladyzhenskaya-Babuska-Brezzi condition for the Stokes problem is established. In [6, 7, 8, 9, 10], a weighted finite element method (FEM) for an approximate solution of elliptic problems with a (ℎ ) rate. singularity has been developed.
In the paper, an -generalized solution of the steady Navier-Stokes equations governing the flow of an incompressible viscous fluid in the rotation form in -shaped domain is defined. We use Picard’s iterative procedure [11] to find a solution of a nonlinear problem. Then, we construct a weighted ifnite element scheme based on the definition of an -generalized solution: 1) the functions of the ifnite element spaces satisfy the mass conservation law in a strong sense — Scott-Vogelius (SV) pair of 2nd order [12]; 2) basis functions are the product of SV functions and weight functions in some degree. This construction allows us better take into account the behavior of the solution in the neighborhood of the singularity point and increase the convergence rate of the approximate solution to the exact one to the first order with respect to the grid step ℎ, i. e. (ℎ) rate in 21, (Ω) norm. derived numerically.
Thus, we develop the numerical method overcomes the so-called pollution efect (see [13]). The same advantage for other hydrodynamics problems was achieved in [14, 15, 16]. The optimal values , ∗ and of the presented weighted FEM using the modern elements of computational technologies were
The paper consists of six sections. Section 2 is devoted to the definition of an -generalized solution. In Section 3, we present the weighted FEM. The iterative procedure for solving the systems of linear algebraic equations is constructed in Section 4. In Section 5, we shaw and discuss the results of computational experiments. Necessary conclusions are made in last section.
‖ ‖ = √ 12 + 22 are the measure and norm of ,
respectively.
Let Ω be a bounded, connected domain in the Euclidean space 2. Denote by Ω̄ and Γ the closure and boundary of Ω, respectively, Ω̄ = Ω ∪ Γ. Let = ( 1, 2) be an element of 2, where = 1 2 and
We write the steady Navier-Stokes equations governing the flow of an incompressible viscous fluid in the convection form: find a velocity field from
= ( ) = ( 1( ), 2( )) and a kinematic pressure = ( ) − ̄ △ + ( ⋅ ∇) + + ∇ = , div = where ̄ > 0 is the kinematic viscosity coeficient (inversely proportional to the Reynolds number), 0, = ( ) = ( 1( ), 2( )) and = ( ) = ( 1( ), 2( )) are given force field in Ω and boundary data on Γ, respectively. Denote by △, ∇ and div the Laplace, gradient and divergence operators in 2, Further, we introduce the necessary notation. Let = ( 1, 2), = ( 1, 2) and – scalar, then 1 2 ⋅ = 1 1 + 2 2, × = (− 2, 1) , curl = − 1 2 + 2 1 .
We have the identity
( ⋅ ∇) = ( curl ) × + ∇ 2.
It follows from the equality ( ⋅ ∇) + ( ⋅ ∇) = ∇( ⋅ ) + ( curl ) × + ( curl ) × and assumption that =
.
StoUkseisngeq(u3a),tiwonitsh: find=a ve+lo21c it2yffieolrdthe system (1), (2) we get the rotation form of the steady Navier and a Bernoulli pressure such that − ̄ △ + ( curl ) × + + ∇ = , div = 0
The system (4), (5) as well as (1), (2) is nonlinear due to the presence of the rotation term ( curl ) × in the momentum equations. The system on (4), (5) and term in particular we linearized by Picard’s procedure (see [11]).
Starting with an initial approximation (0) for which div (0) = 0 in Picard’s iteration constructs a sequence of solutions ( ( ), ( )) by solving the linear system: − ̄ △ ( ) + ( curl ( −1)) × ( ) + ( ) + ∇ ( ) = , div ( ) = 0 ( ) = Γ in on | | 1 1 2 respectively, where ( ) =
2 , | | = 1 + 2, ≥ 0 - integer.
Let 21, (Ω, ) for > 0 be a set of functions from the space 21, (Ω), meets the conditions ∫ 2 ( )
2 Ω̄⧵Ω′ ( ) ≥ 1 > 0, | ( )| ≤ 2( ( ) ) + ∈ Ω , ′ where
= 0, 1 and 2 21, (Ω). Denote by 2, (Ω, ) a positive constant which is not depend on , with the norm of a space
a set of functions from the space 2, (Ω) which subject to conditions (9) (only for = 0) with a norm of a space 2, (Ω). Let 02, (Ω, ) = { ∈ 2, (Ω, ) ∶ ∫ = 0}. = (6) (7) (8) (9) (6) (see [11]). in force. ( ( ), ( )
Note that the initial Bernoulli pressure in (7) need not be specified. If ̄ be a not too small and be a not too large, the steady Navier-Stokes equations (4), (5) have a unique solution ( , ) and the iterates ), = 1, 2, in (7), (8) converge to it as → ∞ for any choice of the arbitrary (0) satisfying Note that for a linearized system (7), (8) the conservation laws of the mass and momentum remain In the article, we consider the special case of a bounded polygon domain Ω. Let Ω be a -shaped domain with one reentrant obtuse corner equals to 3 on the boundary and its vertex coincides with the origin. We define an -generalized solution in e2ach Picard’s iteration of the problem (7), (8) and construct the efective weighted FEM. Thus, we solve the nonlinear problem (4), (5) governing the lfow of a incompressible viscous fluid in the rotation form and show the advantage of our approximate method over the classical approaches in a -shaped domain by the computational simulations.
Let us introduce the notation and define necessary spaces of generalized functions. Denote by ′ Ω = {{ ∈ Ω̄ ∶ ‖ ‖ ≤ < 1, > 0} a part of a -neighborhood of a point (0, 0) contained in Ω̄. Let ′ ‖ ‖, ∈ Ω , , ∈ Ω̄ ⧵ Ω ′
be a weight function.
Denote by 2, (Ω) and 21, (Ω) the spaces of functions ( ) with a bounded norms ( ) = and ‖ ‖ 2, (Ω) = √ √ ∫ 2 ( ) 2
( ) ‖ ‖ 21, (Ω) = ‖ ( )| ( )|‖22(Ω) + ‖ ( )| 1 ( )|‖22(Ω), 2, (Ω, ) ( diferentiable functions with a compact support in (Ω, ) ⊂ 21, (Ω, )) be a closure by
21, (Ω) norm of a set of the infinitelyΩ comply with the conditions (9). We will say ( ) ∈ 21,/2(Γ, ), if exists a function Φ( ) ∈ 21, (Ω, ) such that Φ( )|Γ = ( ) and ‖ ‖ 21,/2(Γ, ) = inf ‖Φ‖ 21, (Ω).
( ( ), ) = ∫ [ 2 ( ) ⋅ + ̄ ∇ ( ) ⋅ ∇( 2 ) + 2 (( curl ( −1)) × ( )) ⋅ ] , ( ( ), ) + ( , ( )) = ( ), ( ( ), ) = 0 hold, where bilinear and linear forms are as follows 20, (Ω, ) ifelds 1/2(Γ, ) and
(Ω, ) on Γ and in Ω, respectively.
We introduce the concept of an -generalized solution for the linearized problem (7), (8). Definition 1. The pair ( ) ∈ 21, (Ω, ) and ( ) ∈ 20, (Ω, ) is called an -generalized solution of the problem (7), (8), where ( ) satisfies a condition (8) on Γ for any pair ∈ 21, (Ω, ) and ∈ 1 2,
and Γ are triangulation nodes subsets for the velocity ( , ( )) = − ∫ ( ) div ( 2 ) , ( ( ), ) = − ∫ ( 2 ) div ( ) , ( ) = ∫ 2 ⋅ and ∈ 2, (Ω, ), ∈ 21,/2(Γ, ), ≥ ≥ 0.
Perform triangulation Υℎ based on the barycentric partition of the elements triangulation ℎ of the domain Ω. Then, we divide each element ∈ ℎ of the quasi-uniform (macroelement) into three triangles ). Let
(finite element), ∈ Υℎ (their common vertex is in the barycenter of the macroelement and be the vertices and midpoints of the sides ∈ Υℎ, respectively. Introduce the √ and ‖ ‖ 21, (Ω) = ‖ 1‖2 21, (Ω) + ‖ 2‖2 21, (Ω)
For the vector field = ( 1, 2) we define sets 2, (Ω, ) and 12, (Ω, ) such that ∈ 2, (Ω, ) and ∈ 21, (Ω, ), respectively, with a bounded norms ‖ ‖ 2, (Ω) = ‖ 1‖22, (Ω) + ‖ 2‖22, (Ω)
for the first set √ for the second one. Similarly, we define the sets of vector 1) 2) notation of sets: ifeld components in = to the node at the appropriate
. = { ∪ }, where Ω and on Γ, respectively; = { } of triangulation nodes for the Bernoulli pressure, where the node an exact match ⋃ ∈Υℎ We denote by Ωℎ =
the totality of the finite elements with sides of order ℎ. Next, we describe the Scott-Vogelius (SV) element pair (see [12]). For the components of the velocity field, we ℎ = { ℎ ∈ (Ω) ∶ ℎ| ∈ 2( ), ∀ ∈ Υℎ}, ℎ = ℎ × ℎ; ℎ = { ℎ ∈ 2(Ω) ∶ ℎ| ∈ 1( ), ∀ ∈ Υℎ, ∫ ℎ = 0}. use polynomials of the second degree ( ℎ), and for the pressure — the first one ( ℎ):
The SV pair has useful feature, namely div ℎ ⊂ ℎ. Next, we represent special basis functions and construct a scheme of the weighted finite element method. To each node ∈ Ω ( ∈ ) we associate the basis function where ∈ ℎ, ( ) = , and ∗ are real parameters. pressure ,ℎ
( ) we seek as a
The spaces ℎ and ℎ for the components of the velocity field and pressure are defined as linear span of the basis functions { }
and { } , respectively. Let 0ℎ be a subspace of ℎ ∶ 0ℎ = { ℎ ∈ ℎ ∶ ℎ( )| ∈ Γ = 0}. The approximate components of the velocity field ,ℎ = ( ,ℎ, 1, ,ℎ, 2) and ( ) ( ) ( ) ( ) = ∗ ( ) ⋅ ( ), ( ( ) = ∗ ( ) ⋅ ( )), = 0, 1, … , ( = 0, 1, …), = 0, 1, … ( ∈ ℎ, ( ) = , , = 0, 1, …); = { 1, = , 0, ≠ , , ∗ 21, (Ωℎ, ), 0ℎ ⊂ 2, (Ωℎ, ), ℎ ⊂ 20, (Ωℎ, ). 1 ( ) ( )
( ) ,ℎ, 1( ) = ∑ 2 ( ), ,ℎ, 2( ) = ∑ 2 +1 ( ), ,ℎ ( ) = ∑ ( ) ( ) ( ) ( ), ( ) where ( ) = − ∗( [ /2]) ̃ , ( ) = − ∗( ) ̃ ( ). The coeficients ( ) and ( ) in (10) are found as a result of solving a system (11), (12) (see below). Let ℎ = ℎ × ℎ, 0ℎ = 0ℎ × 0ℎ and ℎ ⊂ the problem (7), (8) for any pair ℎ ∈ 0ℎ and ℎ ∈ ℎ the equalities
( ) Definition 2. The pair ,ℎ ∈ ℎ and
,(ℎ ) ∈ ℎ is called an approximate -generalized solution of ( ) ℎ ( ,ℎ , ) + ( ℎ
( ) , ,ℎ ) = ( ℎ) and ( )
ℎ ( ,ℎ , ) = 0 hold, where ∈ 2, (Ω, ), ∈
21,/2(Γ, ), ≥ ≥ 0. get a system of linear algebraic equation:
Thus, we construct a weighted FEM to find an -generalized solution for the problem (7), (8). We ( ) ( ) + ( ) = and
( )
= , ( ) , 2 , … , 1
( ) ( ) ( ) ( ) , 3 , …) , ( ) = ( 0 , 1 , 2 , …) and be a vector of values of the linear form where ( ) = ( 0
( ) ( ).
and thus we will approximately solve the original nonlinear problem in rotation form (4), (5): Now, we present an iterative procedure for solving the sequences of systems view (12), = 1, 2, 3, … 1. Let (0) and (0) be an arbitrary vectors such that (,ℎ0) ( )| ∈ Γ = ( ), div (,ℎ0) ( )| ∈ = 0 (for example (,ℎ0) ( )| ∈
(0) = ) and ,ℎ ( )| ∈ Ω 2. Realize the Picard’s procedure = 0, 1, 2, … until the stopping condition is fulfilled: = 0. a) Let 0( ) ∶= ( ) and (0 )
∶= ( ); b) We construct an internal convergent iterative process (see [17]). For = 0, 1, … , − 1 ∶ −1 ( ) ( ) +1 = ( ) + ̂ ( − ( )
− ( )) ( ) +1 = + ̂ −1
( ) +1; (10) (11) (12) that on each ∈ Υℎ ∶ c) Let ( +1) ∶= ( ) and ( +1) ∶= ( )
, where ̂ and ̂ are the preconditioning matrices to and = −1 , respectively.
At first, we build a preconditioner ̂ applying an incomplete LU factorization. We employ the GMRES(5)-method (see [18]). If we have error 0 = ̂ −1( − ) for the problem = , then the Arnoldi procedure will build an orthonormal basis of the subspace: Span{ 0, ̂ −1 0, … , ( ̂ −1 )4 0} .
Further, we construct an auxiliary matrix ̃ to ̂ , which is the weight mass matrix , ∗, ̄ , such ( , ∗, ̄ ) =
1 ̄ ∫ 2( + ∗) ( ) ( ) , , = 0, 1, … .
After that, we define a diagonal matrix ̄ = ̄ , ∗, ̄ , where ̄ , ∗, ̄ (
) = ∑( , ∗, ̄ ) . to determining the vector ⋄ ∶= ̂ −1 we need to find a solution of the internal procedure:
It is known (see [19]), that such diagonal lumping ̄ is a good preconditioner to matrix ̃ . In order 1) 0 = ; 3) ⋄ = . 2) = −1 + ̄ −1( − ̃ −1) ( = 1, … , ); We use the GMRES(5)-method: (Span{̄, ( ̄ −1 ̃ )̄, … , ( ̄ −1 ̃ ) 4 ̄}, ̄ = ̄ −1( − ̃ −1)).
and 2 ( ) = −1 + ℎ, respectively, into elementary squares , where , = 0, … , , ℎ Let Ω = (−1, 1)×(−1, 1)⧵[0, 1]×[−1, 0]. Then we split Ω̄ by the horizontal and vertical lines 1( ) = −1+ ℎ, number. After that, we divide each right corner) into two triangles
by the diagonal (the lower left corner connects to the upper (macroelements). Further, each macroelement is partitioned into three triangles (barycentric partition). Consider a solution ( , ) of nonlinear problem (4), (5) which has a singularity in the vicinity of the reentrant corner = 32 with apex at the origin (0, 0) ∶ 1( 1, 2) = ( 1 + 22)2 ((1 + ) ( 1, 2) ⋅ sin(arctg 2 ) + ( 1, 2) ⋅ cos(arctg 2 2
2( 1, 2) = ( 1 + 2 ) ( ( 1, 2) ⋅ sin(arctg 2 ) − (1 + ) ( 1, 2) ⋅ cos(arctg 2 2 2 2 1 = 2 ,
– even 1 )
) , 1 ) ) , where + + 1
1 ( 1, 2) = ( 1 + 2 ) 2 2
( −1 (1 + )2 ( 1, 2) + ( 1, 2) 2 ( 1, 2) = cos((1 − ) arctg 2 ) − cos((1 + ) arctg 2 )+ sin((1 + ) arctg 21 ) ⋅ cos( )
sin((1 − ) arctg 21 ) ⋅ cos( ) ( 1, 2) = (1 + ) sin((1 + ) arctg 2 ) − (1 − ) sin((1 − ) arctg 2 )+ +cos((1 + ) arctg 2 ) ⋅ cos( ) − cos((1 − ) arctg 2 ) ⋅ cos( ), 1 )
, 1 1 + 1 1 − 1 − 1 1 The error between the generalized solution ℎ
and exact one in the 21(Ω) norm.
The error between an -generalized solution ,ℎ and exact one in the 21, (Ω) norm, for diferent values the given limit values , = 1, 2, oΩf the generalized solution = ∗ = 0, = 1.
(in percentage of their total number), where the errors are more than
N= the given limit values , = 1, 2, of an Ω
- generalized solution = 1.8, ∗ = −0.31, = 0.014. (in percentage of their total number), where the errors , are more than
N=
Numerical experiments were carried out on grids with diferent steps 2 (Ω). It is a typical situation in non-convex polygonal domains.
1 Let = ̄ = 1 and = 0.54448. The pair of functions ( , ) is analytic in Ω̄ ⧵ (0, 0), but ∉ 22(Ω) . The errors of the gen were derived numerically. eralized (classical FEM with = 0, = 1, ∗ = ∗ = 0) and -generalized (presented weighted FEM) solutions were determined using the modulus of the diference between the exact solution , and approximate one at the nodes , i. e. = | ( ) − ℎ, ( )| for the generalized solution and = | ( ) − ,ℎ, ( )| for an -generalized one, where ∈ , = 1, 2, and also in the norms of generalized functions. See Figures 1-2 and Tables 1-4. The optimal values of parameters , ∗ and Ω
1 ℎ 1
The results of computational experiments for the steady Navier-Stokes equations (4), (5) lead to the following conclusions:
1) An approximate -generalized solution by the weighted FEM converges to the exact one with a (ℎ) rate in the 21, (Ω) norm (see Table 2), while the approximate generalized solution by the classical FEM converges to the exact one with a (ℎ0.54) rate in the 12(Ω) norm (see Table 1). In other words, the proposed method suppresses the so-called pollution efect [13].
2) For all values of , and ∗ from the range of optimal values ( ∼ ℎ, ∼ 2 and ∗ ∼ 1 − ) an approximate -generalized solution converges to the exact one with a (ℎ) rate in the 12, (Ω) norm.
3) The number of nodes and their surroundings by using a weighted FEM, in which the values of the absolute errors , , = 1, 2, do not exceed the given values, increases with and is much more then by using the classical FEM (see Tables 3-4) and Figures 1-2.
The reported study was supported by RSF according to the research project No. 21-11-00039. Computational resources for the numerical experiments were provided by the Shared Services Center "Data Center of FEB RAS". respectively. ∗ = −0.35): ) = 140, ) = 280 and ) = 140, ) = 280 for the 1st and 2nd components of ,ℎ , lem, Dokl. Akad. Nauk SSSR 309 (1989) 1318–1320. [1] V. A. Rukavishnikov, On the diferential properties of -generalized solution of Dirichlet prob[2] V. A. Rukavishnikov, On the existence and uniqueness of an -generalized solution of a boundary value problem with uncoordinated degeneration of the input data, Doklady Mathematics 90 [3] V. A. Rukavishnikov, E. V. Kuznetsova, The -generalized solution of a boundary value problem with a singularity belongs to the space 917. doi:10.1134/S0012266109060147.
2, + /2+ +1(Ω, ), Diferential Equations 45 (2009) 913– +2 [4] V. A. Rukavishnikov, S. G. Nikolaev, On the -generalized solution of the Lamé system with corner singularity, Doklady Mathematics 92 (2015) 421–423. doi:10.1134/S1064562415040080. [5] V. A. Rukavishnikov, A. V. Rukavishnikov, Weighted finite element method for the Stokes problem with corner singularity, Journal of Computational and Applied Mathematics 341 (2018) [6] V. A. Rukavishnikov, A. Y. Bespalov, An exponential rate of convergence of the finite element method for the Dirichlet problem with a singularity of the solution, Doklady Mathematics 62 (2000) 266–270. [7] V. A. Rukavishnikov, E. V. Kuznetsova, A finite element method scheme for boundary value problems with noncoordinated degeneration of input data, Numerical Analysis and Applications 2 (2009) 250–259. doi:10.1134/S1995423909030069. [8] V. A. Rukavishnikov, E. I. Rukavishnikova, Numerical method for Dirichlet problem with degeneration of the solution on the entire boundary, Symmetry 11 (2019) 1455. doi:10.3390/ sym11121455. [9] V. A. Rukavishnikov, A. O. Mosolapov, E. I. Rukavishnikova, Weighted finite element method for elasticity problem with a crack, Computers and Structures 243 (2021) 106400. doi:10.1016/ j.compstruc.2020.106400. [10] V. A. Rukavishnikov, O. P. Tkachenko, Dynamics of a fluid-filled curvilinear pipeline, Applied
Mathematics and Mechanics 39 (2018) 905–922. doi:10.1007/s10483-018-2338-9. [11] M. Benzi, H. Golub, J. Liesen, Numerical solution of saddle point problems, Acta Numerica 14 (2005) 1–137. doi:10.1017/S0962492904000212. [12] L. R. Scott, M. Vogelius, Norm estimates for a maximal right inverse of the divergence operator in spaces of piecewise polynomials, Mathematical Modelling and Numerical Analysis 19 (1985) 111–143. [13] H. Blum, The influence of reentrant corners in the numerical approximation of viscous flow problems, volume 30 of Numerical Treatment of the Navier-Stokes Equations, Springer, 1990. doi:110.1007/978-3-663-14004-7-4. [14] V. A. Rukavishnikov, A. V. Rukavishnikov, New approximate method for solving the Stokes problem in a domain with corner singularity, Bulletin of the South Ural State University. Ser. Mathematical Modelling, Programming & Computer Software 11 (2018) 95–108. doi:10.14529/ mmp180109. [15] V. A. Rukavishnikov, A. V. Rukavishnikov, New numerical method for the rotation form of the
Oseen problem with corner singularity, Symmetry 11 (2019) 54. doi:10.3390/sym11010054. [16] V. A. Rukavishnikov, A. V. Rukavishnikov, The method of numerical solution of the one stationary hydrodymics problem in convective form in -shaped domain, Computer Research and Modeling 12 (2020) 1291–1306. doi:10.20537/2076-7633-2020-12-6-1291-1306. [17] J. H. Bramble, J. E. Pasciak, A. T. Vassilev, Analysis of the inexact Uzawa algorithm for saddle point problems, SIAM Journal on Numerical Analysis 34 (1997) 1072–1092. doi:10.1137/ S0036142994273343. [18] Y. Saad, Iterative methods for sparse linear systems, 2nd. ed., SIAM, Philadelphia, 2003. [19] M. Olshanskii, A. Reusken, Analysis of a Stokes interface problem, Numerische Mathematik 103 (2006) 129–149. doi:10.1007/S00211-005-0646-X.