김경종
21 KiB
10.7 Modeling discontinuities as an enriched feature using the extended finite element method
• “Modeling discontinuities as an enriched feature using the extended finite element method,” Section 10.7.1
10.7.1 MODELING DISCONTINUITIES AS AN ENRICHED FEATURE USING THE EXTENDED FINITE ELEMENT METHOD
Products: Abaqus/Standard Abaqus/CAE Abaqus/Viewer
References
• *ENRICHMENT
• *ENRICHMENT ACTIVATION
• “Using the extended finite element method to model fracture mechanics,” Section 31.3 of the Abaqus/CAE User’s Guide
Overview
Modeling discontinuities, such as cracks, as an enriched feature:
• is commonly referred to as the extended finite element method (XFEM);
• is an extension of the conventional finite element method based on the concept of partition of unity;
• allows the presence of discontinuities in an element by enriching degrees of freedom with special displacement functions;
• enables the modeling of discontinuities in the fluid pressure field as well as fluid flow within the cracked element surfaces as in hydraulically driven fracture;
• does not require the mesh to match the geometry of the discontinuities;
• is a very attractive and effective way to simulate initiation and propagation of a discrete crack along an arbitrary, solution-dependent path without the requirement of remeshing in the bulk materials;
• can be simultaneously used with the surface-based cohesive behavior approach (see “Surface-based cohesive behavior,” Section 37.1.10) or the Virtual Crack Closure Technique (see “Crack propagation analysis,” Section 11.4.3), which are best suited for modeling interfacial delamination;
• can be performed using the static procedure (see “Static stress analysis,” Section 6.2.2), the implicit dynamic procedure (see “Implicit dynamic analysis using direct integration,” Section 6.3.2), the low-cycle fatigue analysis using the direct cyclic approach (see “Low-cycle fatigue analysis using the direct cyclic approach,” Section 6.2.7), the geostatic stress field procedure (see “Geostatic stress state,” Section 6.8.2), or coupled pore fluid diffusion/stress analysis (see “Coupled pore fluid diffusion and stress analysis,” Section 6.8.1);
• can also be used to perform contour integral evaluations for an arbitrary stationary surface crack without the need to refine the mesh around the crack tip;
• allows contact interaction of cracked element surfaces based on a small-sliding formulation;
• allows the application of distributed pressure loads to the cracked element surfaces;
• allows the output of some surface variables on the cracked element surfaces;
• allows both material and geometrical nonlinearity; and
• is available only for first-order stress/displacement solid continuum elements, first-order displacement/pore pressure solid continuum elements, and second-order stress/displacement tetrahedron elements.
Modeling approach
Modeling stationary discontinuities, such as a crack, with the conventional finite element method requires that the mesh conforms to the geometric discontinuities. Therefore, considerable mesh refinement is needed in the neighborhood of the crack tip to capture the singular asymptotic fields adequately. Modeling a growing crack is even more cumbersome because the mesh must be updated continuously to match the geometry of the discontinuity as the crack progresses.
The extended finite element method (XFEM) alleviates the shortcomings associated with meshing crack surfaces. The extended finite element method was first introduced by Belytschko and Black (1999). It is an extension of the conventional finite element method based on the concept of partition of unity by Melenk and Babuska (1996), which allows local enrichment functions to be easily incorporated into a finite element approximation. The presence of discontinuities is ensured by the special enriched functions in conjunction with additional degrees of freedom. However, the finite element framework and its properties such as sparsity and symmetry are retained.
Introducing nodal enrichment functions
For the purpose of fracture analysis, the enrichment functions typically consist of the near-tip asymptotic functions that capture the singularity around the crack tip and a discontinuous function that represents the jump in displacement across the crack surfaces. The approximation for a displacement vector function with the partition of unity enrichment is
\mathbf {u} = \sum_ {I = 1} ^ {N} N _ {I} (x) [ \mathbf {u} _ {I} + H (x) \mathbf {a} _ {I} + \sum_ {\alpha = 1} ^ {4} F _ {\alpha} (x) \mathbf {b} _ {I} ^ {\alpha} ],
where N _ { I } ( x ) are the usual nodal shape functions; the first term on the right-hand side of the above equation, { \mathbf { u } } _ { I } , , is the usual nodal displacement vector associated with the continuous part of the finite element solution; the second term is the product of the nodal enriched degree of freedom vector, \mathbf { a } _ { I } , and the associated discontinuous jump function H ( x ) across the crack surfaces; and the third term is the product of the nodal enriched degree of freedom vector, { \bf b } _ { I } ^ { \alpha } , and the associated elastic asymptotic crack-tip functions, F _ { \alpha } ( x ) . The first term on the right-hand side is applicable to all the nodes in the model; the second term is valid for nodes whose shape function support is cut by the crack interior; and the third term is used only for nodes whose shape function support is cut by the crack tip.
Figure 10.7.1–1 illustrates the discontinuous jump function across the crack surfaces, H ( x ) , which is given by
H (\mathrm{x}) = \left\{ \begin{array}{c l} 1 & \text {if} (\mathbf {x} - \mathbf {x} ^ {*}). \mathbf {n} \geq 0, \\ - 1 & \text {otherwise}, \end{array} \right.
text_image
Crack tip n s X* X r θ s n X
Figure 10.7.1–1 Illustration of normal and tangential coordinates for a smooth crack.
where is a sample (Gauss) point, \mathbf { x } ^ { * } is the point on the crack closest to , and is the unit outward normal to the crack at \mathbf { x } ^ { * } .
Figure 10.7.1–1 illustrates the asymptotic crack tip functions in an isotropic elastic material, F _ { \alpha } ( x ) , which are given by
F _ {\alpha} (x) = [ \sqrt {r} \sin \frac {\theta}{2}, \sqrt {r} \cos \frac {\theta}{2}, \sqrt {r} \sin \theta \sin \frac {\theta}{2}, \sqrt {r} \sin \theta \cos \frac {\theta}{2} ],
where ( r , \theta ) is a polar coordinate system with its origin at the crack tip and \theta = 0 is tangent to the crack at the tip.
These functions span the asymptotic crack-tip function of elasto-statics, and \begin{array} { r } { { \sqrt { r } } \sin { \frac { \theta } { 2 } } } \end{array} takes into account the discontinuity across the crack face. The use of asymptotic crack-tip functions is not restricted to crack modeling in an isotropic elastic material. The same approach can be used to represent a crack along a bimaterial interface, impinged on the bimaterial interface, or in an elastic-plastic power law hardening material. However, in each of these three cases different forms of asymptotic crack-tip functions are required depending on the crack location and the extent of the inelastic material deformation. The different forms for the asymptotic crack-tip functions are discussed by Sukumar (2004), Sukumar and Prevost (2003), and Elguedj (2006), respectively.
Accurately modeling the crack-tip singularity requires constantly keeping track of where the crack propagates and is cumbersome because the degree of crack singularity depends on the location of the crack in a non-isotropic material. Therefore, we consider the asymptotic singularity functions only when modeling stationary cracks in Abaqus/Standard. Moving cracks are modeled using one of the two alternative approaches described below.
Modeling moving cracks with the cohesive segments method and phantom nodes
One alternative approach within the framework of XFEM is based on traction-separation cohesive behavior. This approach is used in Abaqus/Standard to simulate crack initiation and propagation. This is a very general interaction modeling capability, which can be used for modeling brittle or ductile fracture. The other crack initiation and propagation capabilities available in Abaqus/Standard are based on cohesive elements (“Defining the constitutive response of cohesive elements using a traction-separation description,” Section 32.5.6) or on surface-based cohesive behavior (“Surface-based
cohesive behavior,” Section 37.1.10). Unlike these methods, which require that the cohesive surfaces align with element boundaries and the cracks propagate along a set of predefined paths, the XFEM-based cohesive segments method can be used to simulate crack initiation and propagation along an arbitrary, solution-dependent path in the bulk materials, since the crack propagation is not tied to the element boundaries in a mesh. In this case the near-tip asymptotic singularity is not needed, and only the displacement jump across a cracked element is considered. Therefore, the crack has to propagate across an entire element at a time to avoid the need to model the stress singularity.
Phantom nodes, which are superposed on the original real nodes, are introduced to represent the discontinuity of the cracked elements, as illustrated in Figure 10.7.1–2. When the element is intact, each phantom node is completely constrained to its corresponding real node. When the element is cut through by a crack, the cracked element splits into two parts. Each part is formed by a combination of some real and phantom nodes depending on the orientation of the crack. Each phantom node and its corresponding real node are no longer tied together and can move apart.
text_image
original nodes phantom nodes crack Ω₀⁺ crack Ω₀⁻ Ωₚ⁺ crack Ω₀⁻
Figure 10.7.1–2 The principle of the phantom node method.
The magnitude of the separation is governed by the cohesive law until the cohesive strength of the cracked element is zero, after which the phantom and the real nodes move independently. To have a set of full interpolation bases, the part of the cracked element that belongs in the real domain, \Omega _ { 0 } , is extended to the phantom domain, \Omega _ { p } . Then the displacement in the real domain, \Omega _ { 0 } , can be interpolated by using the degrees of freedom for the nodes in the phantom domain, \Omega _ { p } . The jump in the displacement field is realized by simply integrating only over the area from the side of the real nodes up to the crack; i.e., \Omega _ { 0 } ^ { + } and \Omega _ { 0 } ^ { - } . This method provides an effective and attractive engineering approach and has been used for simulation of the initiation and growth of multiple cracks in solids by Song (2006) and Remmers (2008). It has been proven to exhibit almost no mesh dependence if the mesh is sufficiently refined.
Modeling hydraulically driven fracture
The cohesive segments method in conjunction with phantom nodes discussed above can also be extended to model hydraulically driven fracture. In this case additional phantom nodes with pore pressure degrees
of freedom are introduced on the edges of each enriched element to model the fluid flow within the cracked element surfaces in conjunction with the phantom nodes that are superposed on the original real nodes to represent the discontinuities of displacement and fluid pressure in a cracked element. The phantom node at each element edge is not activated until the edge is intersected by a crack. The flow patterns of the pore fluid in the cracked elements are shown in Figure 10.7.1–3. The fluid is assumed to be incompressible. The fluid flow continuity, which accounts for both tangential and normal flow within and across the cracked element surfaces as well as the rate of opening of the cracked element surfaces, is maintained. The fluid pressure on the cracked element surfaces contributes to the traction-separation behavior of the cohesive segments in the enriched elements, which enables the modeling of hydraulically driven fracture.
text_image
cracked enriched elements tangential flow normal flow
Figure 10.7.1–3 Flow within a cracked element.
Modeling moving cracks based on the principles of linear elastic fracture mechanics (LEFM) and phantom nodes
Another alternative approach to modeling moving cracks within the framework of XFEM is based on the principles of linear elastic fracture mechanics (LEFM). Therefore, it is more appropriate for problems in which brittle crack propagation occurs. Similar to the XFEM-based cohesive segments method described above, the near-tip asymptotic singularity is not considered, and only the displacement jump across a cracked element is considered. Therefore, the crack has to propagate across an entire element at a time to avoid the need to model the stress singularity. The strain energy release rate at the crack tip is calculated based on the modified Virtual Crack Closure Technique (VCCT), which has been used to model delamination along a known and partially bonded surface (see “Crack propagation analysis,” Section 11.4.3). However, unlike this method, the XFEM-based LEFM approach can be used to simulate crack propagation along an arbitrary, solution-dependent path in the bulk material without the requirement of a pre-existing crack in the model.
The modeling technique is very similar to the XFEM-based cohesive segment approach described above where phantom nodes are introduced to represent the discontinuity of the cracked element when
the fracture criterion is satisfied. The real node and the corresponding phantom node will separate when the equivalent strain energy release rate exceeds the critical strain energy release rate at the crack tip in an enriched element. The traction is initially carried as equal and opposite forces on the two surfaces of the cracked element. The traction is ramped down linearly over the separation between the two surfaces with the dissipated strain energy equal to either the critical strain energy required to initiate the separation or the critical strain energy required to propagate the crack depending on whether the VCCT or the enhanced VCCT criterion is specified.
Modeling low-cycle fatigue crack propagation based on the principles of LEFM
The XFEM-based LEFM approach can also be used to simulate a discrete crack growth subjected to subcritical cyclic loading in a low-cycle fatigue analysis using the direct cyclic approach (“Low-cycle fatigue analysis using the direct cyclic approach,” Section 6.2.7). The fracture energy release rates at the crack tips in the enriched elements are calculated based on the above mentioned modified VCCT technique. The onset and crack growth are characterized by using the Paris law, which relates the relative fracture energy release rates to crack growth rates as illustrated in Figure 10.7.1–4. This approach has been used to model progressive delamination under a sub-critical cyclic loading along a known and partially bonded surface (see “Low-cycle fatigue criterion” in “Crack propagation analysis,” Section 11.4.3). However, unlike this method, the XFEM-based LEFM approach can be used to simulate fatigue crack propagation along an arbitrary, solution-dependent path in the bulk material.
line
| G | da/dN |
|---|---|
| G_thresh | Low |
| G_pl | High |
| G_C | High |
Figure 10.7.1–4 Fatigue crack growth governed by the Paris law.
Using the level set method to describe discontinuous geometry
A key development that facilitates treatment of cracks in an extended finite element analysis is the description of crack geometry, because the mesh is not required to conform to the crack geometry. The level set method, which is a powerful numerical technique for analyzing and computing interface motion, fits naturally with the extended finite element method and makes it possible to model arbitrary crack growth without remeshing. The crack geometry is defined by two almost-orthogonal signed distance functions, as illustrated in Figure 10.7.1–5. The first, \phi , describes the crack surface, while the second, \psi , is used to construct an orthogonal surface so that the intersection of the two surfaces gives the crack front. \mathbf { n } ^ { + } indicates the positive normal to the crack surface; \mathbf { m } ^ { + } indicates the positive normal to the crack front. No explicit representation of the boundaries or interfaces is needed because they are entirely described by the nodal data. Two signed distance functions per node are generally required to describe a crack geometry.
text_image
crack surface (φ = 0) orthogonal surface (ψ = 0) m+ n+ crack front (intersection of ψ and φ)
Figure 10.7.1–5 Representation of a nonplanar crack in three dimensions by two signed distance functions \phi and \psi . .
Defining an enriched feature and its properties
You must specify an enriched feature and its properties. One or multiple pre-existing cracks can be associated with an enriched feature. In addition, during an analysis one or multiple cracks can initiate in an enriched feature without any initial defects. However, multiple cracks can nucleate in a single enriched feature only when the damage initiation criterion is satisfied in multiple elements in the same
time increment. Otherwise, additional cracks will not nucleate until all the pre-existing cracks in an enriched feature have propagated through the boundary of the given enriched feature. If several crack nucleations are expected to occur at different locations sequentially during an analysis, multiple enriched features can be specified in the model. Enriched degrees of freedom are activated only when an element is intersected by a crack. Only stress/displacement or displacement/pore pressure solid continuum elements can be associated with an enriched feature.
Input File Usage: *ENRICHMENT
Abaqus/CAE Usage: Interaction module: Special→Crack→Create→XFEM
Defining the type of enrichment
You can choose to model an arbitrary stationary crack or a discrete crack propagation along an arbitrary, solution-dependent path. The former requires that the elements around the crack tips are enriched with asymptotic functions to catch the singularity and that the elements intersected by the crack interior are enriched with the jump function across the crack surfaces. The latter infers that crack propagation is modeled with either the cohesive segments method or the linear elastic fracture mechanics approach in conjunction with phantom nodes. However, the options are mutually exclusive and cannot be specified simultaneously in a model.
| Input File Usage: | Use the following option to specify a crack propagation analysis (default): *ENRICHMENT, TYPE=PROPAGATION CRACKUse the following option to specify an analysis with stationary cracks: *ENRICHMENT, TYPE=STATIONARY CRACK |
| Abaqus/CAE Usage: | Use the following input to specify a crack propagation analysis:Interaction module: crack editor: toggle on Allow crack growthUse the following input to specify an analysis with stationary cracks:Interaction module: crack editor: toggle off Allow crack growth |
Assigning a name to the enriched feature
You must assign a name to an enriched feature, such as a crack. This name can be used in defining the initial location of the crack surfaces, in identifying a crack for contour integral output, in activating or deactivating the crack propagation analysis, and in generating cracked element surfaces.
Input File Usage: *ENRICHMENT, NAME=name
Abaqus/CAE Usage: Interaction module: Special→Crack→Create: XFEM: Name: name
Identifying an enriched region
You must associate the enrichment definition with a region of your model. Only degrees of freedom in elements within these regions are potentially enriched with special functions. The region should consist of elements that are presently intersected by cracks and those that are likely to be intersected by cracks as the cracks propagate.