450 likes | 503 Views
Complex Variable Finite Element Methods for Fracture Mechanics Analysis. Harry Millwater, Ph.D. David Wagner, MS Jose Garza, MS Andrew Baines, BS Kayla Lovelady, BS Carolina Dubinsky, MS Thomas Ross, MS Sivirt Student Seminary, April 2013. Overview.
E N D
Complex Variable Finite Element Methods for Fracture Mechanics Analysis Harry Millwater, Ph.D. David Wagner, MS Jose Garza, MS Andrew Baines, BS Kayla Lovelady, BS Carolina Dubinsky, MS Thomas Ross, MS Sivirt Student Seminary, April 2013
Overview • Development of the complex variable finite element method (ZFEM) • Derivatives wrt shape, material prop, load and others • Implemented into the Abaqus finite element code • Verification against known solutions • Applications to fracture mechanics • Future efforts
Finite Difference Method Im Perturb along the imaginary axis Finite Differencing Re h Forward Differencing Determining h is problematic
Complex Taylor Series Expansion Im Perturb along the imaginary axis Subset of Fourier Differentiation h Re h can be “very” small ~ 10-30
Example: Finite Difference CTSE No subtraction of nearly equal numbers h can be as small as desired Exact derivative can be recovered as h becomes small
Finite Element Analysis • Shape sensitivities “difficult” to compute using finite difference method • Mesh must be perturbed enough to “see” the perturbation but not too large to distort elements h must be kept small
Complex Variable Advantage Perturb nodal coordinates in imaginary axis only! Minimal mesh distortion issues!
Computational Issues Standard finite element EXCEPT complex nodal coordinates Stiffness matrix is now complex->COMPLEX solver required Apply an imaginary displacement to the nodal coordinates to represent a shape (domain) change Shape change is arbitrary Sensitivities for entire displacement/strain/stress available Imaginary displacement
Example: 1st Der. of σθ w.r.t. Ri Analytic Solution CTSE UTSA Matlab code • A. Voorhees, H.R. Millwater, R.L. Bagley, “Complex Variable Methods for Shape Sensitivity of Finite Element Models,” Finite Elem. Anal. Des., 47 (2011) 1146–1156, doi:10.1016/j.finel.2011.05.003
L2 Norm – Error controlled by Mesh L2 norm reduces as mesh is refined
Current Effort – Implementation into Abaqus Abaqus user element implementation (uel) 6 dof/node (3 real, 3 imag) Abaqus cannot solve complex stiffness matrix, represent as real n x n complex matrix solved as 2n x 2n real matrix
Abaqus Example Plane Strain ½ P ½ P Real nodes (physical coordinates) Imaginary nodes (perturbation in physical coordinates Abaqus UEL (user element) Plane stress/strain Element has 8 nodes (4 real, 4 imag)
Abaqus UEL Example Code COMPLEX(KIND=r8),DIMENSION(MCRD ,NNODE/2 ) :: zCOORDS COMPLEX(KIND=r8),DIMENSION(NNODE,NNODE)::zAMATRX !Complex Stiffness Matrix !Construct complex nodal coordinate array zCOORDS = CMPLX(COORDS(1:2,1:3), COORDS(1:2,NNODE/2+1:NNODE)*h) f1 = zCOORDS(1,2)*zCOORDS(2,3) - zCOORDS(1,3)*zCOORDS(2,2) f2 = zCOORDS(1,3)*zCOORDS(2,1) - zCOORDS(1,1)*zCOORDS(2,3) f3 = zCOORDS(1,1)*zCOORDS(2,2) - zCOORDS(1,2)*zCOORDS(2,1) b1 = zCOORDS(2,2) - zCOORDS(2,3) b2 = zCOORDS(2,3) - zCOORDS(2,1) b3 = zCOORDS(2,1) - zCOORDS(2,2) c1 = zCOORDS(1,3) - zCOORDS(1,2) c2 = zCOORDS(1,1) - zCOORDS(1,3) c3 = zCOORDS(1,2) - zCOORDS(1,1) Area = half * (f1 + f2 + f3) zAMATRX = MATMUL(Bt, MATMUL(Cthick,zB_ip) ) * Area
Obtaining Derivatives • Shape sensitivity – input imaginary coordinates to represent shape change. (All Imag nodes not perturbed have a value of zero) Nodal Coordinates 1 (0,0) 2 (1,0) 3 (0,1) 4 (1,1) 5 (0,0) 6 (0,0) 7 (0,1)h 8( 0,1)h ½ P ½ P 7 3 8 4 L h = 10-10 1 5 2 6
Obtaining Derivatives • Load sensitivity – Apply perturbation in loading to IM nodes. Applied Loads 3 (0,1/2 P) 4 (0,1/2 P) 7 (0,1/2)h 8( 0,1/2)h ½ P ½ P 7 3 8 4 Imaginary coordinates all zero h = 10-10 1 5 2 6
Obtaining Derivatives • Material sensitivity – Apply perturbation to constitutive matrix. ½ P ½ P 7 3 8 4 h = 10-10 Imaginary coordinates all zero 1 5 2 6
Abaqus 2D Verification Analytical Abaqus
Weight Function Development Calculation of the partial derivative of the crack opening displacement with respect to crack length required • Standard research approaches: • Assume 3-4 term approximations to du/da or approximate the weight fn. directly. Use multiple reference solutions to solve for Mi
Perturbation of Crack Length Crack tip element can be perturbed in the imaginary domain – no perturbation of real mesh Perturb a no. of elements around the crack tip Perturbation of Crack Length in Imaginary domain
Accuracy: Infinite Array of Cracks Crack Opening Displacement 67 – 181 nodes on crack line Accurate calculation of the crack opening displacement with respect to crack length all along the crack line
Infinite Array Weight Function Accuracy controlled by the user: high order expansions possible 3-term[reference] 7-term [WCTSE] 3-term[WCTSE]
Accuracy: Crack from a Hole WCTSE ±0.15% WCTSE ±0.06% WCTSE weight function consistently better than published weight fns. D. Wagner, and H. Millwater, “2D Weight Function Development using a Complex Taylor Series Expansion Method,” Engng Fract Mech 86 (2012), 23-37
Strain Energy Release Rate Double Cantilever Beam (DBC) Energy release rate (Exact) = 0.800 20 Contours 8 4 2
Strain Energy Release Rate Double Cantilever Beam (DBC) Exact = 0.8 Exact = 0.800
Implementation into Abaqus **[…] ** ** ** **===================== Real Mesh ======================= *NODE, NSET=rnodes 1, 0.0d0, 0.0d0 2, 2.0d0, 0.0d0 3, 1.2d0, 1.0d0 4, 0.2d0, 1.0d0 5, 1.0d0, 0.0d0 6, 1.6d0, 0.5d0 7, 0.6d0, 1.0d0 8, 0.1d0, 0.5d0 ** ** ** ** ** ** ** ** ** ** *ELEMENT, TYPE=CPE8, ELSET=relems 1, 1, 2, 3, 4, 5, 6, 7, 8 ** **===================== Material ======================== *SOLID SECTION, ELSET=relems, MATERIAL=aluminum 1d0 *MATERIAL, NAME=aluminum *ELASTIC 100d9, 0.3d0 **======== Node Sets, Element Sets, and Surfaces ======== *NSET, NSET=rbotedge 1 , 5, 2 ** ** **================ Boundary Conditions ================== *BOUNDARY rbotedge, 2 1, 1 ** ** **==================== Load Step ======================== *STEP, NAME=Step-1 *STATIC , DIRECT 1d0, 1d0, 1d-05, 1d0 **[…] **[…] *USER ELEMENT, TYPE=U28, NODES=16, COORDINATES=2, UNSYMM, PROPERTIES=5, IPROPERTIES=2, VARIABLES=108 1, 2 **===================== Real Mesh ======================= *NODE, NSET=rnodes 1, 0.0d0, 0.0d0 2, 2.0d0, 0.0d0 3, 1.2d0, 1.0d0 4, 0.2d0, 1.0d0 5, 1.0d0, 0.0d0 6, 1.6d0, 0.5d0 7, 0.6d0, 1.0d0 8, 0.1d0, 0.5d0 **==================== Complex Mesh ===================== *NODE, NSET=inodes 101, 0.0d0, 0.0d0 102, 0.0d0, 0.0d0 103, 1.0d0, 0.0d0 104, 0.0d0, 0.0d0 105, 0.0d0, 0.0d0 106, 0.5d0, 0.0d0 107, 0.5d0, 0.0d0 108, 0.0d0, 0.0d0 *ELEMENT, TYPE=U28, ELSET=zelems 101, 1, 2, 3, 4, 5, 6, 7, 8, 101,102,103,104,105,106,107,108 **===================== Material ======================== ** ** ** *UEL PROPERTY, ELSET=zelems 100d9, 0.3d0, 1d-9, 1d0, 1d0, 3, 0 **======== Node Sets, Element Sets, and Surfaces ======== *NSET, NSET=rbotedge 1, 5, 2 *NSET, NSET=ibotedge 101,105,102 **================ Boundary Conditions ================== *BOUNDARY rbotedge, 2 1, 1 ibotedge, 2 101, 1 **==================== Load Step ======================== *STEP, NAME=Step-1, UNSYMM=YES *STATIC, DIRECT 1d0, 1d0, 1d-5, 1d0 **[…] Abaqus Input File Changes Needed Imaginary nodes UEL definition UEL properties
Multicomplex Mathematics Bi-complex numbers representation : Matrix representation of bi-complex numbers: Tri-complex numbers representation :
Bicomplex Analysis Higher Order Sensitivities Thick walled cylinder
Progressive Fracture • Construct a 3rd order Taylor series of the strain energy using tricomplex elements in front of the crack tip. • Predict the crack path along the max energy release • Progress the crack, remesh, repeat
Progressive Fracture • Construct a 3rd order Taylor series of the strain energy using tricomplex elements. Predict the crack path along the max energy release • Progress the crack, remesh, repeat
Progressive Fracture Similar results to Franc2D with far fewer steps
ZFEM User Element Features • “Easy” to Use • Only a Single Input File Needed • Single User Subroutine File • No Abaqus Configuration Needed • Analytic Extensions of Built-in Elements • Identical Response of Real values • 2D (6,8 noded), 3D (8, 15, 20) • Modular Fortran 95 Implementation • Extensible • Very Good Performance abaqus job=mycomplexmodel user=zfem
Future Interests • High order sensitivities: • Multicomplex mathematics element library (complete for linear elastic statics) • Extension to non-linear materials • Plasticity enhancement (in progress) • Composites, anisotropic materials • 3D progressive fracture (in progress) • Structural dynamics (in progress) • Thermoelastic analysis • Expanded element library • Plates and shells
Acknowledgements • Efficient Sensitivity Methods for Probabilistic Lifing and Engine Prognostics, Pat Golden, AFRL/RXLMN, Aug. 2007-Sep. 2010 • Efficient Finite Element-based 3D Fracture Mechanics Crack Growth Analysis using Complex Variable Sensitivity Methods, DoD PETTT, Sep. 2010 - Aug. 2011 • Implementation of Complex Variable Finite Element Methods in Abaqus, DOD PETTT, Sep. 2011- Aug. 2012 • Enhanced Fracture Mechanics Crack Growth Analysis using Complex Variable Sensitivity Methods, AFOSR (David Stargel), May 2011-2014