Two-Level Finite Element Approximation for Oseen Viscoelastic Fluid Flow

In this paper, a two-level finite element method for Oseen viscoelastic fluid flow obeying an Oldroyd-B type constitutive law is presented. With the newly proposed algorithm, solving a large system of the constitutive equations will not be much more complex than the solution of one linearized equation. The viscoelastic fluid flow constitutive equation consists of nonlinear terms, which are linearized by taking a known velocity b(x), and transforms into the Oseen viscoelastic fluid flow model. Since Oseen viscoelastic fluid flow is already linear, we use a two-level method with a new technique. The two-level approach is consistent and efficient to study the coupled system which contains nonlinear terms. In the first step, the solution on the coarse grid is derived, and the result is used to determine the solution on the fine mesh in the second step. The decoupling algorithm takes two steps to solve a linear system on the fine mesh. The stability of the algorithm is derived for the temporal discretization and obtains the desired error bound. Two numerical experiments are executed to show the accuracy of the theoretical analysis. The approximations of the stress tensor, velocity vector, and pressure field are P1-discontinuous, P2-continuous and P1-continuous finite elements respectively.


Introduction
In nature, most fluids are non-Newtonian, and many researchers have investigated the behaviour of non-Newtonian fluid extensively.Non-Newtonian fluid has numerous applications in different sectors such as in the natural substances, in food processing, and biological fluid.Some familiar examples of the natural substances are magma, lava, gums and so on.In the food processing area: butter, cheese, jam, ketchup, soap, yogurt, and in the biological fluid flow sector: blood, saliva, mucus, synovial fluids, are all examples of non-Newtonian fluid.Viscoelastic fluid is one of the essential non-Newtonian fluid which satisfies both the viscous property of the fluid and elastic estate of solid.Over the last few decades, the developments of viscoelastic fluid research have achieved significant progress, but the rheological properties of the viscoelastic fluid is not similar as Newtonian fluid.To describe the viscoelastic behaviour of the liquid, it was a significant challenge to formulate a suitable constitutive equation.In 1950, James G. Oldroyd [1] introduced a new paradigm to study the dilute solution of polymeric molecular behaviour.After that, many models have been developed to study the viscoelastic fluid such as Phan-Thien-Tanner, Larson, Maxwell model and so on.
The Oseen fluid flow model for Newtonian fluid can be obtained from the Navier-Stokes equation.By fixing the velocity, u(x) = b(x), of the nonlinear term of the momentum equation, the Navier-Stokes equation can be transformed into a linear system.In the viscoelastic fluid flow model, under a creeping flow assumption, the inertia term u • ∇u is neglected in the momentum equation.So in the viscoelastic fluid flow model, the non-linearity occurs only in the constitutive equation [2].Thus, no simplification takes place in the momentum equation, and the whole model equation consists of three unknowns which are: scalar pressure field p, velocity vector field u, and symmetric stress tensor σ.
There are two usual approaches to solve the viscoelastic fluid flow models: the Streamline Upwind Petrov Galerkin (SUPG) method and the discontinuous Galerkin (DG) approximation.In 1973, Reed and Hill first proposed the discontinuous Galerkin technique in [3].To solve the neutron transport equation, Lesaint and Raviart discussed the discontinuous method in [4] for hyperbolic PDEs.DG methods have become popular due to its computational flexibility, ability to incorporate physical properties, element-wise conservative and implementable on an unstructured mesh.Fortin and Fortin [5,6] first introduced the DG method for viscoelastic fluid flow.Barrenger and Sandri [7] first proposed error analysis for the steady-state case for viscoelastic fluid flow model using the DG method.For the steady-state problem, Najib and Sandri [8] introduced the semi-decoupled scheme by using fixed point technique.Also, some researchers used the Streamline Upwind Petrov Galerkin (SUPG) approximation [9] to deal with the constitutive equation.For example, in [10,11] Ervin and Miles analysed an implicit Euler time discretization and a SUPG discretization for the constitutive equation.To solve Oseen viscoelastic fluid flow, Lee et al. used domain decomposition method in [12], the defect correction process at high Weissenberg number in [13], two-level stabilized mixed finite element method in [14], stabilized Lagrange-Galerkin method for the nonlinear scheme in [15].To discretize the elliptic problems, the two-level method was first introduced by J. Xu in [16].Layton-Tobiska [17] used this method for solving Navier-Stokes equation and the same technique is studied in [18][19][20][21][22][23][24][25][26].
In this paper, we consider the two-level method to investigate the Oseen viscoelastic fluid flow for error estimation.To solve a multi-model, the two-level approach is familiar with its computational flexibility.Since in the Oseen viscoelastic fluid flow model, the constitutive equation is already linear where the non-linearity vanishes because of creeping flow.This new feature allows us to consider a two-level approach differently.In the first step, we solve the momentum equation and constitutive equation on a coarse mesh, and in the second step, we decoupled the momentum equation and constitutive equation with respect to velocity and stress.Using the coarse mesh solution, we determine the solution on the fine mesh.The implementation in the temporal discretization of the proposed two-level approach is as follows: we use two families of partitions T h and T H (h << H) to subdivide our domain into triangles with sides of length h and H respectively and the corresponding finite element spaces are denoted by (X h , Q h , S h ) and (X H , Q H , S H ). The well-posedness of the scheme is derived and obtains the desired error bound.The optimal convergence order is verified by considering a true solution.The exclusive feature of the two-level method for Oseen viscoelastic fluid flow is illustrated experimentally by considering a 4:1 abrupt contraction channel flow.
The remaining part of the paper is organized as follows: in Section 2, the Oseen viscoelastic fluid flow model is introduced.In Section 3, the finite element approximation and the variational formulation is discussed.We propose the algorithm for the two-level method, investigate the existence and uniqueness of finite element solution also perform the error analysis in Section 4. The results of the two numerical simulations illustrate in Section 5. Finally, a short conclusion is discussed in Section 6.

Model Equations
In this section, we describe the Oseen viscoelastic fluid flow model.

Model Problem
Let Ω be an open, bounded, simply connected domain in R 2 with the Lipschitz continuous boundary Γ.Consider the Johnson-Segalman model problem where σ denotes the polymeric stress tensor, u is the velocity vector, p is the scalar pressure field and λ is the Weissenberg number defined as the ratio of average relaxation time of polymer to an external given time which is typical length to a typical velocity in the flow.The total stress tensor + σ, consists of Newtonian part and viscoelastic part [27].For notational convenience, we use Einstein's convention of summation and denote differentiation with comma as ∂u ∂x i is written as u ,i and ∂u ∂t is written u t .For a scalar function p, gradient of p is a vector ∇p with (∇p) i = p ,i ; for a vector u, gradient of u is a tensor ∇u with (∇u) ij = u i,j .For a vector function u, [28].We assume that p has zero mean value over Ω. α is a number such that 0 < α < 1, which may be considered as the fraction of viscoelastic viscosity, and f is the body force.The deformation tensor, D(u), and the vorticity tensor, W(u), are given by where g a (σ, ∇u) is defined by We make the assumption for b:

The Weak Derivative and Finite Element Discretization
We introduce some notation first.The L 2 (Ω) and L p (Ω) norms are denoted by the inner product (•, •), and • L p , with the special cases of L 2 (Ω) and L ∞ (Ω) norms being written as • and • ∞ .For m ∈ N, we denote the norm associated with the Sobolev space W m,p (Ω) by • W m,p , with the special case W m,2 (Ω) being written as H m (Ω) with the norm • m and seminorm | • | m .In order to introduce a variational formulation, we set the spaces X, Q, S, V for velocity u, pressure p, and stress σ respectively.Velocity Space : In order to formulate the corresponding weak form of ( 7)- (10), multiply the Equations ( 7), ( 8) and ( 9) by the test function τ, v, and q respectively.Integrating by parts and applying the divergence theorem gives us: To find (u, It is well known that the velocity and pressure spaces, X and Q, satisfy the inf-sup (or LBB) condition [29] inf where C is a positive constant independent of h and H. Applying the inf − sup condition in the weak form (11)- (13) gives us, find (u, σ) : (σ, Let T µ is a uniformly regular triangulation of Ω such that Ω = {∪K : K ∈ T µ }, µ = h, H and µ = max K∈T µ µ K .Assume that there exist positive constants c 1 , c 2 such that c 1 µ µ K c 2 ρ K , where µ K is the diameter of K, ρ K is the diameter of the greatest ball included in K.The classical Taylor-Hood FE are used for the approximation in space of (X, Q) : P2-continuous in velocity u, P1-continuous in pressure p, and P1-discontinuous approximation is considered for the stress σ.Due to the hyperbolic nature of the constitutive equation, a stabilization technique is needed for the finite element simulation of viscoelastic flows.Streamline upwinding and discontinuous Galerkin method are the commonly used discretization techniques to handle this problem.We use the discontinuous Galerkin method for approximating the stress.The corresponding FE spaces are defined as where P i (K), i = 1, 2, denotes the space of polynomials of degree ≤ i on K µ ∈ T µ .It is well known [29] that the Taylor-Hood pair (X µ , Q µ ) satisfies the discrete inf-sup condition.For discontinuous stress, we need to use an upwinding technique introduced in [7]: for this, we define where ∂K µ is the boundary of K µ ∈ T µ and n is the outward unit normal to ∂K µ , and The term ((b • ∇)σ, τ) is approximated by means of an operator B µ on (X µ , S µ , S µ ) is defined by The discontinuous Galerkin finite element approximation of ( 11)-( 13) is as follows: to find Applying the inf-sup condition For error estimates, the standard approximation results are as follows: Let σh and ũh be the L 2 projections of σ ∈ S h , and and on Γ: The internal boundaries are: summing over all K and defining Γ h := {∪∂K, K ∈ T h } \ Γ, i.e. the internal boundaries yields:

Two-Level Algorithm's Existence and Uniqueness, Error Analysis of Oseen Viscoelastic Fluid Flow Model
In this section, we recall the two-level method for the steady state viscoelastic fluid flow and propose the two-grid method for Oseen viscoelastic fluid flow model.We analyse the existence and uniqueness of the model problem and derive the error analysis of the two-level process of Oseen viscoelastic fluid flow model.

Two-Level Method for Steady State Viscoelastic Fluid Flow Model
First, we recall the traditional two-level method for steady state viscoelastic fluid flow model.
Step 1.To estimate the solution on nonlinear coarse mesh: Step 2. To estimate the solution on linear fine mesh:

Two-Level Method for Oseen Viscoelastic Fluid Flow
Step 1.To find the solution on nonlinear coarse mesh of the problem: Step 2. To find the solution on linear fine mesh of the problem: (σ For 0 < α < 1, there exist C 0 , C 0 and H 0 such that if (1)-( 4) admits a solution (σ, ) and there exists a constant C, independent of H, such that: Furthermore, there exists an open set ε such that (σ For M satisfying 1 − 2λMd > 0, and f ∈ H −1 (Ω), there exists a unique solution (σ h , u h ) ∈ (S h × X h ) of the Equations (38)-(40).
Proof of Theorem 2. Multiplying the Equation (42) by 2α and adding with (41) as follows Now we will show that the left hand side of (44) is continuous and coercive on Hence, which shows the continuity of the two-level algorithm.For coercivity, using (18) we have which shows the coercivity of the two-level algorithm.Hence, the well-posedness of the two-level finite element scheme for Oseen viscoelastic fluid flow is obtained.

Error Analysis
In this subsection, we derive the error analysis for the solution (σ µ , u µ ) of the proposed two-level-method.Theorem 3. If (u h , σ h ) satisfies a unique solution of step 2 of the two-level methods of Oseen viscoelasic fluid then the following estimate holds Proof of Theorem 3. To find (σ, u) ∈ (S h × V h ) satisfies for all (τ h , v h ) ∈ (S h × V h ) as follows: The inf − sup condition in step 2 of the algorithm can be written as: for all and Subtracting ( 53) from ( 51) and ( 54) from (52) as follows: and Adding and subtracting σh , ũh in the Equations ( 55) and (56) which yields Also, Since σh and ũh are orthogonal projection of σ and u also, Then Equation ( 58) becomes Now let τ h = σh − σ h , v h = ũh − u h in Equations (57) and (59) respectively.From Equation (57 and from the Equation (59) we get By using Equation ( 60) and ( 18) we get where using Young's inequality Now combining all the estimates of Σ Substituting Equation (69) in Equation ( 62) By the triangle inequality, we have From Equation (61) Using the triangle inequality Using the result of ( 24), Theorem 1 and substitution of (75) in ( 72) And Hence, Choose appropriately, we have which completes the proof of the error analysis of the two-level method for Oseen viscoelastic fluid flow.
Remark 1.The error estimate of the two-level method for Oseen viscoelastic fluid flow presented in Theorem 3 is a rough error, and the estimation of optimal convergence order is an open question.The scaling of h = H 3/2 guarantees the linear convergence.In our understanding, the two-level method for viscoelastic fluid flow model is not optimal, which is reflected by the papers [24,25].In both articles, authors didn't get the optimal error order for the multi-level and two-level method for viscoelastic fluid flow model and the consequences for Oseen viscoelastic fluid flow model have the same phenomenon.In the computational sense, the two-level method is flexible, and the implementation is easy.

Numerical Tests
In this section, two numerical experiment results are presented to illustrate the exclusive features of the proposed two-level method for Oseen viscoelastic fluid flow.First numerical simulation is known as an analytical solution test.A well-known approach for the validation of the theoretical analysis is to show the computation of the errors and the order of convergence for the proposed scheme.In the exact solution test, the computation of the errors and the order of convergence are presented for one-level and two-level algorithms with P2 − P1 − P1dc finite elements for Oseen viscoelastic fluid flow.The second numerical test is known as the benchmark problem called 4:1 abrupt contraction channel flow, a prototypical problem to show the viscoelastic flow behaviour [30,31].In 4:1 contraction channel flow, we demonstrate the graphical representation of streamlines, pressure oscillation, and reentrant corner flow behaviour of the steady-state viscoelastic fluid flow model and Oseen viscoelastic fluid flow model for one-level and two-level algorithms.The computations carried out for the solutions under the fine mesh converge to a solution computed using the standard coarse mesh.The one-level and two-level methods for the viscoelastic fluid flow model are well studied where the flow behaviour and patterns of the streamlines are regular in shape.We compare the behaviour of the streamlines and patterns of the contour for the one-level and two-level methods for Oseen viscoelastic with the stationary viscoelastic fluid flow model which appear in the similar precision.The comparison between these two models and two methods give us a strong overview of the proposed two-level schemes accuracy.According to the theoretical analysis, continuous piecewise quadratic elements P2 are used for velocity, continuous piecewise linear elements P1 are used for pressure, and discontinuous piecewise linear elements P1dc are used for stress.All the numerical tests are performed by using the public domain software Freefem++ [32].The figures and graphs are drawn by MATLAB and Tecplot 360 software package.

Analytic Solution Test
The computational domain for analytical solution is considered as a unit square Ω = [0, 1] × [0, 1], and the parameters λ, α, a in the equations are specified as 1.0, 0.5 and 0 [25,33,34], respectively.Hence, the true solution is given by [33] In the numerical computation b(x) is considered as an exact solution of u.Considering the global domain Ω as a unit square which is divided into identical squares to create the triangular meshes.
For convenience the denotations of one-level errors are as follows The denotations for two-level errors are as follows The standard way to verify the theoretical analysis accuracy is the computation of errors and the order of convergence with the analytical solution.The order of convergence is represented by Order = log(Error/Error1) log(h/h ) .The term "Order" denotes the global rate of convergence.The notation h and h denote the mesh size with the global errors "Error" and "Error1".The objective of the representation of the following two tables is to verify the theoretical analysis accuracy by computing the order of the convergence for the one-level and two-level method.
In the one-level method for Oseen viscoelastic fluid flow model, we observe that the convergence order for velocity u, stress σ, and pressure p achieve more than 2 order for L 2 -norm while velocity in H 1 -norm obtained convergence order more than 1.5.For the two-level method, the convergence order for velocity u in L 2 -norm is nearly second order, but for stress σ, the order is nearly 1.5 for L 2 -norm.The convergence order for pressure p in L 2 -norm is nearly 1.5 order.For velocity u in H 1 -norm the convergence order is 1.5.The results also illustrate for one-level, and two-level method for Oseen viscoelastic fluid flow with P2 − P1 − P1dc finite element, achieve the desired error estimation with exact convergence order.The computation of errors and order of convergence validate the theoretical analysis accuracy for the relation H ≈ h 2/3 .Moreover with increases of the value of λ has a small effect in the error order appear in a similar precision for the one-level and two-level method.

4:1 Contraction Channel Flow
In this subsection, we examine a benchmark problem for viscoelastic fluid flow which is known as 4:1 contraction channel flow.The study of this classic reference problem is used to investigate the performance of various two-dimensional finite element, finite volume, hybrid, and spectral viscoelastic solvers.It has various interesting applications in polymeric industries.The fluid behaviour of contraction channel flow was studied widely by many authors in [35][36][37], also extensively used to show the behaviour of the streamlines and contour patterns.Here the case of planar flow through a contraction geometry with a ratio of 4:1 concerning upstream and downstream channel widths are considered.The domain is constructed such that the channel lengths are sufficiently long for fully developed Poiseuille flow at both the inflow and outflow boundaries.We demonstrate several graphical representations of viscoelastic fluid flow model and Oseen viscoelastic fluid flow behaviour for one-level and two-level method with reentrant corner flow, streamlines behaviour, pressure gradient flow behaviour.The main focus of this experiments is to illustrate the proposed algorithms efficiency, and accuracy is comparing with the well studied steady-state viscoelastic fluid flow model's behaviour of streamlines and flow patterns.
For 4:1 contraction channel flow of Oseen viscoelastic fluid, we first determine the one-level solution of u 1 and u 2 associated with one-level contraction domain.For two-level contraction domain, the solution of uh 1 and uh 2 are determined in a similar way.We use those solutions in Oseen one-level and two-level computation as b 1 , b 2 , bh 1 and bh 2 , respectively.The computations are performed on a uniformly refined mesh, shown in Figure 1.The method of structured mesh constructs the domain for one-level H and for two-level h = H 3/2 .The one-level and two-level domain geometry consist of the vertices (0, 0), (8, 0), (8, 0.25), (4, 0.25), (4, 1), (0, 1).The inflow and outflow boundaries are represented by ∂Γ in = {(x, y) : x = 0, 0 y 1} and ∂Γ out = {(x, y) : x = 8, 0 y 0.25} for both one-level and two-level method [27,30,31].For the velocity, For stress, on ∂Γ in , No slip boundary conditions are imposed for the velocity on the solid walls of the contraction, and symmetry condition is imposed on the bottom of the computational domain.Besides, the physical parameters Re, α, λ, and a are chosen as 1, 8/9, 0.7 and 1, respectively.The computation carried out using P2 − P1 − P1dc finite elements for one-level and two-level method.The size of mesh for the 4:1 contraction channel are 0.448783.For the one-level method, number of elements are 1408, number of nodes for P2 is 2945, P1 is 769, P1dc is 4224.For the two-level method, number of elements are 12,496, number of nodes for P2 is 25,373, P1 is 6439, P1dc is 37,488.Figure 2 illustrates the graphical representation of the horizontal and vertical velocity near re-entrant corner along the vertical line x = 4.0265 [27].We demonstrate a one-level and two-level method for steady state viscoelastic and Oseen viscoelastic fluid flow model.From the figure, we observe that the velocity near re-entrant corner is almost continuous for one-level and two-level method.The graphical representation illustrates the flow behaviour of the horizontal velocity near re-entrant corner for one-level and two-level method is almost continuous for Oseen viscoelastic fluid and similar to the viscoelastic fluid flow.For the vertical velocity, the flow behaviour achieves a little different flow path but appear continuously.The graphical representation of the velocity near re-entrant corner shows the accuracy of the two-level scheme.
In the Figures 3 and 4, we represent the streamlines and magnitude for steady state viscoelastic fluid flow and Oseen viscoelastic fluid flow for one-level and two-level algorithms.The figures demonstrate that the flow enters through the inflow boundary ∂Γ in = {(x, y) : x = 0, 0 y 1} and flow out through the outflow boundary ∂Γ out = {(x, y) : x = 8, 0 y 0.25}.As expected, the flow behaviour of the viscoelastic type always creates a vortex in the corner [27,34].Oseen viscoelastic fluid flow is a reduced linear model of the viscoelastic fluid flow model.The accuracy, behaviour of the streamlines, flow characteristics and the shape of the contour is examined by the Figure 3.In Figure 3, the flow pattern and nature of the streamlines are presented for the steady-state viscoelastic model and Oseen viscoelastic fluid flow model with the one-level algorithm.The contour is drawn by using the values generated for v.The comparison of both figures for one-level method illustrates that the reduced Oseen viscoelastic fluid flow model's shape of the contour is regular and obtains similar accuracy.The flow rate and behaviour of the streamlines are also almost same.In the Figure 4, the fluid flow behaviour of the two-level algorithm is presented for the steady-state viscoelastic model and Oseen viscoelastic fluid flow model.The contour is drawn for the two-level method with uh where the red colour in the outflow channel shows the maximum flow speed.In the inflow boundary, the flow speed is less, which is represented by the contour colour and the slower flow is in the corner section where a vortex form.The shape of contour, the behaviour of the streamlines and flow rate obtain similar accuracy for the two-level method for steady state viscoelastic fluid flow and Ossen viscoelastic fluid flow.All the figures appear with similar precision, and the patterns of the streamlines are regular which validate the numerical computation's accuracy.The physically valid results give the relevance of the numerical methods.
In Figure 5, we represent the pressure contour for one-level and two-level methods for Oseen viscoelastic fluid.The figures reveal the similar pattern and no pressure oscillation which also confirm the two-level process accuracy.

Conclusions
In this contribution, the two-level decoupled method for Oseen viscoelastic fluid flow is presented and implemented for theoretical analysis with numerical experiments.Oseen viscoelastic fluid flow model is a reduced linear model of the viscoelastic fluid flow model.The proposed numerical scheme is reasonable for Oseen type problem where the non-linearity vanishes.To solve a multi-model, two-level approach is familiar with its computational flexibility.The method developed herein reduces the computational cost and implementation is easy.The solution on the coarse grid is derived in the first step, and the result is used to determine the solution on the fine mesh in the second step.The well-posedness of the scheme ensures the validation of the temporal discretization.The scaling of h = H 3/2 guarantees the linear convergence and the optimality of the convergence order is an open question.The numerical experiments support the theoretical analysis accuracy.The optimal convergence order is obtained by considering an analytical solution of the model problem.The graphs and figures of the viscoelastic fluid flow model and Oseen viscoelastic fluid flow model for the one-level and two-level methods appear in similar precision, which also validate the numerical ) a is related to material parameter with a ∈ [−1, 1].For the case, a = 1, the Johnson-Segalman model reduces to the well-known Oldroyd-B model.The steady state viscoelastic fluid flow model (1)-(4) can be transformed into Oseen viscoelastic fluid flow model considering a known velocity b(x) in the nonlinear terms of the constitutive equation instead of unknown velocity u(x) to make the constitutive equation linear.The Oseen viscoelastic fluid flow model as follows:

) 4 . 3 .
Existence and Uniqueness of the Finite Element Solution Theorem 1. (Existence and uniqueness of the steady state viscoelastic fluid flow model step 1)

Figure 1 .
Figure 1.The geometry of contraction mesh for one-level (up) and two-level (down).

Figure 2 .Figure 3 .
Figure 2. Horizontal velocity (left), vertical velocity (right) near re-entrant corner for steady state viscoelastic and Oseen viscoelastic fluid flow one-level and two-level method.

Figure 4 .
Figure 4. Streamlines and magnitude of velocity contours u for two-level method: steady state viscoelastic fluid flow model (left), Oseen viscoelastic fluid flow model (right).

Figure 5 .
Figure 5. Representation of pressure contour p: for one-level (left) and two-level (right) method for Oseen viscoelastic fluid flow.