Numerical Analysis and Comparison of Three Iterative Methods Based on Finite Element for the 2D/3D Stationary Micropolar Fluid Equations

In this paper, three iterative methods (Stokes, Newton and Oseen iterative methods) based on finite element discretization for the stationary micropolar fluid equations are proposed, analyzed and compared. The stability and error estimation for the Stokes and Newton iterative methods are obtained under the strong uniqueness conditions. In addition, the stability and error estimation for the Oseen iterative method are derived under the uniqueness condition of the weak solution. Finally, numerical examples test the applicability and the effectiveness of the three iterative methods.


Introduction
The micropolar fluid equations can be used to describe the flow behavior of polymeric fluids. Based on the classical Navier-Stokes equations, Eringen [1] first proposed the micropolar fluid equations in 1965. The micropolar fluid equations postulate that every particle not only has mass and velocity but also has microinertia and microrotation, and can support not only stress, but couple stress. Micropolar fluids have been widely used in modern industry, biology, engineering and other fields. The dynamic behavior of fluid flow, such as particle suspensions, liquid crystals, lubrication and animal blood can be described by micropolar fluid equations [2,3].
Many scholars are devoted to study the well-posedness of solutions for the micropolar fluid equations. Galdi and Rionero [4] discussed the existence and uniqueness of weak solutions for the initial boundary value problem of the micropolar fluid equations. Rojasmedar and Boldrini [5] proved the global existence of weak solutions by the Galerkin method. Łukaszewicz [6,7] proved the local existence of strong solutions and the global existence of weak solutions by using the linearization method and fixed point theorems. Dong and Chen [8] proposed the regularity criteria of weak solutions for the 3D micropolar fluid equations.
In the past decades, the finite element method (FEM) has been widely used to solve the fluid dynamics equations. He and Li [9] introduced the three iterative methods for Navier-Stokes equations, and derived stability conditions for different iterative algorithms. Dong et al. [10] discussed and analyzed the three classical iterative methods for MHD equations. Based on the asymptotic behavior of the perturbation solutions of the steady Navier-Stokes equations, He [11] proposed the Euler implicit/explicit iterative scheme and proved the corresponding stability condition. In order to improve the calculation efficiency of Navier-Stokes equations, Xu [12], Xu and He [13] considered the two-level methods based on the iterative solutions. Layton [14,15] proposed and analyzed the two-level method to solve the nonlinearity. Huang et al. [16] proposed the two-level stabilized finite element method based on local Gauss integration. Readers can refer to [17][18][19][20] and the references cited to further understand the finite element iterative methods.
Recently, the FEM has been applied to solve the micropolar fluid equations. Ortega-Torres and Rojas-Medar [21] employed the fully discrete penalty finite element method to solve micropolar fluid equations, and proved optimal error estimates of linear velocity, pressure, and angular velocity. Nochetto et al. [22] proposed the first-order semi-implicit fully discrete finite element method. Jiang and Yang [23] proposed some projection methods and analyzed the stability and error estimates for the classical first-order projection scheme. Maimaiti and Liu [24] proposed first-order and second-order pressure-modified projection methods, and analysed the stability of time semi-discrete and fully discrete schemes.
According to the above literature summarization, the finite element method and three iterative methods have been successfully adopted to simulate Navier-Stokes equations, and some numerical achievements about the unsteady micropolar fluid equations are obtained. Numerical analysis and comparison of finite element iterative methods for the stationary micropolar fluid equation have not been reported in the open literature. It is well known that the stationary micropolar fluid equations are a strong coupled nonlinear system. This system contains two nonlinear terms, and velocity u, pressure p and angular velocity ω are coupled together. When the viscous effect in the micropolar fluid equations is weak, which means the nonlinear effect is strong, the main challenges of numerically solving it are how to construct efficient and stable iterative algorithms, rigorous theoretical analysis and application in practical problems. Based on the above difficulties, the three iterative methods are considered.
In this paper, three iterative methods are presented to simulate the stationary micropolar fluid equations for the different viscosity effect. The uniform stability and convergence of the proposed methods are also analysed. The main conclusions are reported in Theorems 5-9. In the numerical example section, we compared the effectiveness of the three iterative methods. Meanwhile, the classical benchmark problems and the bearing lubrication problem with inhomogeneous boundary conditions are simulated.
The rest of this paper is organized as follows. In Section 2, some necessary theories of Sobolev spaces are introduced. In Section 3, the stability and error estimation of the finite element discrete problem is provided. In Section 4, the main part of this paper, the stability and error estimation of the three iterative methods under different uniqueness conditions are analyzed and compared. In Section 5, some numerical tests are provided to illustrate the correctness of the theoretic analysis and contrast the effectiveness of the proposed methods. In Section 6, the main conclusions of this paper are presented.

Preliminaries
In this paper, we assume that Ω is a convex polygonal/polyhedral domain with Lipschitz continuous boundary in R d , d = 2 or 3. We consider the stationary micropolar fluid equations as follows [1,6]: where u is the velocity, p is the pressure, ω is the angular velocity, f and g are the external force terms, ν is the kinematic Newtonian viscosity, ν r is the dynamic microrotation viscosity, and the positive constants c a , c d , c 0 are new viscosities connected with the asymmetry of the stress tensor, respectively. In order to simplify the calculation, we set ν 0 = ν + ν r , c 1 = c a + c d , c 2 = c 0 + c a − c d > 0, and introduce the general definition of the curl operator as In particular, when d = 2, let For simplicity, let consider the following homogeneous boundary conditions: Let W k,r (Ω) be the standard Sobolev spaces for all nonnegative integers k and r with norm · k,r . As usual, the Hilbert space H k (Ω) = W k,2 (Ω) are equipped with norm v k when r = 2. To obtain the weak formulation of (1)-(2), we introduce the following Sobolev spaces For simplicity, we employ the product space W(Ω) = X × X with the usual graph where ·, · denotes duality product between the function spaces H 1 (Ω) d and its dual. For the purpose of the subsequent analysis, the following bilinear and trilinear forms are introduced: The following estimates about the trilinear form b(·, ·, ·) are classical [18,20,25].
where N is a constant, γ 0 is a positive constant that may dependent on Ω. The following Sobolve spaces inequalities are used frequently With the above notations, the variational formulation of (1)-(2) reads: find ((u, ω), p) ∈ W × M such that According to the above estimates, H'older inequality, Sobolev inequalities and the following inequalities the following properties of A(·, ·) and B(·, ·, ·) can be derived: where C max = max{ν 0 + 2 √ 2ν r γ 0 , c 1 + c 2 d + 2 √ 2ν r γ 0 + 4ν r γ 2 0 }, C min = min{ν 0 , c 1 − c 2 d} and λ = √ 2N. Note that the bilinear form d(·, ·) is continuous on W × M and satisfies the classical inf-sup condition [6,18]: Next, let us define the following norms: For the sake of convenience in writing, we set Łukaszewicz in [6] proved the existence and uniqueness of the weak solution of (1) and (2) by considering the auxiliary elliptical linear problem and adopting the Leray-Schauder principle. Next, we give a new proof of well-posedness by using the Banach fixed point theorem.
Theorem 1. If ν, ν r , c a , c d and c 0 satisfy the uniqueness condition then, there exists a unique solution to (9) and satisfies Proof. Let us define the function space For given (u, ω) ∈ Z, we consider the following linear problem: From the saddle-point theory [18], the problem (17) possesses a unique solution ((w, φ), p) ∈ W × M. Taking (v, ψ) = (w, φ) and q = p, using (11) and (12) we can drive the (16).
Next, let us prove the uniqueness condition (15). From the problem (17), we can define the mapping Making a difference for i = 1 and i = 2 in (18) and By (11), (13) and (15), we have According to the Banach fixed point theorem, we know that Φ has a fixed point in Z, which is the solution of problem (9). This completes the proof.

Finite Element Approximation
In this section, let T h be a regular and quasi-uniform partition of Ω into triangles for d = 2 or tetrahedras for d = 3, and h = max K∈T h h k is a real positive parameter tending to 0.
Next, we employ the following conforming finite element pairs to approximate the velocity, angular velocity and pressure.
whereb is the bubble function, and P 1 (K) denotes the space of polynomials of a degree equal to 1 on K. For convenience, let W h = X h × X h . In this case, (W h , M h ) satisfies the following assumptions [6,17,18].

Assumption 1.
There exists a positive constant β 0 that only depends on Ω such that Assumption 2. There exists a mapping R h : and a mapping Q h : With the above discrete spaces, the finite element approximation of (9) reads: find Similarly, let us define the discrete analogue of space V as Furthermore, we have the following discrete estimates [9,17]: Next, we prove the existence and uniqueness, stability and convergence of the discrete problem (26). Lemma 1. The trilinear term B(·, ·, ·) satisfies the following estimate Proof. The estimate can be easily derived by (4)- (8) and H'older inequality.
Proof. The existence proof is obvious. Next, let us prove the stability and uniqueness.
Proof. The proof is put in "Appendix B".

Three Iterative Methods
Considering the effect of different viscosity coefficients, the three iterative methods for micropolar fluid equations are proposed. Stokes iterative method: Newton iterative method: Oseen iterative method: The initial value ( Next, we establish the stability and error estimates of the three iterative methods. For simplicity, let (e n ,

Lemma 2.
Under the conditions of Theorem 1, Assumptions 1 and 2, the initial value satisfies the initial errors (e 0 , ξ 0 ) and η 0 satisfy the following bounds (36) and applying (11), the first term of (37) can be deduced. Similarly, choosing (v, ψ) = (A h u 0 h , A h ω 0 h ) and q = 0, the second term of (37) can be derived.

Theorem 5. Under the conditions of Theorem 1, Assumptions 1 and 2, if
and (e m , ξ m ), η m satisfy the following bounds for all m ≥ 0.
then (u m h , ω m h ) defined by the Newton iterative method satisfies and (e m , ξ m ), η m satisfy the following bounds for all m ≥ 0.

Theorem 7. Under the conditions of Theorem 1, Assumptions 1 and 2, if
then (u m h , ω m h ) defined by the Oseen iterative method satisfies (e n , ξ n ) and η n satisfy the following bounds for all m ≥ 0.
In the next theorem, another expression of the error estimates will be given, which can be controlled by (u n h − u n−1 h , ω n h − ω n−1 h ).

Theorem 8. Under the conditions of Theorems 5-7, there hold
and (e m , ξ m ), η m (m ≥ 1)) satisfy for the Stokes and Oseen iterative methods, and for the Newton iterative method.

Remark 1.
The error factor σ is not easy to compute during program implementation. Hence, we can adopt the error estimate from Theorem 8. In this case, the error can be controlled by the iterative error and the mesh size h.

Theorem 9.
Under the conditions of Theorems 3-7, the optimal error estimates of the three methods satisfy for the Stokes iterative method with 0 < σ < 2/5; for the Newton iterative method with 0 < σ < 5/11; for the Oseen iterative method with 0 < σ < 1.

Numerical Examples
In this section, four numerical examples are presented to test the correctness of the theoretical results and the validity of the proposed methods. In the first two examples, the three iterative methods are implemented for different viscosity coefficients. Meanwhile, the triangular cavity problem and bearing lubrication problem are also simulated in the last two examples. The P1b − P1 − P1b finite element pair is used to approximate u, p and ω. We choose the public domain finite element software FreeFem++ [26] to implement the algorithms. For numerical implementations, the iterative tolerance is set as 1.0 × 10 −6 .
When the viscosity coefficients ν = ν r = 0.5 × 10 −3 , only the Oseen iteration can maintain convergence. The results are agreed with the ones predicted by the theoretical analysis. Since the results of 2D and 3D are similar, we present only the 2D results here. Please refer to the Appendix C for 3D results.

Driven Cavity Flow
In this example, the classical benchmark problem, known as driven cavity flow is considered. Let the domain Ω = [0, 1] 2 , c a = c d = 0.5, c 0 = 1, f = 0, g = 0. The boundary condition for velocity reads u = (1, 0) on y = 1, and u = (0, 0) on the other three boundaries. Meanwhile, the angular velocity satisfies the homogeneous Dirichlet boundary condition.
In Table 4, we compare the applicability of the three iterative methods under different viscosity coefficients. As the table shows, the Oseen iteration is the most flexible method. This is consistent with our theoretical analysis. Thus, the Oseen iterative method is adopted. In Figures 2-4, the velocity streamlines, angular velocity contour lines of ω = (0, 0, ω 3 ) and pressure isobars for three different viscosity coefficients ν = ν r = 0.5, ν = ν r = 0.005, ν = ν r = 0.001 are displayed.   It can be seen that when the viscosity is small, the velocity streamlines, angular contour lines and pressure isobars satisfy symmetry. With the viscosity decreases, the main vortex center moves towards the geometric center of the cavity. Meanwhile, the secondary vortices appear near the corner of the cavity. The flow property of micropolar fluids is almost identical with the Navier-Stokes equations in [27].

Triangular Cavity Flow
Motivated by [28], the flow of micropolar fluids inside an isosceles triangular cavity is considered. Figure A1 shows the schematic diagram of flow geometry and the computational mesh. Similarly, we only consider the Oseen iterative method in this example. Set c a = c d = 0.5, c 0 = 1, f = 0, g = 0, the boundary conditions are as follows: Similar to square cavity flow, when the hydrodynamic viscosity decreases, the velocity streamlines, angular velocity contour lines and pressure isobars no longer satisfy symmetry. In summary, the flow property is still exhibited as laminar flow. Please refer to the Appendix D for Figures A2-A4.

The Bearing Lubrication Problem
In ordinary operating conditions, the lubricant liquid can be considered as fluid suspension. Thus the micropolar fluid equations can be used to simulate the lubrication problem [29]. In this example, the application of micropolar fluids in a non-concentric bearing lubrication problem is considered. The fluid region is a ring domain between the outer boundary Γ 1 with radius r 1 and the inner boundary Γ 2 with radius r 2 .
In this example, we only present the numerical experimental results for the Newton iterative method since it is faster than the Stokes and Oseen iterative methods. Figures 5-8 show the evolution diagram of velocity components (horizontal velocity, vertical velocity), angular velocity and pressure at three different rotational angular velocities ω r . The numerical results show that the velocity components, angular velocity and pressure increase with the increase in the rotational angular velocity ω r . The pressure increase is the most obvious. Therefore, the bearing is capable of supporting higher loads. The angular velocity contour lines tend to revolve around the inner circle. This means that the higher the rotational angular velocity ω r , the stronger the micropolarity effect of the fluid.

Conclusions
In this paper, three iterative methods are analyzed and compared for the micropolar fluid equations. The strong stability conditions for the Stokes and Newton iterations can be determined by the parameter σ, i.e., 0 < σ < 2/5 and 0 < σ < 5/11, respectively. Furthermore, the uniqueness condition 0 < σ < 1 for the Oseen iterative method is derived. Theoretical results and numerical simulations show that the iterative algorithm for the Stokes iteration is the simplest for larger viscosity coefficients, and the Newton iteration is the most efficient method for relatively large viscosity coefficients. The Oseen iterative method is the most flexible method because its stability and convergence condition is the weakest. Furthermore, the simplified bearing lubrication problem with inhomogeneous boundary conditions is considered and some interesting physical phenomena are observed. In addition, we notice that the results of pressure in the third example are not ideal. In future research, we will consider to improve pressure stability by adding stability terms.
Author Contributions: Formal analysis, X.X. and D.L.; Methodology, X.X. and D.L. All authors have read and agreed to the published version of the manuscript.

Acknowledgments:
The authors sincerely thank the editor and referees for their valuable comments and suggestions which helped us to improve the quality of this paper.

Conflicts of Interest:
The authors declare no conflict of interest.