Application of the Fictitious Domain Method for Navier-Stokes Equations

Almas Temirbekov,Zhadra Zhaksylykova,Yerzhan Malgazhdarov and Syrym Kasenov

1Al-Farabi Kazakh National University,Almaty,Kazakhstan

2Abai Kazakh National Pedagogical University,Almaty,Kazakhstan

3Sarsen Amanzholov East Kazakhstan University,Ust-Kamenogorsk,Kazakhstan

Abstract:To apply the fictitious domain method and conduct numerical experiments,a boundary value problem for an ordinary differential equation is considered.The results of numerical calculations for different values of the iterative parameter τ and the small parameter ε are presented.A study of the auxiliary problem of the fictitious domain method for Navier-Stokes equations with continuation into a fictitious subdomain by higher coefficients with a small parameter is carried out.A generalized solution of the auxiliary problem of the fictitious domain method with continuation by higher coefficients with a small parameter is determined.After all the above mathematical studies,a computational algorithm has been developed for the numerical solution of the problem.Two methods were used to solve the problem numerically.The first variant is the fictitious domain method associated with the modification of nonlinear terms in a fictitious subdomain.The model problem shows the effectiveness of using such a modification.The proposed version of the method is used to solve two problems at once that arise while numerically solving systems of Navier-Stokes equations:the problem of a curved boundary of an arbitrary domain and the problem of absence of a boundary condition for pressure in physical formulation of the internal flow problem.The main advantage of this method is its universality in development of computer programs.The second method used calculation on a uniform grid inside the area.When numerically implementing the solution on a uniform grid inside the domain,using this method it’s possible to accurately take into account the boundaries of the curved domain and ensure the accuracy of the value of the function at the boundaries of the domain.Methodical calculations were carried out,the results of numerical calculations were obtained.When conducting numerical experiments in both cases,quantitative and qualitative indicators of numerical results coincide.

Keywords:Fictitious domain method;Navier-Stokes equations;difference schemes;approximation;computational algorithm;numerical experiment

1 Introduction

At the beginnig we will give a literature review on the application of the fictitious domain method.Currently,there are several methods for the numerical solution of boundary value problems in complex geometric domains,such as the method of curved grids and the fictitious domain method.The construction of curved grids for the numerical solution of problems requires the transformation of the equation into curved coordinates,which has a more complex form than the original equations.When constructing curved grids,various requirements are imposed on difference equations,which makes the construction of curved grids a difficult mathematical task.Therefore,for the numerical solution of a wide class of problems of mathematical physics in an arbitrary domain,it is effective to use the fictitious domain method[1,2].The fundamental works on the application of the fictitious domain method include the works of Vabischevich et al.[1-8].Works[1,2]are devoted to the fictitious domain method in the numerical solution of problems of mathematical physics in complex areas.The fictitious domain method is based on the transition to a problem in a regular area that entirely contains the original one.The issues of substantiation of such an approach at the differential level in the study of boundary value problems for elliptic and parabolic equations,eigenvalue problems are considered.Modifications of well-known iterative methods are constructed to solve grid problems that arise when using the fictitious domain method.The possibilities of the fictitious domain method are illustrated by examples of solving problems of ideal and viscous incompressible fluid,filtration under a hydraulic structure.In[3,4],the fictitious domain method is applied to elliptic equations.Error estimates are obtained for solving problems using the fictitious domain method.Lagrange multipliers based on methods of fictitious domains are applied.Papers[5-8]are devoted to the application of the fictitious domain method for solving computational fluid dynamics problems.

The fictitious domain method for the Navier-Stokes equations is the subject of works by Bugrov et al.[9].Smagulov et al.[10]proposed the fictitious domain method for hydrodynamic equations in multiply connected domains.

In[11],a variation of the fictitious domain method for the Navier-Stokes equations for a viscous incompressible fluid in the velocity and pressure variables was proposed.Moreover,a constant value was imposed on pressure at the boundary of the auxiliary domain.This condition allows obtaining the Dirichlet problem for the Poisson equation for pressure which admits the development and application of effective iterative numerical methods for its solution.

The fictitious domain method is utilized in solving problems of computational fluid dynamics such as modeling the motion of particles in a fluid flow[12],the motion of an incompressible fluid[13],the two-phase Stokes problem including the surface tension force[14],the problem of flow around moving or deformable bodies fluid flow[15-17]and others.

In[18],the spectral element/HP method was applied to perform direct numerical simulation(DNS)for one row of turbine blades T106A using open source software Nektar++.The main goal of the current research is to perform preliminary researches for uniform,stable flow around an aerodynamic profile using the Nektar++solver for 2D Navier-Stokes equations for incompressible flow.

In[19],the fictitious domain method was applied to simulate the interaction of liquid particles with the Navier slip boundary condition.Numerical experiments show that the method works well for an anisotropic particle in a flow with the Navier slip boundary condition.

In[20],fictitious domain method with aH1-penalty for the Stokes problem with a Dirichlet boundary condition is studied.Work[21]is devoted to the application of the fictitious domain method in the numerical simulation of a pulse oscillation converter.In[22],fictitious domain method with a distributed Lagrange multiplier was studied for parabolic problems of a jump-like type with moving boundaries.

In[23],the influence of rotation and hydrostatic initial stresses under the influence of an electromagnetic field incident on the outer surface of a semiconductor medium was investigated.

Papers[24-27]are devoted to the study of fictitious domain method for problems with discontinuous coefficients.Elliptic type equations with strongly varying coefficients are considered and investigated in[26-27].Equations of this type are obtained using fictitious domain method.A special method is used for the numerical solution of an elliptic equation with strongly varying coefficients.For the obtained problem,the convergence rate estimation theorem of the developed numerical algorithm is proved.Based on obtained estimates,a computational algorithm was developed and numerical calculations were performed to illustrate the effectiveness of the proposed method.Work[28]is devoted to the development of new exact and numerical solutions for the 1-dimensional integro-differential Ito equation using the methods of 1-expansion and finite differences,respectively.Trigonometric,hyperbolic and rational solutions are successfully presented.The stability and accuracy of the obtained numerical simulation are discussed.The presented graphical comparison shows that the exact and numerical solutions almost coincide with each other.Analytical and numerical solutions of the generalized Benjamin-Bona-Mahoney equation(GBBM)are investigated in[29].The exact solution is obtained analytically,while numerical solutions are demonstrated using some methods,namely adaptive moving mesh and uniform mesh methods.The exact solution is presented in the form of convergent power series.Finite differences are also used to discretize the BBM equation.

The works[30-32]are devoted to the application of a parallel computational algorithm for fictitious domain method.In[30],a parallel computational algorithm is constructed using the fictitious domain method for the three-dimensional Helmholtz equation.Work[31]is devoted to modeling turbulent flow in a channel using fictitious domain method.The paper[32]describes the possibilities of using the fictitious domain method for biomechanics problems.

The work[33]is aimed at studying the effect of rotation on the general model of generalized thermo-microstriction equations for a homogeneous isotropic elastic semi-spatial solid whose surface is subjected to thermal shock.Comparisons are made with the results in the presence and absence of rotation,as well as in the presence and absence of microextension constants between the two theories.

In[34],the influence of variable thermal conductivity,which depends on temperature,is considered in the context of photothermal diffusion(PTD).The PTD process is applied using the theory of thermoelasticity under chemical action.The presented model describes the interaction between elasticthermal-plasma waves based on the properties of the material of a semiconductor elastic medium.The Laplace transform is used to solve control equations in one dimension of a thin circular plate.Complete solutions in the time domain are observed using the numerical approximation method.Physical fields with some comparisons are presented analytically and graphically.

The paper[35]presents the results of a study of the Navier-Stokes equations in numerical modeling in two-connected domains.Two methods are used to solve this equation numerically.In the first method,the condition of unambiguity of pressure is used for numerical solutions of difference equations for the functions of current and velocity vortex.To obtain a numerical solution of the elliptic equation,which is obtained for current functions is found as the sum of two problems.The second alternative method for solving the difference problem,fictitious domain method is considered.In this method,there is no need to satisfy the condition of one-digit pressure,so it is easy to implement.

In[36],the augmented domain method is considered for numerical simulation of the flow of a viscous incompressible fluid in complex geometric domains.The problem is considered in a discretely given two-connected domain with a curved boundary.Spline interpolation of curved boundaries is performed.A monotonic finite difference scheme and an algorithm for numerical realization of the Navier-Stokes equation for a viscous incompressible fluid are developed.Numerical results are obtained for different numbers of grid nodes.

Glowinski et al.[5,37-39]considers a family of domain decomposition methods based on the explicit use of the Lagrange multiplier defined on the actual boundary.The proposed technique associated with true boundary conditions is common for modeling inviscid incompressible potential flows.According to the proposed method,the original differential problem is reduced to an optimal control problem with a saddle point,and the iterative method of conjugate gradients is used for its numerical implementation.

In[40],the fictitious domain method without Lagrange multipliers is used for numerical simulation of the suspension concentration.The boundary is tracked to delimit the domains occupied by liquid and solid particles.The particle indicator function is constructed using the Heaviside function.Further,the Heaviside function is approximated using the hyperbolic tangent in a small neighborhood(proportional to the grid step) around the boundary.Such a continuous function ensures smooth change,improves numerical accuracy and reliability.

Let us dwell on the advantages and disadvantages of the existing variants of the fictitious domain methods.The method proposed and implemented in[5,37-39]allows satisfying the boundary conditions on the actual boundary using variational principles.Two meshes are constructed for the finite element solution in the extended domain:triangulation over the entire domain and a curvilinear mesh on the actual boundary.Then,the desired variables are matched on the curvilinear boundary at each iteration of the conjugate gradient method.In some cases,this process may behave worse than the solution of the original problem on a non-uniform mesh consistent with a curvilinear boundary.The results of numerical calculations show convergence to the solution“in the average”due to the use of variational principles since the functional containing the main equation and the boundary condition is minimized.

Therefore,in this paper,the authors propose a method that allows constructing a homogeneous difference scheme in the entire extended domain which is a convenient tool in terms of programming automation.At the same time,a reasonable continuation of the main equation coefficients leads to the convergence of the solution to the desired solution in the original domain,which is confirmed by mathematically proven statements and the results of numerical calculations.

In this work,two methods are applied for the numerical solution of the formulated problem.The first one is the fictitious domain method associated with the modification of nonlinear terms in a fictitious subdomain.The model problem shows the efficiency of using such a modification.The proposed variation of the method is used to solve two problems at once that arise in the numerical solution of the Navier-Stokes equations:the problem of the curvilinear boundary of an arbitrary domain and the problem of the absence of a boundary condition for pressure in the physical formulation of the internal flow problem.

The second method used the calculation on a uniform mesh inside the domain.The solution on a uniform mesh inside the domain makes it possible to accurately take into account the boundaries of the curved domain and ensures the accuracy of the function value on the boundaries of the domain in the numerical implementation.

2 One-Dimensional Problem

Fig.1 and Tab.1 show the results of numerical calculations for various values of the iterative parameterτand small parameterε.Methodological calculations were carried out using two schemes.The calculation results demonstrate good convergence of the auxiliary problem solution of the fictitious domain method to the solution of the main problem in the original domain.In this case,the maximum deviation is 0,00003472.

Table 1:Results for different values of the parameter ε

3 Multivariable Problem

Formulation of the problem.Let us consider an initial-boundary value problem for an unsteady flow of a viscous incompressible fluid in a bounded domain Ω ⊂R2with a curved borderS.The problem is reduced to solving the system of nonlinear Navier-Stokes equations in the velocity-pressure variables:

4 Numerical Calculations

We consider the numerical solution of the auxiliary problem(13)-(16)to illustrate the advantages of the proposed approach.

Let us take a curvilinear channel with solid boundaries(Fig.2)as the domain under consideration.D1denotes the fictitious domain and Ω denotes the physical domain.

We use the finite difference method and the scheme of splitting by physical processes[42]for the numerical implementation of the auxiliary problem(13)-(16).The integration domain is covered by the so-called MAC-grid[43],in which the nodes for determining the pressure are located inside a rectangular grid cell,and the nodes for determining the velocity component are located on its faces.The scheme taking into account the sign is used when approximating the convective terms.Thus,an approximation with the second order of accuracy in space and the first order in time is provided.

The prescribed values of pressure and zero values of the tangential component of the fluid flow velocity at the inlet and outlet of the computational domain were set in the numerical implementation.At the «solid» boundaries,the pressure values are given as linear functions and the tangential component of the velocity is equal to zero.

Obviously,the sum of the equations corresponding to Stages I and III gives the original equation of motion(13),and Eq.(39)corresponding to Stage II is obtained from(40)by applying the divergence operator to the last one taking into account the continuity equation.

The following physical interpretation of the given splitting scheme is proposed.At Stage I,it is assumed that the transfer of momentum (momentum per unit mass) is carried out only due to convection and diffusion.The velocity field thus obtained does not satisfy the incompressibility condition,in general.In this paper,the implicit scheme is used at Stage I in contrast to the classical version of the method of splitting into physical processes[42].A similar approach was substantiated and numerically implemented in[44].Iterative schemes were proposed for solving auxiliary grid Navier-Stokes equations.This approach improves the condition number of the system of linear equations matrix and makes it possible to speed up the convergence of the solution.

At Stage II,the pressure field is found from the solution of the Poisson equation based on the found intermediate velocity field and taking into account the solenoidality condition of the velocity vector.At this stage,numerical methods for solving grid elliptic equations with Dirichlet boundary conditions for pressure were used.

At Stage III,it is assumed that the transfer is carried out only due to the pressure gradient,and the convection and diffusion are absent.

In the numerical implementation,grids containing 50×20,100×40,150×60,200×80 nodes and the following dimensionless parameter values were used:

0＜xi ＜2,0＜yi ＜1,τ=0.001,μ=1/Re=0.001,ε=10-5.

The upper and lower solid curved boundaries are described by the equations

y1=0.2-0.1 cos(2πx)

y2=0.8-0.1 cos(2πx)

Currently,there are several methods for the numerical solution of boundary value problems in complex geometric domains,such as the method of curved grids and fictitious domain method.The construction of curved grids for the numerical solution of problems requires the transformation of the equation into curved coordinates,which has a more complex form than the original equations.Therefore,for the numerical solution of a wide class of problems of mathematical physics in an arbitrary domain,it is effective to use the fictitious domain method.

The solution of this problem is implemented in two ways.In the first case,the fictitious domain method by the leading coefficients was used.The main advantage of this method is its versatility in the development of computer programs for the numerical simulation of a wide class of problems of mathematical physics.In the second case,a calculation on a was used.Methodological calculations have been carried out,and the results of numerical calculations have been obtained.The quantitative and qualitative indicators of the numerical results coincide for both cases when conducting numerical experiments.

The boundary of the physical domain is smeared when solving the problem by the fictitious domain method,and therefore the solutions at the boundaries may differ from the boundary condition,although it gives reliable results of the flow.In addition,the fictitious domain method is easily implemented.But since the problem is ill-conditioned at the first stage,an implicit scheme was used where the boundary conditions of the integer iteration step were used in the calculations.

Fig.3 shows the vector field of the problem solved using the fictitious domain method,respectively,Fig.4 shows the velocity module of the problem solved using the dummy domain method.Curved lines can be seen from the velocity modulus,Poiseuille currents are formed corresponding to the curvature of the boundary.

Fig.5 shows the vector field of the problem solved on a consistent grid,respectively,Fig.6 shows the velocity module of the problem solved on a consistent grid.In Figs.4 and 6,one can see the qualitative and quantitative correspondences of the formed Poiseuille flows when solving the problem by different methods.

Numerical implementation of the solution on a uniform grid inside the domain makes it possible to accurately take into account the boundaries of the curved domain and ensures the accuracy of the function value on the domain boundaries.The only drawback of this approach is the lack of a specific algorithm,and one has to come up with separate conditions for determining the boundary and boundary nodes of the grid for each problem,especially when boundary conditions of the second kind are imposed.

Tab.2 shows the maximum values of the velocity components on the indicated sections of the channel and their discrepancies.It follows from the table that there are small deviations which indicate a quantitative coincidence of the indicators.

Table 2:The value of the velocity components on the sections

Table 2:Continued

It can be seen from Figs.7-9 that the lower curved boundary relative to the computation domain is convex,and the upper one,on the contrary,is concave in the area of the section X=0.5 and X=1.5.This affects the profile values of the velocity component U at the lower boundary to a greater extent than at the upper boundary since they flow around the hill at the lower boundary and fall into the concavity pit at the upper boundary(Fig.5).The values of the velocity component V are positive everywhere since the fluid flow moves in a positive direction(upward).Also,it can be seen that the values of the velocity component U2 have negative values in the border areas of the upper boundary which indicates the formation of vortex motions in the cavities.

In Fig.8,on the contrary,the values of the velocity component U at the upper boundary are greater than at the lower one.Vortex motions are observed near the lower boundary from the presented graph of the velocity component U2.The values of the velocity component V are negative since the flow moves in a decreasing direction relative to the y axis(downward).The results presented in Fig.7 are similar to those obtained in Fig.5.

The calculations used a uniform grid with dimensions of 50 × 20,100 × 40,150 × 60,200 ×80.A numerical experiment was carried out on a modern personal computer with the following characteristics:Intel(R)Core(TM)i9-10900F CPU@2.80GHz,RAM 32 GB.

Thus,in order to apply the fictitious domain method and conduct numerical experiments,the boundary value problem for an ordinary differential equation is first considered.The results of numerical calculations for different values of the iterative parameterτand the small parameterεare presented.After successful application of the fictitious domain method for an ordinary differential equation,a more complex problem of applying the fictitious domain method for the Navier-Stokes equation in natural variables is considered.

Further,the research of the auxiliary FDM problem at the differential level for the Navier-Stokes equations with continuation into a fictitious subdomain by the higher coefficients with a small parameter is carried out.Methods of a priori estimates are used for the mathematical study of the problems under consideration.A generalized solution of the auxiliary FDM problem with continuation by higher coefficients with a small parameter is determined.The fictitious domain method is used to solve many problems of computational fluid dynamics.Currently,there are several methods for the numerical solution of boundary value problems in complex geometric domains,such as the method of curved grids and fictitious domain method.The construction of curved grids for the numerical solution of problems requires the transformation of the equation into curved coordinates,which has a more complex form than the original equations.Therefore,for the numerical solution of a wide class of problems of mathematical physics in an arbitrary domain,it is effective to use the fictitious domain method.

Thus,in this paper,two methods are applied for the numerical solution of the formulated problem.The first one is the fictitious domain method associated with the modification of nonlinear terms in a fictitious subdomain.The model problem shows the efficiency of using such a modification.The proposed variation of the method is used to solve two problems at once that arise in the numerical solution of the Navier-Stokes equations:the problem of the curvilinear boundary of an arbitrary domain and the problem of the absence of a boundary condition for pressure in the physical formulation of the internal flow problem.The main advantage of this method is its versatility in the development of computer programs.

The second method used the calculation on a uniform mesh inside the domain.The solution on a uniform mesh inside the domain makes it possible to accurately take into account the boundaries of the curved domain and ensures the accuracy of the function value on the boundaries of the domain in the numerical implementation.Tabs.3 and 4 present the computation time and number of iterations obtained with the use of the two methods for different grid nodes.It can be seen from the presented tables that the calculation time and the number of iterations are noticeably longer in case of a numerical solution on a uniform grid inside the domain.This is due to the fact that the conditions for belonging of the calculation node to the computational area are checked at each stage of the algorithm.This affects the counting time and the number of iterations.

Table 3:Counting time for different grid configurations

Table 4:Number of iterations for different grid nodes

Acknowledgement:Authors thank those who contributed to write this article and give some valuable comments.

Funding Statement:This research is funded by the Science Committee of the Ministry of Education and Science of the Republic of Kazakhstan(Grant No.AP09058430).

Conflicts of Interest:The authors declare that they have no conflicts of interest to report regarding the present study.