2D Stationary Crack Analysis: J-integral and Stress Intensity Factor Calculation

We talked about calculating J-integral and K1 using XFEM in the last post. XFEM is a powerful method as it eliminates the need to remesh crack-tip regions. However, there are some limitations. For example, material behaviour is assumed to be linearly elastic. The available fracture criteria are valid only for cracks in homogeneous linear elastic materials. The analysis is assumed to be quasi-static. Other analysis types are not available in ANSYS. Therefore, it is still pertinent to explore the common method where KSCON command is used to model stress concentration at the crack tip. It is worth mentioning that KSCON command does not support 3D modelling. For 3D problems, cube elements with mid-point nodes can be used. By collapsing one face of the cube, the cube degenerates to a prism element and a similar effect as that of KSCON can be achieved. 

In this example, we are going to solve a similar 2D problem as the one in the last post. A single-edge cracked plate is under uniform tension. As the crack lies along the symmetric axis, only half of the plate is modelled. The results of this example were compared with those obtained using XFEM and good agreement was observed. With this common method, it is possible to calculate J-integral and stress intensity factors even when the material behaviour is elastic-plastic. 

111000

Fig.1. Mesh refinement at the crack tip.

The APDL is attached here (annotations are in lowercase):

FINISH
/CLEAR
/TITLE,STATIONARY CRACK ANALYSIS
/PREP7
WIDTH=10
LENGTH=20
CRLENGTH=5			!crack length
NUM_CONTOUR=10		!the number of contours
ET,1,183
KEYOPT,1,3,2		!plane strain
MP,EX,1,3E4
MP,NUXY,1,0.3
MP,DENS,1,1
!define keypoints:
SELTOL,1E-8		!set the tolerance
K,1
K,2,CRLENGTH,0
K,3,WIDTH,0
K,4,WIDTH,LENGTH/2
K,5,0,LENGTH/2
A,1,2,3,4,5
!define crack tip singularity:
KSCON,2,WIDTH/420,1,14,0.75
ESIZE,WIDTH/120
AMESH,ALL
!boundary condition: top face	
NSEL,S,LOC,Y,LENGTH/2
SF,ALL,PRES,-100
ALLSEL
LSEL,S,,,2
NSLL,S,1
DSYM,SYMM,Y,0		!symmetric boundary
ALLSEL
NSEL,S,LOC,X,WIDTH
NSEL,R,LOC,Y,0
D,ALL,ALL
ALLSEL

/SOLU
ANTYPE,0
TIME,1
DELTIM,0.1,1E-1,0.2
OUTRES,ALL,ALL
!local coordinate system:
LOCAL,11,0,0,0,0,,		!the origin of this local coordinate system does not have to be the crack tip, but its axis direction is important
!select nodes located along the crack front and define it as crack front/tip node component
NSEL,S,LOC,X,CRLENGTH
NSEL,R,LOC,Y,0
CM,CRACK_TIP_NODE_CM,NODE
ALLSEL
! Define J-integral calculation:
CINT,NEW,1
CINT,TYPE,JINT
CINT,CTNC,CRACK_TIP_NODE_CM
CINT,NORM,11,2		!define the normal direction of crack plane
CINT,NCON,NUM_CONTOUR
CINT,SYMM,ON   !symmetrical crack
! Define K calculation:
CINT,NEW,2       ! initiate a new calculation as #2
CINT,TYPE,SIFS   ! specify stress-intensity factor calculations 
CINT,CTNC,CRACK_TIP_NODE_CM
CINT,NORM,11,2
CINT,NCON,NUM_CONTOUR
CINT,SYMM,ON   !symmetrical crack
OUTRES,CINT,ALL
SOLVE
FINISH
/POST1
SET,LAST,LAST
/OUT
/COM ****** RESULTS ******
/COM
/COM ****** PRINT NODAL RESULTS ******
/COM
/COM >>> JINTEGRAL
/COM
PRCINT,1,,JINT
/COM
/COM >>> MODE 1 STRESS INTENSITY FACTOR
/COM
PRCINT,2,,K1
/COM
/COM >>> MODE 2 STRESS INTENSITY FACTOR
/COM
PRCINT,2,,K2
/COM
PLNSOL,S,Y
/CVAL,1,0,500,1000,2000,4000,5000,6000,10000 
/REPLOT
!export J-integral results:
*DIM,JDATA,,NUM_CONTOUR
*DO,I,1,NUM_CONTOUR
*GET,J_CONTOUR,CINT,1,CTIP,2,CONTOUR,I,DTYPE,JINT
JDATA(I)=J_CONTOUR
*ENDDO
*CREATE,JOUTPUT,MAC
*CFOPEN,JRESULTS,TXT
*VWRITE,JDATA(1)
(F10.6)
*CFCLOS
*END
JOUTPUT

JKresult

Fig.2. The results of J-integral and K1.

111001

Fig.3. Nodal stress in Y direction.

One thought on “2D Stationary Crack Analysis: J-integral and Stress Intensity Factor Calculation

  1. Pingback: 3D Stationary Crack Analysis: J-integral and Stress Intensity Factor Calculation | Xutao Sun

Leave a comment