CN111611738B - Interface problem simulation method based on stable generalized finite element - Google Patents
Interface problem simulation method based on stable generalized finite element Download PDFInfo
- Publication number
- CN111611738B CN111611738B CN202010436741.XA CN202010436741A CN111611738B CN 111611738 B CN111611738 B CN 111611738B CN 202010436741 A CN202010436741 A CN 202010436741A CN 111611738 B CN111611738 B CN 111611738B
- Authority
- CN
- China
- Prior art keywords
- interface
- gfem
- function
- sgfem
- finite element
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 120
- 238000004088 simulation Methods 0.000 title claims abstract description 28
- 239000011159 matrix material Substances 0.000 claims abstract description 69
- 230000008569 process Effects 0.000 claims description 12
- 238000012545 processing Methods 0.000 claims description 11
- 238000013459 approach Methods 0.000 claims description 7
- 238000004422 calculation algorithm Methods 0.000 claims description 6
- 230000008878 coupling Effects 0.000 claims description 5
- 238000010168 coupling process Methods 0.000 claims description 5
- 238000005859 coupling reaction Methods 0.000 claims description 5
- 239000006185 dispersion Substances 0.000 claims description 4
- 230000006872 improvement Effects 0.000 claims description 4
- 238000012886 linear function Methods 0.000 claims description 4
- 239000012141 concentrate Substances 0.000 claims description 3
- 238000012885 constant function Methods 0.000 claims description 3
- 238000006467 substitution reaction Methods 0.000 claims description 3
- 239000000758 substrate Substances 0.000 claims description 3
- 230000010354 integration Effects 0.000 description 7
- 241000894007 species Species 0.000 description 5
- 238000012937 correction Methods 0.000 description 4
- 238000012360 testing method Methods 0.000 description 4
- 239000000463 material Substances 0.000 description 3
- 230000008901 benefit Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 239000002131 composite material Substances 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000008030 elimination Effects 0.000 description 2
- 238000003379 elimination reaction Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 239000000835 fiber Substances 0.000 description 2
- 239000012530 fluid Substances 0.000 description 2
- 238000002156 mixing Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 239000007787 solid Substances 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007423 decrease Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000013467 fragmentation Methods 0.000 description 1
- 238000006062 fragmentation reaction Methods 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 230000001788 irregular Effects 0.000 description 1
- 239000007788 liquid Substances 0.000 description 1
- 238000005192 partition Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000000518 rheometry Methods 0.000 description 1
- 230000006641 stabilisation Effects 0.000 description 1
- 238000011105 stabilization Methods 0.000 description 1
- 230000036962 time dependent Effects 0.000 description 1
- 230000005514 two-phase flow Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
The invention provides an interface problem simulation method based on stable generalized finite elements, which comprises the following steps: constructing GFEM with a collection enrichment strategy according to interface problem GFEM/XFEM; modifying an enrichment function in the GFEM enrichment strategy to construct a novel SGFEM basic function; introducing a unilateral distance function, and adjusting an SGFEM basic function; and solving the interface problem based on the adjusted SGFEM basic function to complete the simulation of the interface problem. According to the interface problem simulation method based on the stable generalized finite element, an SGFEM is provided for the interface problem, so that the precision of the traditional GFEM is improved, the condition number of the stiffness matrix is greatly reduced, and the numerical simulation of the interface problem is more accurate and robust. The method is suitable for grids and interfaces at any relative position.
Description
Technical Field
The invention relates to the technical field of engineering computational mechanics software application, in particular to an interface problem simulation method based on stable generalized finite elements.
Background
Interfacial problems occur in many engineering applications and physical phenomena such as fluid-solid coupling, bi-material, multiphase flow, porous media, three-dimensional fiber composites, and free boundary problems, among others [1]Brenner H.Interfacial transport processes and rheology.Elsevier,2013 ] [2]Li Z,Ito K.The immersed interface method: numerical solutions of PDEs involving interfaces and irregular domains, siam,2006. Solutions to these problems are typically discontinuous at the interface. Numerical simulation techniques are indispensable for complex interface problems. Taking the two-dimensional ellipse and parabolic interface problem as an example for the description herein, the three-dimensional problem can be solved with similar ideas.
Order theIs a bounded area, +.>Is its boundary. The interface Γ divides Ω into two regions Ω 0 And omega 1 When dealing with smooth interface problems, i.e. Γ is smooth, corner points are not contained, as shown in fig. 1. Consider the elliptical interface problem:
the solution of this problem satisfies the interface condition:
[u] Γ =0, on Ω (2)
At->Upper part
Wherein the method comprises the steps ofIs boundary->Is the unit external normal vector of->Is the unit normal vector at the interface, [ u ]] Γ ＝u 1 -u 0 Representing the jump of the physical quantity of the two areas at the interface. a (P) is a slicing constant function,
wherein a is 0 ，a 1 Is interface constant, satisfies 0 < ζ 0 ≤a i ≤ζ 1 < ++i=0, 1 andit is evident that a (P) is discontinuous at the interface Γ. Condition (2) shows that the true solution is continuous at the interface, but +.>Is discontinuous at the interface, so the smoothness of true solutions is low, which presents challenges to the numerical simulation method.
Approximation space of conventional finite element method (finite element method, FEM)Is as follows:
wherein I is h Is index set of grid nodes, phi i As a finite element medium shape function, if the mesh is triangular,is a piecewise linear function, if quadrilateral, +.>Is a piecewise bilinear function. Conventional FEM, in the face of such interface problems, must either mesh match or refine the mesh along the interface to effectively simulate the problem [ 3] due to the lack of polynomial function in describing the discontinuity characteristics]Lin T，Lin Y，Zhang X.Partially penalized immersed finite element methods for elliptic interface problems.SIAM J.Numer.Anal.，2015，53：1121-1144.Particularly in the complex engineering problem of the interface evolving along with time, the grid division burden is heavy, and the requirement of engineering calculation is difficult to adapt.
In recent decades, a finite element method based on a non-fitting grid has achieved significant progress, and the grid division of the method is simple and the grid is fixed, and the position of the grid is irrelevant to an interface. Generalized or extended finite element method [ 4]]Babuska I，Banerjee U，Osborn J，Survey of meshless and generalized finite element methods：a unified approach，Acta Numer，2003，12∶1-125.[5]FriesT-P，Belytschko T.The extended/generalized finite element method：an overview of the method and its applications.Internat.J.Numer.Methods Engrg.，2010，84：253-304.[6]Babuska I, melenk j.the partition of unity finite element method, international j.numer methods engrg.1997, 40:727-758 (generalized/extended FEM, GFEM/XFEM) belongs to one of the non-fitting finite element methods, which couples special functions capable of expressing true and resolved local information on the basis of a fixed regular grid, thereby improving the approximation accuracy of complex non-smooth problems or meeting the special approximation requirements for specific problems, which are called enrichment functions. GFEM/XFEM has attracted considerable interest over the past 20 years and has been successfully applied to a number of engineering problems such as material modeling, crack propagation and solid-liquid interactions, etc. [4, 5]]. For simplicity, GFEM is used to denote GFEM/XFEM in the following description. Based on the deep research of GFEM on the construction theory of unit-shaped functions, complex problems with any internal features (such as cracks, holes and the like) and external features (such as concave angles, edges and the like) can be solved on simple and area-independent finite element grids [4, 5]]. GFEM has long been applied to solve fixed and mobile interface problems [7]Kergrene K，Babuska I，Banerjee U.Stable Generalized Finite Element Method and associated iterative schemes：application to interface problems.Comput.Methods Appl.Mech.Engrg.，2016，305：1-36.[8]Chessa J，Belytschko T.An extended finite element method for two-phase fluids.J.Appl.Mech.，2003，70：10-17.[9] N，Cloirec M，Cartraud P，Remacle J.A computational approach to handle complex microstructure geometries.Comput.Methods Appl.Mech.Engrg，2003，192/>：3163-3177.[10]Sauerland H，Fries T-P.The stable XFEM for two-phase flows.Comput.&Fluids，2013，87：41-49.[11]Cheng K，Fries T-P，Higher-order XFEM for curved strong and weak discontinuities，Internat.J.Numer.Methods Engrg.，2010，82：564-590.[12]Zhang Q，Babuska I.A stable generalized finite element method (SGFEM) of degree two for interface problems.Comput.Methods Appl.Mech.Engrg.，2020，363：112889。
Two conventional GFEM's are currently in common use: topological generalized finite element method (topological GFEM) and modified generalized finite element method (corrected GFEM). The approximation space of the topological GFEM is
Wherein the enrichment function D is a distance function with respect to the interface Γ, defined as D (P) =dist (P, Γ),is an index set of enriched nodes, defined as +.>e s Representing a grid cell. [11]In medium analysis topology GFEM, they are in [11 ] because of the poor approximation structure caused by the mixing unit (blending elements)]The correction GFEM is proposed to approximate the space.
Wherein is defined asζ is a ramp function, defined as +.>
Currently, these GFEM have several major difficulties in solving interface problems: first, GFEM may lack stability with a stiffness matrix having a condition number much greater than that of FEM [13]Babuska I and U.Banerjee U.Stable generalized finite element method,Comput.Methods Appl.Mech.Engrg.,2011 (201-204): 91-111; second, GFEM may lack robustness, condition numbers may become ill-conditioned as interfaces approach cell boundaries [14]Zhang Q,Babuska I,Banerjee U.Robustness in Stable Generalized Finite Element Methods (SGFEM) applied to Poisson problems with crack singularities, comp.methods appl.mech.engrg.,2016, 311:476-502. The lack of stability and robustness makes the condition number of the stiffness matrix of GFEM large, resulting in catastrophic rounding errors in the elimination, or the convergence speed can become very slow when solving the linear system of equations iteratively. Furthermore, because of the introduction of enrichment functions other than finite element shape functions, many of the standard operations of finite elements cannot be applied to GFEM, such as mass concentration [15]Thomee V.Galerkin finite element methods for parabolic problems.Springer,1984.[16]Menouillard T,R othor J, combescure a, et al efficiency explicit time stepping for the eXtended Finite Element Method (X-FEM). International. J. Methods Engrg,2006, 68 (9): 911-939 [17]Elguedj T,Gravouil A,Maigre H.An explicit dynamics extended finite element method.Part 1 ]: mass lumping for arbitrary enrichment functions. Comput. Methods appl. Mech. Engrg,2009, 198 (30-32): 2297-2317.
Disclosure of Invention
The invention provides an interface problem simulation method based on stable generalized finite elements, which aims to overcome the technical defect that the stability and robustness loss exist when the existing GFEM solves the interface problem, so that disastrous rounding errors in a primordial elimination method are caused, or the convergence speed becomes very slow when a linear equation set is solved iteratively.
In order to solve the technical problems, the technical scheme of the invention is as follows:
the interface problem simulation method based on the stable generalized finite element comprises the following steps:
s1: constructing GFEM with a collection enrichment strategy according to interface problem GFEM/XFEM;
s2: modifying an enrichment function in the GFEM enrichment strategy to construct a novel SGFEM basic function;
s3: introducing a unilateral distance function, and adjusting an SGFEM basic function;
s4: and solving the interface problem based on the adjusted SGFEM basic function to complete the simulation of the interface problem.
In the scheme, the SGFEM is provided for the interface problem, so that the precision of the traditional GFEM is improved, the condition number of the rigidity matrix is greatly reduced, and the numerical simulation of the interface problem is more accurate and robust. The method is suitable for grids and interfaces at any relative position.
In the scheme, the novel SGFEM developed by the invention can improve the precision and condition number of the traditional GFEM/XFEM, so that the numerical simulation of the interface problem is more accurate and robust, and the method has high application value in solving the problems in engineering application and physical phenomena because the interface problem widely exists in the engineering application and physical phenomena.
The step S1 specifically includes:
order theIs a bounded area, +.>Is the boundary thereof, the interface Γ divides Ω into two regions Ω 0 And omega 1 When the problem of smooth interface is solved, i.e. Γ is smooth, no corner points are contained; consider the elliptical interface problem:
at->Upper part
The solution of this problem satisfies the interface condition:
[u] Γ =0, on Ω (2)
Wherein the method comprises the steps ofIs boundary->Is the unit external normal vector of->Is the unit normal vector at the interface, [ u ]] Γ ＝u 1 -u 0 Representing the jump of the two area physical quantities at the interface; a (P) is a slicing constant function,
wherein a is 0 ，a 1 Is interface constant, satisfies 0 < ζ 0 ≤a i ≤ζ 1 < ++i=0, 1 andit is evident that a (P) is discontinuous at the interface Γ; condition (2) shows that the true solution is continuous at the interface, but +.>The interface is discontinuous, so that the smoothness of true solution is low, and challenges are brought to a numerical simulation method; approximation space of general finite element method FEM +.>The expression is as follows:
wherein I is h Is index set of grid nodes, phi i As a finite element medium shape function, if the mesh is triangular,is a piecewise linear function, if quadrilateral, +.>The function of (2) is a piecewise bilinear function; the generalized or extended finite element method GFEM/XFEM belongs to one of non-fitting finite element methods, and on the basis of a fixed regular grid, the coupling can express special functions of true-solution local information, so that the approximation precision of complex non-smooth problems is improved or special approximation requirements on specific problems are met, and the special functions are called enrichment functions; thus, a topological generalized finite element method GFEM with a collective enrichment strategy was constructed.
Wherein the approximation space of the topological GFEM is
Wherein the enrichment function D is a distance function with respect to the interface Γ, defined as D (P) =dist (P, Γ),is an index set of enriched nodes, defined as +.>e s Representing a grid cell; correcting the topology GFEM, wherein the corrected GFEM approximation space is as follows:
wherein is defined asζ is a ramp function, defined as +.>
Wherein, the step S2 is specifically to replace the enrichment function D of GFEM by adding a simple local processing step of subtracting the interpolation term, thereby constructing a novel SGFEM basic function.
In the step S2, a simple local processing step of subtracting the interpolation term is added, which specifically includes:
for smooth interface problems, useInstead of D, wherein->Representing a finite element interpolation function of the continuous function f based on the region omega; GFEM based on this improvement is known as SGFEM [13, 14]The method comprises the steps of carrying out a first treatment on the surface of the The approximation space of SGFEM is defined as follows:
wherein the method comprises the steps ofAs in topology GFEM (4); subtracting their finite element interpolation for enrichment functions, which mitigates their correlationOne of the measures of linear dependence of finite element functions [12-14]。
The step S3 specifically includes:
unlike the distance function D (P) in conventional GFEM: =dist (P, Γ) (Γ represents interface), introducing a new unilateral distance functionThe concrete steps are as follows:
wherein one side of Γ is D (P), and the other side is 0; the one-sided distance function is simpler than the expression of the distance function D, so both numerical integration and algorithms are simplified. Moreover, this one-sided distance function facilitates the demonstration of the optimal convergence of SGFEM [12]. The approximation space of SGFEM using the single-sided distance function is defined as follows:
the adjusted SGFEM basic function can keep optimal convergence, and the scale condition number of the rigidity matrix can be the same order as that of a finite element method.
Wherein, in the step S3, in order to facilitate the display of the result, the rigidity matrix scale condition number SCN of the GFEM is defined:
let A be FEM, GFEM, or SGFEM stiffness matrix, D be a diagonal matrix, defined asThe dimensional condition number of a is defined as follows:
K：＝κ 2 (DAD) (9)
wherein kappa is 2 The condition number of the symmetric matrix is the ratio of the maximum eigenvalue to the minimum nonzero eigenvalue of the symmetric matrix; stiffness matrix condition number K of FEM FEM ＝O(h -2 ) The method comprises the steps of carrying out a first treatment on the surface of the Stiffness moment of formula (8)Array scale condition number K SGFEM The method is the same order as the finite element and has robustness, which solves the robustness difficulty of the stability of the conventional GFEM.
Wherein the method further comprises the steps of:
s5: and processing the adjusted SGFEM basic function based on a quality concentration strategy, so as to achieve the purposes of reducing the quality matrix number and simplifying the algorithm.
In the step S5, the quality concentration policy specifically includes:
a typical parabolic interface problem is adopted to illustrate the quality concentration strategy:
the interface conditions are satisfied:
processing the parabolic interface problem (9) with reference to a discretization process of the interface problem (1), wherein the time dispersion adopts a backward differential format: instant node o=t 0 ＜t 1 ＜...＜t N Time interval [0, T =t instead of time interval]Is called asFor the time step, let { t } n ：t n ＝nτ，0≤n≤N}，/>The parabolic interface problem (9) is therefore based on a finite dimensional space/>The fully discrete format of (a) is as follows:
backward Euler-Galerkin format:
find outMake->
Crank-Nicolson Galerkin format:
find outSo that
Wherein the method comprises the steps ofRewriting (17) and (18) into linear equations, respectively
And
wherein the method comprises the steps ofFor stiffness matrix->For the quality matrix->Is->Is a substrate of (a);
the order of the shape function is rearranged such that the quality matrix M of SGFEM has the form:
wherein M is 11 In connection with the standard finite element FE section, M 22 In relation to the enriched fraction of SGFEM, M 12 Related to the inner product of standard FE and enriched fraction; m is M 11 Matrix size n 2 ，M 22 Matrix size of Cn, n 2 Cn, are constants; the quality centralization strategy is to use a diagonal matrix for MTo approximate substitution, wherein:
in other words the first and second phase of the process,centralizing each row of elements in M to be arranged on the diagonal position of the row; the method adopts a quality matrix which only concentrates the relevant part of the finite element in M, namely only for M in M 11 Concentrating, wherein other parts are kept unchanged, and the concentrated mass matrix is as follows:
wherein,is a diagonal matrix with a diagonal element M 11 Row sums of the corresponding rows.
Wherein, in practical application, because M 11 Matrix size n 2 ，M 22 The matrix size of (C) is Cn, C is a constant, and n is n when the discretization size is large 2 > Cn, i.e., the portion of the enrichment distance function, i.e., the number of enrichment nodes near the interface is one-dimensional less than the number of finite element nodes.
Compared with the prior art, the technical scheme of the invention has the beneficial effects that:
according to the interface problem simulation method based on the stable generalized finite element, an SGFEM is provided for the interface problem, so that the precision of the traditional GFEM is improved, the condition number of the stiffness matrix is greatly reduced, and the numerical simulation of the interface problem is more accurate and robust. The method is suitable for grids and interfaces at any relative position.
Drawings
FIG. 1 is a schematic view of an interface problem area;
FIG. 2 is a schematic flow chart of the method of the present invention;
FIG. 3 is a schematic diagram of enrichment node formats of different GFEM species;
FIG. 4 is a schematic diagram of an example of a cell intersecting an interface and Gaussian integration rules on each triangle;
fig. 5 shows the energy relative error EE for different interface coefficient ratios, left: (ρ=10), right: (ρ=100);
fig. 6 is a scale condition number SCN, left, for different interface coefficient ratios: (ρ=10), right: (ρ=100);
fig. 7 is when δ=1.88×10 -2 Schematic of enrichment schemes of different GFEM species;
FIG. 8 is a stiffness matrix scale condition number;
fig. 9 is a quality concentration strategy: left: 12849 non-zero points before concentration; right: after concentration, 4529 non-zero points.
Detailed Description
The drawings are for illustrative purposes only and are not to be construed as limiting the present patent;
for the purpose of better illustrating the embodiments, certain elements of the drawings may be omitted, enlarged or reduced and do not represent the actual product dimensions;
it will be appreciated by those skilled in the art that certain well-known structures in the drawings and descriptions thereof may be omitted.
The technical scheme of the invention is further described below with reference to the accompanying drawings and examples.
Example 1
As shown in fig. 2, the interface problem simulation method based on the stable generalized finite element comprises the following steps:
s1: constructing GFEM with a collection enrichment strategy according to interface problem GFEM/XFEM;
s2: modifying an enrichment function in the GFEM enrichment strategy to construct a novel SGFEM basic function;
s3: introducing a unilateral distance function, and adjusting an SGFEM basic function;
s4: and solving the interface problem based on the adjusted SGFEM basic function to complete the simulation of the interface problem.
In a specific implementation process, the invention provides the SGFEM aiming at the interface problem, which not only improves the precision of the traditional GFEM, but also greatly reduces the condition number of the rigidity matrix, so that the numerical simulation of the interface problem is more accurate and robust. The method is suitable for grids and interfaces at any relative position.
In a specific implementation process, the novel SGFEM developed by the invention can improve the precision and condition number of the traditional GFEM/XFEM, so that the numerical simulation of the interface problem is more accurate and robust, and the interface problem widely exists in engineering application and physical phenomena, so that the method has high application value in solving the problems in engineering application and physical phenomena.
More specifically, the step S1 specifically includes:
order theIs a bounded area, +.>Is the boundary thereof, the interface Γ divides Ω into two regions Ω 0 And omega 1 When the problem of smooth interface is solved, i.e. Γ is smooth, no corner points are contained; consider the elliptical interface problem:
at->Upper part
The solution of this problem satisfies the interface condition:
[u] Γ =0, on Ω (2)
At->Upper part
Wherein the method comprises the steps ofIs boundary->Is the unit external normal vector of->Is the unit normal vector at the interface, [ u ]] Γ ＝u 1 -u 0 Representing the jump of the two area physical quantities at the interface; a (P) is a fragmentation constantThe function of the function is that,
wherein a is 0 ，a 1 Is interface constant, satisfies 0 < ζ 0 ≤a i ≤ζ 1 < ++i=0, 1 andit is evident that a (P) is discontinuous at the interface Γ; condition (2) shows that the true solution is continuous at the interface, but +.>The interface is discontinuous, so that the smoothness of true solution is low, and challenges are brought to a numerical simulation method; approximation space of general finite element method FEM +.>The expression is as follows:
wherein I is h Is index set of grid nodes, phi i As a finite element medium shape function, if the mesh is triangular,is a piecewise linear function, if quadrilateral, +.>The function of (2) is a piecewise bilinear function; the generalized or extended finite element method GFEM/XFEM belongs to one of non-fitting finite element methods, and on the basis of a fixed regular grid, the coupling can express special functions of true-solution local information, so that the approximation precision of complex non-smooth problems is improved or special approximation requirements on specific problems are met, and the special functions are called enrichment functions; thus, build up with collectionTopology generalized finite element method GFEM incorporating enrichment strategy.
More specifically, the approximation space of the topological GFEM is
Wherein the enrichment function D is a distance function with respect to the interface Γ, defined as D (P) =dist (P, Γ),is an index set of enriched nodes, defined as +.>e s Representing a grid cell; correcting the topology GFEM, wherein the corrected GFEM approximation space is as follows:
wherein is defined asζ is a ramp function defined as
More specifically, the step S2 is specifically to replace the enrichment function D of GFEM by adding a simple local processing step of subtracting the interpolation term, so as to construct a new SGFEM basis function.
More specifically, in the step S2, a simple local processing step of subtracting the interpolation term is added, specifically:
for smooth interface problems, useInstead of D, wherein->Representing a finite element interpolation function of the continuous function f based on the region omega; GFEM based on this improvement is known as SGFEM [13, 14]The method comprises the steps of carrying out a first treatment on the surface of the The approximation space of SGFEM is defined as follows:
wherein the method comprises the steps ofAs in topology GFEM (4); for finite element interpolation where enrichment functions subtract out them, this is one of the approaches to mitigate their linear dependence on finite element functions [12-14 ]]。
More specifically, the step S3 specifically includes:
unlike the distance function D (P) in conventional GFEM: =dist (P, Γ) (Γ represents interface), introducing a new unilateral distance functionThe concrete steps are as follows:
wherein one side of Γ is D (P), and the other side is 0; the one-sided distance function is simpler than the expression of the distance function D, so both numerical integration and algorithms are simplified. Moreover, this one-sided distance function facilitates the demonstration of the optimal convergence of SGFEM [12]. The approximation space of SGFEM using the single-sided distance function is defined as follows:
the adjusted SGFEM basic function can keep optimal convergence, and the scale condition number of the rigidity matrix can be the same order as that of a finite element method.
Wherein, in the step S3, in order to facilitate the display of the result, the rigidity matrix scale condition number SCN of the GFEM is defined:
let A be FEM, GFEM, or SGFEM stiffness matrix, D be a diagonal matrix, defined asThe dimensional condition number of a is defined as follows:
K：＝κ 2 (DAD) (9)
wherein kappa is 2 The condition number of the symmetric matrix is the ratio of the maximum eigenvalue to the minimum nonzero eigenvalue of the symmetric matrix; stiffness matrix condition number K of FEM FEM ＝O(h -2 ) The method comprises the steps of carrying out a first treatment on the surface of the Rigidity matrix scale condition number K of equation (8) SGFEM The method is the same order as the finite element and has robustness, which solves the robustness difficulty of the stability of the conventional GFEM.
More specifically, the method further comprises the steps of:
s5: and processing the adjusted SGFEM basic function based on a quality concentration strategy, so as to achieve the purposes of reducing the quality matrix number and simplifying the algorithm.
In the step S5, the quality concentration policy specifically includes:
the quality concentration strategy is illustrated by a typical parabolic interface problem [15 ]:
the interface conditions are satisfied:
processing the parabolic interface problem (9) with reference to a discretization process of the interface problem (1), wherein the time dispersion adopts a backward differential format: instant node 0=t 0 ＜t 1 ＜...＜t N Time interval [0, T =t instead of time interval]Is called asFor the time step, let { t } n ：t n ＝nτ，0≤n≤N}，/>The parabolic interface problem (9) is therefore based on a limited dimensional space +.>The fully discrete format of (a) is as follows:
backward Euler-Galerkin format:
find outMake->
Crank-Nicolson Galerkin format:
find outSo that
Wherein the method comprises the steps ofWill (17) and (18)Rewritten as a linear system of equations, respectively
And
wherein the method comprises the steps ofFor stiffness matrix->For the quality matrix->Is->Is a substrate of (a);
the order of the shape function is rearranged such that the quality matrix M of SGFEM has the form:
wherein M is 11 In connection with the standard finite element FE section, M 22 In relation to the enriched fraction of SGFEM, M 12 Related to the inner product of standard FE and enriched fraction; m is M 11 Matrix size n 2 ，M 22 Matrix size of Cn, n 2 Cn, are constants; the quality centralization strategy is to use a diagonal matrix for MTo approximate substitution, wherein:
in other words the first and second phase of the process,by focusing each row of elements in M on the diagonal of the row [15]]The method comprises the steps of carrying out a first treatment on the surface of the The quality concentration is a standard operation of time-dependent equations, which can greatly improve the computational efficiency. Since the mass matrix of GFEM contains the fraction of the enrichment function, it is similar to (15). Therefore, the standard practice of the finite element method cannot be applied to the quality of GFEM. There are many literature studies on this problem, such as [16, 17]. SGFEM uses modified enrichment function +.>Unlike normal GFEM, therefore [16, 17]The progress of the study in (c) cannot be applied to SGFEM. The method adopts a quality matrix which only concentrates the relevant part of the finite element in M, namely only for M in M 11 (15) Concentrating, wherein other parts are kept unchanged, and the concentrated mass matrix is as follows:
more specifically, the method comprises the steps of,is a diagonal matrix with a diagonal element M 11 Row sums of the corresponding rows. Studies have shown that mass concentration strategies can be applied very efficiently to SGFEM.
More specifically, in practical applications, M 11 The other parts are not concentrated, but this does not have a great influence. Because of M 11 Matrix size n 2 ，M 22 The matrix size of (C) is Cn, C is a constant, and n is n when the discretization size is large 2 > Cn, i.e. the part of the enrichment distance function, i.e. the number of enrichment nodes near the interface is less one-dimensional than the number of finite element nodes, the enrichment part is extremely small or even negligible compared to the finite element partNot counted.
In a specific implementation process, the invention provides the SGFEM aiming at the interface problem, which not only improves the precision of the traditional GFEM, but also greatly reduces the condition number of the rigidity matrix, so that the numerical simulation of the interface problem is more accurate and robust. The method is suitable for grids and interfaces at any relative position. The innovation herein includes the following three points: two stabilization techniques, including the use of modified enrichment functions (minus their interpolation) and unilateral distance functions, one mass concentration strategy applicable to SGFEM.
In the specific implementation process, the invention develops a novel SGFEM which can improve the precision and condition number of the traditional GFEM/XFEM, so that the numerical simulation of the interface problem is more accurate and robust, and the interface problem widely exists in engineering application and physical phenomena, such as fluid-solid coupling, bi-material, multiphase flow, porous medium, three-dimensional fiber composite material, free boundary problem and the like, so the invention has high application value in solving the problems in engineering application and physical phenomena. In addition, the method has the advantages of strong algorithm innovation, high calculation precision, clear and reasonable software architecture, strong embeddability and expansibility, and has important value for upgrading the current finite element commercial software program.
Example 2
More specifically, on the basis of example 1, solution analysis was performed in combination with specific interface problems:
consider the elliptical interface problem (1) on the area Ω= (0, 1) × (0, 1). Consider a circular interface Γ whose equation is expressed asWherein (1)>The true solution u of the structure (1) is as follows:
wherein (r, θ) is a center of (x) 0 ，y 0 ) Polar coordinates of (c). It can be demonstrated that u satisfies the right-hand term f, g of the interface condition (2), (1) can be obtained by substituting true solutions into equation (1). Two sets of parameters a (P) were tested: (1) a, a 0 =1 and a 1 ＝10；(2)a 0 =1 and a 1 =100. Interface coefficient ratio is defined asThus, the coefficient ratios ρ=10 and ρ=100 for the two sets of coefficients. Region Ω= [0,1]×[0，1]Discretizing with uniform square grid with side length +.>n＝2 j+1 J=1, 2,..7. Obviously, the grid is extremely simple and independent of the location of the interface. We will test and compare the following methods:
FEM: a standard bilinear finite element method (3) based on the uniform square grid;
GFEMTopo: topology GFEM (4) based on enrichment function D;
GFEMcor: with enrichment function ζ (D-D (P) i ) A) corrected GFEM (5);
SGFEM: based on a single-sided distance functionAnd a finite element interpolation modified SGFEM (8).
As shown in fig. 3, the enriched node formats (n=32) of different species of GFEM are shown.
More specifically, the test true solution and the approximation solution u obtained by the method h The energy norm relative error (EE) and the dimensional condition number SCN (9) of the stiffness matrix a are defined as follows:
next, a numerical integration formula is described, using a standard 4 x 4 gaussian integration rule for cells that do not intersect the interface; for the cell intersecting the interface, a straight line is used to connect the interface intersection and the cell boundary, and the cell is decomposed into 4-6 sub-triangles, each of which uses a 12-point Gaussian integration rule, as shown in FIG. 4. This rule is simple and sufficient integration accuracy can be observed in numerical experiments.
In the implementation, fig. 5 and 6 show the relative error of energy norms and condition number of the above-mentioned different methods when h takes different values for the interface coefficient ratio ρ=10 or 100, respectively. As can be seen from fig. 5, when the interface coefficient ratio ρ=10, the convergence order of FEM, topology GFEM is O (h 0.5 ) Correcting the convergence order of the GFEM and the SGFEM to be O (h); when the interface coefficient ratio ρ=100, the convergence order of SGFEM and modified SGFEM is still O (h), the convergence order of FEM and topology GFEM is improved, but the optimum O (h) is not reached, and in fact, when h is sufficiently small, the convergence of FEM and topology GFEM is still O (h 0.5 ). It can be seen from FIG. 6 that both FEM and SGFEM are of the same order O (h -2 ) The increase, which is the SCN order of a typical FEM, is also the condition number O (h in the case of topology and correction GFEM -2 ) This does not seem to be different from SGFEM, but we will see the behavior of each method in the following robustness test.
More specifically, robustness (robustness), which is an important performance index of GFEM, is discussed, and any interface position and any meshing situation may occur in practical application, and no meshing is related to the interface, which is an important advantage of GFEM. We will test the robustness of different GFEM species, i.e. how the condition number changes as the linear interface approaches the cell boundary.
The region Ω= [0,1] × [0,1] was discretized with a uniform grid of 16×16, the interface locations were represented by the square y=0.375+δ, and fig. 7 shows the enriched node schemes for different GFEM species.
Let δ=0.0375×2 -j (j=1, 2,.,. 20), grid parametersAt this time->Ranging from 0.3 to 5.72X10 -7 . The stiffness matrix condition number change for FEM, topology GFEM and SGFEM will be observed as δ decreases, interface Γ approaching grid line y=0.375.
When ρ=10 and ρ=100, the dimensional condition number of the stiffness matrix is as shown in fig. 8, and it is apparent that the condition number of the correction GFEM varies drastically when the interface is close to the cell boundary, meaning that the correction GFEM is not robust. FEM, topology GFEM and SGFEM are robust, their condition numbers remain smooth as delta changes. From numerical experiments, it can be seen that SGFEM is the only method to achieve both the best convergence order and stability and robustness in the GFEM presented herein.
More specifically, consider the parabolic interface problem (10) on the area Ω= (0, 1) × (0, 1). The curved interface expression and the scope of investigation are consistent with the previous case, and the true solution u is constructed as follows:
two sets of parameters a (P) were still tested: (1) Alpha 0 =1 and a 1 ＝10；(2)a 0 =1 and a 1 =100. The ratio of them is ρ=10 and ρ=100, respectively. The finite element dispersion scheme is as described above; make the time interval [0, T ]]＝[0，1]Time stept n N, both time-discrete schemes backward euler (11) and CN (12) are considered. Fig. 9 shows the distribution of non-zero elements of the SGFEM quality matrix when the mesh subdivision is 32 x 32.
In the specific implementation process, the final step of approximate solution is calculatedL of (2) 2 Errors and energy errors, the number of scale conditions of the mass matrix was also tested. Tables 1 to 4 +.>Represents L 2 Norm relative error, ++>Representing the energy norm relative error. As can be seen from a combination of the data in the table below, the errors before and after the mass concentration are substantially identical, with an order of magnitude reduction in SCN after concentration for the condition number SCN of the mass matrix. The effectiveness of the mass concentration approach of the present invention can thus be found.
Table 1ρ=10 relative error and condition number of backward Euler format
Table 2ρ=100 relative error and condition number of backward Euler format
Table 3ρ=10 relative error and condition number of the C rank-Nicolson format
/>
Table 4ρ=100 relative error and condition number of the C rank-Nicolson format
It is to be understood that the above examples of the present invention are provided by way of illustration only and not by way of limitation of the embodiments of the present invention. Other variations or modifications of the above teachings will be apparent to those of ordinary skill in the art. It is not necessary here nor is it exhaustive of all embodiments. Any modification, equivalent replacement, improvement, etc. which come within the spirit and principles of the invention are desired to be protected by the following claims.
Claims (2)
1. The interface problem simulation method based on the stable generalized finite element is characterized by comprising the following steps of:
s1: according to interface problem GFEM/XFEM, constructing GFEM with a collection enrichment strategy, specifically:
order theIs a bounded area, +.>Is the boundary thereof, the interface Γ divides Ω into two regions Ω 0 And omega 1 When the problem of smooth interface is treated, the gamma is smooth and does not contain corner points; consider the elliptical interface problem:
on Ω (1)/(1)>At->Upper part
The solution of the elliptical interface problem satisfies the interface condition:
[u] Γ =0, on Ω (2)
At->Upper part
Wherein the method comprises the steps ofIs boundary->Is the unit external normal vector of->Is the unit normal vector at the interface, [ u ]] Γ ＝u 1 -u 0 Representing the jump of the two area physical quantities at the interface; a (P) is a slicing constant function,
wherein a is 0 ，a 1 Is interface constant, satisfies 0 < ζ 0 ≤a i ≤ζ 1 < ++i=0, 1 andit is evident that a (P) is discontinuous at the interface Γ; interface condition (2) indicates that the true solution is continuous at the interface, but +.>The interface is discontinuous, so that the smoothness of true solution is low, and challenges are brought to a numerical simulation method; approximation space of general finite element method FEM +.>The expression is as follows:
wherein I is h Is index set of grid nodes, phi i As a finite element medium shape function if the mesh is triangularThe shape of the product is that,is a piecewise linear function, if quadrilateral, +.>The function of (2) is a piecewise bilinear function; the generalized or extended finite element method GFEM/XFEM belongs to one of non-fitting finite element methods, and on the basis of a fixed regular grid, the coupling can express special functions of true-solution local information, so that the approximation precision of complex non-smooth problems is improved or special approximation requirements on specific problems are met, and the special functions are called enrichment functions; therefore, constructing a topological generalized finite element method GFEM with a set enrichment strategy;
wherein,
the approximation space of the topological GFEM is
Wherein the enrichment function D is a distance function with respect to the interface Γ, defined as D (P) =dist (P, Γ), and abbreviated as D in equations (4) - (6),is an index set of enriched nodes, define +.> e s Representing a grid cell; correcting the topology GFEM, wherein the corrected GFEM approximation space is as follows:
wherein is defined asζ is a ramp function, defined as +.>
S2: modifying an enrichment function in the GFEM enrichment strategy to construct a novel SGFEM basic function, specifically, replacing the enrichment function D of the GFEM by adding a simple local processing step of subtracting an interpolation term, thereby constructing the novel SGFEM basic function; wherein,
the local processing step of adding a simple subtraction interpolation term is specifically:
for smooth interface problems, useInstead of D, wherein->Representing a finite element interpolation function of the continuous function f based on the region omega; GFEM based on this improvement is called SGFEM; the approximation space of SGFEM is defined as follows:
subtracting their finite element interpolation for the enrichment function, which is one of the approaches to mitigate its linear dependence on the finite element function;
s3: introducing a unilateral distance function to adjust an SGFEM basic function, wherein the unilateral distance function specifically comprises the following steps:
introducing a new unilateral distance functionThe concrete steps are as follows:
wherein one side of Γ is D (P), and the other side is 0; in the formula (8)Abbreviated as +.>The approximation space of SGFEM using the single-sided distance function is defined as follows:
the adjusted SGFEM basic function keeps optimal convergence, and the scale condition number of the rigidity matrix is the same as that of the finite element method;
to facilitate the presentation of the results, the stiffness matrix scale condition number SCN of GFEM is defined:
let A be the stiffness matrix of FEM, GFEM, or SGFEM, D be a diagonal matrix, defined asThe dimensional condition number of a is defined as follows:
K：＝κ 2 (DAD) (9)
wherein kappa is 2 The condition number of the symmetric matrix is the ratio of the maximum eigenvalue to the minimum nonzero eigenvalue of the symmetric matrix; stiffness matrix condition number K of FEM FEM ＝O(h -2 ) The method comprises the steps of carrying out a first treatment on the surface of the Rigidity matrix scale condition number K of (8) SGFEM The method is same as the finite element and has robustness, which solves the difficulty of robustness of the stability of the conventional GFEM
S4: solving an interface problem based on the adjusted SGFEM basic function to complete the simulation of the interface problem;
s5: the adjusted SGFEM basic function is processed based on a quality centralization strategy, so as to achieve the purposes of reducing the quality matrix number and simplifying the algorithm, wherein,
the quality concentration strategy specifically comprises the following steps:
a typical parabolic interface problem is adopted to illustrate the quality concentration strategy:
the interface conditions are satisfied:
processing the parabolic interface problem by referring to the discretization process of the elliptical interface problem (1), wherein the time dispersion adopts a backward differential format: with node 0=t 0 ＜t 1 ＜…＜t N Time interval [0, T =t instead of time interval]Is called asFor the time step, let { t } n ：t n ＝nτ，0≤n≤N}，/>The parabolic interface problem is therefore based on the finite dimensional space +.>The fully discrete format of (a) is as follows:
backward Euler-Galerkin format:
find outMake->
Crank-Nicolson Galerkin format:
find outSo that
Wherein the method comprises the steps ofRewriting (11) and (12) into a linear system of equations, respectively
And
wherein the method comprises the steps ofFor stiffness matrix->In the form of a quality matrix,is->Is a substrate of (a);
the order of the shape function is rearranged such that the quality matrix M of SGFEM has the form:
wherein M is 11 In connection with the standard finite element FE section, M 22 In relation to the enriched fraction of SGFEM, M 12 Related to the inner product of standard FE and enriched fraction; m is M 11 Matrix size n 2 ，M 22 Matrix size of Cn, n 2 Cn is a constant; the quality centralization strategy is to use a diagonal matrix for MTo approximate substitution, wherein:
in other words the first and second phase of the process,centralizing each row of elements in M to be arranged on the diagonal position of the row; the method adopts a quality matrix which only concentrates the relevant part of the finite element in M and only aims at M in M 11 Concentrating, wherein other parts are kept unchanged, and the concentrated mass matrix is as follows:
wherein,is a diagonal matrix with a diagonal element M 11 Row sums of the corresponding rows.
2. The interface problem simulation method based on stable generalized finite element according to claim 1, wherein in practical application, because of M 11 Matrix size n 2 ，M 22 The matrix size of (C) is Cn, C is a constant, and n is n when the discretization size is large 2 The number of enrichment nodes near the interface is one-dimensional less than the number of finite element nodes.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010436741.XA CN111611738B (en) | 2020-05-21 | 2020-05-21 | Interface problem simulation method based on stable generalized finite element |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010436741.XA CN111611738B (en) | 2020-05-21 | 2020-05-21 | Interface problem simulation method based on stable generalized finite element |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111611738A CN111611738A (en) | 2020-09-01 |
| CN111611738B true CN111611738B (en) | 2024-03-29 |
Family
ID=72195706
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010436741.XA Active CN111611738B (en) | 2020-05-21 | 2020-05-21 | Interface problem simulation method based on stable generalized finite element |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111611738B (en) |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106446432A (en) * | 2016-09-30 | 2017-02-22 | 云翼超算（北京）软件科技有限公司 | Optimal transportation meshless method for solving large deformation of material |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US10061878B2 (en) * | 2015-12-22 | 2018-08-28 | Dassault Systemes Simulia Corp. | Effectively solving structural dynamics problems with modal damping in physical coordinates |
-
2020
- 2020-05-21 CN CN202010436741.XA patent/CN111611738B/en active Active
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106446432A (en) * | 2016-09-30 | 2017-02-22 | 云翼超算（北京）软件科技有限公司 | Optimal transportation meshless method for solving large deformation of material |
Non-Patent Citations (5)
| Title |
|---|
| "A rigorous and unified mass lumping scheme for higher-order elements";Yongtao Yang等;《Computer Methods in Applied Mechanics and Engineering》;第319卷;第491-512页 * |
| "A stable generalized finite element method (SGFEM) of degree two for interface problem";Qinghui Zhang等;《Computer Methods in Applied Mechanics and Engineering》;第363卷;第1-19页 * |
| "Stable generalized finite element method (SGFEM) for parabolic interface problems";Pengfei Zhu等;《Journal of Computational and Applied Mathematics》;第367卷;第1-17页 * |
| Qinghui Zhang等."A stable generalized finite element method (SGFEM) of degree two for interface problem".《Computer Methods in Applied Mechanics and Engineering》.2020,第363卷第1-19页. * |
| 基于改进型XFEM的裂纹分析并行软件实现;王理想;文龙飞;王景焘;田荣;;中国科学:技术科学(11);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN111611738A (en) | 2020-09-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Talischi et al. | Polygonal finite elements for topology optimization: A unifying paradigm | |
| Carrillo et al. | A finite-volume method for nonlinear nonlocal equations with a gradient flow structure | |
| CN108563843B (en) | Method for updating disturbance area of steady compressible flow | |
| Gressier et al. | Positivity of flux vector splitting schemes | |
| Liu et al. | Stress optimization of smooth continuum structures based on the distortion strain energy density | |
| Xing et al. | A finite element‐based level set method for structural optimization | |
| CN116629079B (en) | Method and device for constructing mixed finite element space and solving linear elastic mechanical problem | |
| Esser et al. | An extended finite element method applied to levitated droplet problems | |
| Cao et al. | A new primal-dual weak Galerkin method for elliptic interface problems with low regularity assumptions | |
| CN111611738B (en) | Interface problem simulation method based on stable generalized finite element | |
| Arami Fadafan et al. | Moving particle semi-implicit method with improved pressures stability properties | |
| CN115203997A (en) | Dot matrix-entity composite structure topology optimization method based on multivariate design | |
| Caseiro et al. | A systematic development of EAS three-dimensional finite elements for the alleviation of locking phenomena | |
| Sasikala et al. | On an efficient octic order sub-parametric finite element method on curved domains | |
| Sheshadri | An analysis of stability of the flux reconstruction formulation with applications to shock capturing | |
| Qiu et al. | Stress-based evolutionary topology optimization via XIGA with explicit geometric boundaries | |
| Askar et al. | Analytical and Numerical Results for the Diffusion-Reaction Equation When the Reaction Coefficient Depends on Simultaneously the Space and Time Coordinates | |
| Remacle et al. | Transient Mesh Adaptation Using Conforming and Non Conforming Mesh Modifications. | |
| Abdi | Evolutionary topology optimization of continuum structures using X-FEM and isovalues of structural performance | |
| Ghoneim | A meshfree interface-finite element method for modelling isothermal solutal melting and solidification in binary systems | |
| Ren et al. | A meshfree method for Signorini problems using boundary integral equations | |
| Zheng et al. | Dual surface based approach to block decomposition of solid models | |
| Chen et al. | Topology Optimization Based on SA-BESO | |
| Zhu et al. | ROBUST AND EFFICIENT MIXED HYBRID DISCONTINUOUS FINITE ELEMENT METHODS FOR ELLIPTIC INTERFACE PROBLEMS. | |
| Chakib et al. | Multiscale computational method for nonlinear heat transmission problem in periodic porous media |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |