A method for dealing with viscosity terms in discontinuous finite element discretization based on internal penalty method

Abstract

A method for processing discontinuous finite element discrete viscosity terms based on an internal penalty method belongs to the field of computer-aided engineering. On the basis of the internal penalty method, the present invention converts the Dirichlet boundary into a Newman boundary by introducing a penalty coefficient when solving an elliptic equation with a Dirichlet boundary condition. It is extended to the boundary of an internal unit, and by introducing a bilinear operator as a penalty function, information exchange between adjacent units is established. Based on a high-order precision discontinuous finite element format, an improved internal penalty method is proposed, a new lifting operator is defined, and the tensor is eliminated, which is used to solve the compressible Navier-Stokes equations and numerically simulate fluid behavior. This method eliminates the derivative of the gradient of the viscosity term to the variable, and directly calculates the lifting operator using the expression of the viscosity term, thereby improving the computational efficiency; it can be conveniently extended to other types of equations, such as heat conduction equations, structural mechanics equations, and electromagnetic field equations.

Description

Technical Field

[0001] The invention belongs to the field of high-precision calculation of numerical simulation, and in particular relates to an improved internal penalty discontinuous finite element method for solving compressible fluid. Background Art

[0002] When the discontinuous finite element format was first designed, it could only solve hyperbolic conservation equations with first-order partial derivatives, such as the Euler equation. Subsequently, many scholars proposed different discretization schemes for problems with second-order derivatives or higher-order derivatives. When solving the compressible Navier-Stokes equations, Bassi and Rebay introduced variable gradients as intermediate variables to reduce the original second-order equations to two first-order systems, and then discretized them separately according to the discontinuous finite element framework. The disadvantage of this approach is that it increases the number of equations and has been proven to be incompatible. Another widely used viscous term discretization method is that Cockburn and Shu proposed to introduce variable gradients as intermediate variables, and when constructing variable fluxes and viscous fluxes on the unit interface, a stability coefficient is added to ensure the stability of the format.

[0003] Another discretization method for elliptic equations is called the interior penalty method (IP). However, when the IP method is used to calculate multivariable equations, such as the Navier-Stokes equations, it is necessary to calculate the Jacobian matrix of the viscosity term with respect to the gradient of the variable. This matrix is ​​a fourth-order tensor with a very complex form, which is not conducive to generalization to other equations. Summary of the invention

[0004] Based on the internal penalty method, the present invention converts the Dirichlet boundary into the Newman boundary by introducing the penalty coefficient when solving the elliptic equation with the Dirichlet boundary condition. It is extended to the boundary of the internal unit and establishes the information exchange between adjacent units by introducing a bilinear operator as the penalty function. Based on the high-order precision discontinuous finite element format, an improved internal penalty method is proposed, a new lifting operator is defined, and the tensor is eliminated, which is used to solve the compressible Navier-Stokes equations and numerically simulate the fluid behavior.

[0005] In order to solve the above technical problems, the specific technical solution of the method for processing the discontinuous finite element discrete viscosity term based on the internal penalty method of the present invention is as follows:

[0006] Step 1: Based on the compressible Navier-Stokes equations, it can be written as:

[0007]

[0008] in, and denote the convection term and the viscous term respectively, represents the time partial derivative, represents the spatial gradient operator, U is the conserved variable of the compressible fluid, is the corresponding gradient;

[0009] Step 2: Define sticky items is an intermediate variable Break the Navier-Stokes equations into two equations:

[0010]

[0011] in, is the intermediate variable, For sticky pairs The Jacobian matrix of , : is the double dot product of the tensor;

[0012] Step 3: Multiply the equation in step 2 by the first test function To test whether the solution of the discrete equation satisfies the partial differential equation being solved:

[0013]

[0014] Then perform partial integration:

[0015]

[0016] Due to the discontinuous finite element, the variable is discontinuous on the surface of the element. represents the flux of the variable on the unit surface, and when it is substituted into the above formula, we get

[0017]

[0018] n is the external normal vector of the unit surface, dV is the volume differential element, dS is the surface differential element, is the transpose of the tensor; Represents the dyadic vector of a tensor;

[0019] Step 4: Then perform integration by parts, we get

[0020]

[0021] in The lifting operator on the surface is used as a penalty term.

[0022] Step 5: Multiply the formula in step 2 by the second test function Perform partial integration

[0023]

[0024] is the inviscid flux, calculated according to step 10, is the viscous flux, calculated according to step 11.

[0025] Step 6: To eliminate the intermediate variables in step 5 In step 4, Substitute this into step 5, and we get

[0026]

[0027] Step 7: Since the discontinuous finite element variables are discontinuous, it is necessary to define the flux of the variables on the surface

[0028]

[0029] Step 8: Define the internal penalty factor

[0030]

[0031] {·} and [·] are introduced to represent the average and step of the left and right values ​​of the surface respectively:

[0032]

[0033] [U]=(U + -U - )

[0034] Step 9: Substitute the internal penalty factor obtained in step 8 Substitute into step 6 and obtain the internal penalty discontinuous finite element discretization formula

[0035]

[0036] Step 10: Stickless Flux from Step 5 Calculation using Riemann flux

[0037]

[0038] Maximum eigenvalue

[0039] λ max ={‖u·n‖+c}

[0040] u and c represent the local fluid velocity and sound speed, respectively, and ‖·‖ indicates the absolute value.

[0041] Step 11: Viscous Flux from Step 5 Take

[0042]

[0043] The penalty factor is C IP =6, p is the algebraic order of the discrete polynomial, hf is the scale of the surface, and is obtained using the following formula:

[0044]

[0045] Among them, V + and V - represents the volume of the left and right units, S f is the area of ​​the current front.

[0046] The beneficial effects of the present invention are as follows:

[0047] (1) Eliminate intermediate variables and improve computational efficiency.

[0048] (2) Just calculate the convective flux, viscous flux, and lifting operator according to the variables and their gradient values ​​at the current integration point.

[0049] (3) The derivative of the gradient of the viscosity term with respect to the variable is eliminated, and the lifting operator is directly calculated using the expression of the viscosity term, which can be easily extended to other types of equations, such as heat conduction equations, structural mechanics equations, and electromagnetic field equations. BRIEF DESCRIPTION OF THE DRAWINGS

[0050] Figure 1 Schematic diagram of the spatial solution of Coulter flow in Example 2. DETAILED DESCRIPTION

[0051] In order to better understand the purpose, structure and function of the present invention, a method for processing discontinuous finite element discrete viscosity terms based on the internal penalty method of the present invention is further described in detail below in conjunction with the accompanying drawings.

[0052] Step 1: Based on the compressible Navier-Stokes equations, it can be written as:

[0053]

[0054] in, and denote the convection term and the viscous term respectively, represents the time partial derivative, represents the spatial gradient operator, U is the conserved variable of the compressible fluid, is the corresponding gradient;

[0055] Step 2: Define sticky items As the intermediate variable, the Navier-Stokes equations are decomposed into two equations:

[0056]

[0057] in, is the intermediate variable, For sticky pairs The Jacobian matrix of , : is the double dot product of the tensor;

[0058] Step 3: Multiply the equation in step 2 by the first test function To test whether the solution of the discrete equation satisfies the partial differential equation being solved:

[0059]

[0060] Then perform partial integration:

[0061]

[0062] Due to the discontinuous finite element, the variable is discontinuous on the surface of the element. represents the flux of the variable on the unit surface, and when it is substituted into the above formula, we get

[0063]

[0064] represents the flux of the variable on the unit surface, n is the external normal vector of the unit surface, dV is the volume differential element, dS is the surface differential element, is the transpose of the tensor; Represents the dyadic vector of a tensor;

[0065] Step 4: Then perform integration by parts, we get

[0066]

[0067] in The lifting operator on the surface is used as a penalty term.

[0068] Step 5: Multiply the formula in step 2 by the second test function Perform partial integration

[0069]

[0070] is the inviscid flux, calculated according to step 10, is the viscous flux, calculated according to step 11.

[0071] Step 6: To eliminate the intermediate variables in step 5 In step 4, Substitute this into step 5, and we get

[0072]

[0073] Step 7: Since the discontinuous finite element variables are discontinuous, it is necessary to define the flux of the variables on the surface

[0074]

[0075] Step 8: Define the internal penalty factor

[0076]

[0077] {·} and [·] are introduced to represent the average and step of the left and right values ​​of the surface respectively:

[0078]

[0079] [U]=(U + -U - )

[0080] Step 9: Substitute the internal penalty factor obtained in step 8 into step 6 to obtain the internal penalty discontinuous finite element discretization formula

[0081]

[0082] Step 10: Stickless Flux from Step 5 Calculation using Riemann flux

[0083]

[0084] Maximum eigenvalue

[0085] λ max ={‖u·n‖+c}

[0086] u and c represent the local fluid velocity and sound speed, respectively, and ‖·‖ indicates the absolute value.

[0087] Step 11: Viscous Flux from Step 5 Take

[0088]

[0089] The penalty factor is C IP =6, p is the algebraic order of the discrete polynomial, h f is the scale of the surface, and is obtained using the following formula:

[0090]

[0091] Among them, V + and V - represents the volume of the left and right units, S f is the area of ​​the current front.

[0092] The accuracy of the solution of this method is verified by using Coulter flow. Coulter flow describes the flow between two flat plates, where the upper plate moves at a constant speed and the bottom plate is stationary. The upper wall temperature T1 = 0.85, the lower wall temperature T0 = 0.8, the height H = 2, the length is 2H, the Reynolds number is 1000, and the Mach number Ma ∞ =0.1; the schematic diagram of its spatial solution is as follows Figure 1 shown.

[0093] Its exact solution is:

[0094]

[0095] Where (u, v) is the velocity field, p is the pressure field, and T is the temperature field. According to the calculation formula in step 9, the density L2 error and numerical accuracy are obtained for different polynomial orders for U, where p = 1, 2, 3, indicating that the exact solution is simulated using linear, quadratic, and cubic polynomial distributions, which can theoretically reach the second, third, and fourth order, respectively.

[0096] The calculation is carried out with the uniform flow field as the initial condition, and the boundary conditions are taken as accurate values. The accuracy of the format in solving the Navier-Stokes equations is verified by the exact solution of the example. Table 1 shows that all order formats have achieved the designed accuracy.

[0097] Table 1 Coulter flow accuracy test

[0098]

[0099] It is to be understood that the present invention is described by some embodiments, and it is known to those skilled in the art that various changes or equivalent substitutions may be made to these features and embodiments without departing from the spirit and scope of the present invention. In addition, under the teachings of the present invention, these features and embodiments may be modified to adapt to specific circumstances and materials without departing from the spirit and scope of the present invention. Therefore, the present invention is not limited by the specific embodiments disclosed herein, and all embodiments falling within the scope of the claims of this application are within the scope of protection of the present invention.

Claims

1. A method for processing viscosity terms in discontinuous finite element discretization based on internal penalty method, characterized in that: The following steps are involved: Step 1: Based on the compressible Navier-Stokes equations, it can be written as: in, and denote the convection term and the viscous term respectively, represents the time partial derivative, represents the spatial gradient operator, U is the conserved variable of the compressible fluid, is the corresponding gradient; Step 2: Define sticky items is an intermediate variable Break the Navier-Stokes equations into two equations: in, is the intermediate variable, For sticky pairs The Jacobian matrix of , : is the double dot product of the tensor; Step 3: Multiply the equation in step 2 by the first test function Test whether the solution to a discretized equation satisfies the partial differential equation being solved: Then perform partial integration: Due to the discontinuous finite element, the variable is discontinuous on the surface of the element. represents the flux of the variable on the unit surface, and when it is substituted into the above formula, we get n is the external normal vector of the unit surface, dV is the volume differential element, dS is the surface differential element, is the transpose of the tensor; Represents the dyadic vector of a tensor; Step 4: Then perform integration by parts, we get in The lifting operator on the surface is used as a penalty term; Step 5: Multiply the formula in step 2 by the second test function Perform partial integration in, is the inviscid flux, is the viscous flux; Step 6: To eliminate the intermediate variables in step 5 In step 4, Substitute this into step 5, and we get Step 7: Since the discontinuous finite element variables are discontinuous, it is necessary to define the flux of the variables on the surface Step 8: Define the internal penalty factor The introduced {·} and [·] represent the average and step of the left and right values ​​of the surface respectively: [U]=(U + -IN - ) Among them, the superscript + represents the value on the left side of the face, and the superscript - represents the value on the right side of the face; Step 9: Substitute the internal penalty factor obtained in step 8 Substitute into step 6 and obtain the internal penalty discontinuous finite element discretization formula 2. The method for processing viscosity terms of discontinuous finite element discretization based on the internal penalty method according to claim 1, characterized in that: The inviscid flux Calculation using Riemann flux Maximum eigenvalue λ max ={‖u·n‖+c} u and c represent the local fluid velocity and sound speed, respectively, and ‖·‖ indicates the absolute value; The viscous flux Take Among them C IP is the penalty factor, p is the algebraic order of the discrete polynomial, h f is the scale of the surface, and is obtained using the following formula: Among them, V + and V - represents the volume of the left and right units, S f is the area of ​​the current front.

Images