fea-fracture-mechanics

name	fea-fracture-mechanics
description	FEA fracture mechanics — singular crack-tip elements, J-integral domain method, SIF from FEA (K_I/K_II/K_III), VCCT, cohesive zone model, contour integral ABAQUS/ANSYS.
metadata	{"priority":7,"promptSignals":{"phrases":["FEA fracture mechanics","crack tip element","J-integral FEA","stress intensity factor FEA","VCCT","cohesive zone model","contour integral ABAQUS"],"minScore":3}}

FEA Fracture Mechanics — Complete Skill

Crack-Tip Stress Singularity

Near crack tip, stresses follow r^(-1/2) singularity: σ_ij = K_I / √(2πr) × f_ij(θ) [Mode I]

Standard isoparametric elements cannot capture this — they assume polynomial displacement fields

Singular Crack-Tip Elements

Quarter-Point Elements (QPE)

Concept: move mid-side nodes of 8-node serendipity quad (or 20-node hex) to quarter-point positions (0.25L from crack tip) This produces r^(1/2) displacement variation → r^(-1/2) strain/stress singularity

2D degenerate triangle: Collapse one side of 8-node quad; mid-side nodes at quarter-point Results in r^(-1/2) singularity in all directions from collapsed edge

Construction:

• Ring of QPE (wedge/spider) around crack tip
• Next ring: standard elements
• Fan mesh: 8–16 elements around crack tip typical

Rules:

• QPE size r_tip ≈ a/8 to a/15 (a = crack length)
• Mesh at crack faces: compatible, no penetration (contact or duplicate nodes)

K Extraction from QPE Mesh

Displacement correlation method: K_I = E'/(8π^0.5) × Σ[w_j × (v(r_j, π) - v(r_j, -π)) / r_j^0.5]

E' = E (plane stress); E' = E/(1-ν²) (plane strain) v = crack-opening displacement (COD) Extrapolate K from multiple r-values; use 2-4 points; linear extrapolation to r → 0

More precise: use 3-parameter fit: K = α₀ + α₁√r + α₂r

J-Integral Domain Method

Domain (Virtual Crack Extension) Form

J = ∫∫_A [(W δ₁ⱼ - σᵢⱼ ∂uᵢ/∂x₁) ∂q/∂xⱼ] dA

q(x) = smooth weight function: q = 1 at crack tip, q = 0 on outer contour of domain A W = strain energy density

Advantages over contour integral:

• Domain integral is less sensitive to mesh quality
• Can use standard elements (no QPE needed) if domain is large enough
• Implemented in ABAQUS, ANSYS, MSC Nastran as "contour integral" request

Multiple contours (J_1, J_2, J_3): all should converge to same value; path-independence check Deviation > 5% between contours → mesh too coarse or near strong singularity

Separating K_I, K_II, K_III from J

J = (K_I² + K_II²)/E' + K_III²/(2μ) [combined mode]

Interaction integral (M-integral): K_I = E'/2 × M^(I) (auxiliary Mode I field) Decomposes mixed-mode J into individual K contributions

Virtual Crack Closure Technique (VCCT)

Concept: energy release rate G equals work to close crack from current state G_I = F_y^tip × v_tip / (2 × b × Δa)

F_y^tip = nodal force at crack tip node (mode I) v_tip = crack face opening behind tip b = element width (3D: out-of-plane thickness element behind tip) Δa = crack advance step

Advantages: simple; no re-meshing; works with standard elements Limitations: assumes self-similar crack growth; mode I only (unless extended)

Implementation (ABAQUS): DEBOND with VCCT or *FRACTURE CRITERION, TYPE=VCCT ANSYS: CINT command with VCCT option; results as G_I, G_II, G_III per element

Cohesive Zone Model (CZM)

Traction-separation law: Linear softening: T = K × δ (before damage); 0 (after) Bilinear: T rises to T_max at δ_0, drops to zero at δ_f Exponential (PPR, van den Bosch): physical for polymer fracture

Key parameters: G_Ic = area under T-δ curve = ½ × T_max × δ_f [J/m²] T_max = cohesive strength ≈ σ_ultimate of interface δ_f = critical opening ≈ 2 × G_Ic / T_max

Mixed-mode (Benzeggagh-Kenane criterion): G_total = G_Ic + (G_IIc - G_Ic) × (G_II/G_total)^η

Application: delamination in composites; adhesive joint failure; fatigue crack propagation

ABAQUS: COH3D8 / COH2D4 cohesive elements; or surface-based cohesion ANSYS: INTER204/205 interface elements; CZM material model TBDATA

Extended FEM (XFEM) in Commercial Codes

ABAQUS XFEM: *ENRICHMENT, TYPE=STATIONARY/PROPAGATION No conforming mesh to crack; crack tracked by level sets K_I, K_II, K_III extracted from J via interaction integral Automatic crack propagation in static (quasi-static) steps

ANSYS XFEM: SMART (Separating, Morphing, Adaptive, Remeshing Technology) Automatic crack growth direction from maximum hoop stress criterion Paris law for fatigue cycle-by-cycle crack growth

SIF Validation Benchmarks

Always validate against known solutions:

• Edge crack in tension: K_I = F_I × σ√(πa); F_I from Tada/Rooke-Cartwright
• Center crack: F_I = sec(πa/W)^0.5 (Isida) for W/a correction
• Compact Tension: K_I = (P/BW^0.5) × f(a/W); ASTM E399 expression Target: < 2% difference from handbook for identical geometry

3D vs. 2D Analysis

Plane stress: thin plates, free surface condition; K plane stress ≈ K 2D Plane strain: valid when B > 2.5(K_Ic/σ_y)² — gives lower bound K_Ic 3D: required for through-thickness variation (surface flaw, semi-elliptical crack) K varies along crack front: K = K(s) where s = position along crack front

Fatigue Crack Growth in FEA

Paris law integration: da/dN = C × (ΔK)^m K from FEA at each crack increment → cycle integration → life N_f NASGRO (AFGROW, FRANC3D): more accurate with threshold and R-ratio effects

Output

Provide: K_I, K_II, K_III [MPa√m] at key crack locations, J [kJ/m²], G_I, G_II, G_III [J/m²], path independence check (J contours), comparison to handbook solution, recommended element size at crack tip r_tip [mm], mesh type (QPE ring, XFEM, CZM), failure criterion (K_I vs. K_Ic, J vs. J_Ic), fatigue life N_f [cycles].