abaqus材料子程序
abaqus材料子程序
各向同性材料损伤本构模型
SUBROUTINE
UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,
+ RPL,DDSDDT,DRPLDE,DRPLDT,
+
STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMN AME,
+
NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEW DT,
+
CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC )
INCLUDE 'ABA_PARAM.INC'
CHARACTER*80 CMNAME
DIMENSION STRESS(NTENS),STATEV(NSTATV),
+ DDSDDE(NTENS,NTENS),DDSDDT(NTENS),
+ DRPLDE(NTENS),STRAN(NTENS),DSTRAN(NTENS), + TIME(2),PREDEF(1),DPRED(1),PROPS(NPROPS), + COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3)
DIMENSION STRANT(6),TSTRANT(4),PT(1)
DIMENSION OLD_STRESS(6)
DIMENSION DOLD_STRESS(6),D_STRESS(6)
DIMENSION
C(6,6),CD(6,6),DSTRESS(6),BSTRESS(6),ROOT(3),
+ DFMNDE(6),DDMDE(6),DCDDM(6,6),ATEMP1(6), ATEMP2(6)
PARAMETER
(ZERO=0.D0,ONE=1.D0,TWO=2.D0,FOUR=4.D0,HALF = 0.5D0) C start
C IF (NPROPS.LT.2) THEN
C WRITE(7,*) '** ERROR: UMAT REQUIRES *NPROPS=2'
C STOP
C EN
D IF
E11 =PROPS(1)
V12 =PROPS(2)
G12 =PROPS(1)/TWO/(ONE+PROPS(2))
C Critical values of stresses
XT=PROPS(3)
XC=PROPS(4)
XS=PROPS(5)
GX=PROPS(6) !Fracture energy in matrix ETA=0.001
C Current strain
DO I = 1, NTENS
STRANT(I) = STRAN(I) + DSTRAN(I)
END DO
C Stiffness
DO I = 1, 6
DO J = 1, 6
C(I,J)=ZERO
END DO
END DO
ATEMP = (1+V12)*(1-TWO*V12)
C(1,1) = E11*(1-V12)/ATEMP
C(2,2) = E11*(1-V12)/ATEMP
C(3,3) = E11*(1-V12)/ATEMP
C(1,2) = E11*V12/ATEMP
C(1,3) = E11*V12/ATEMP
C(2,3) = E11*V12/ATEMP
C(4,4) = G12
C(5,5) = G12
C(6,6) = G12
DO I = 2, 6
DO J = 1, I-1
C(I,J) = C(J,I)
END DO
END DO
C Critical values of strains
XET=XT/(C(1,1)-2*V12*C(1,2))
XEC=XC/(C(1,1)-2*V12*C(1,2))
XES=XS/C(4,4)
DMOLD = STATEV(1)
C Strain initiation criterion
A11 = STRANT(1)**TWO+STRANT(2)**TWO+STRANT(3)**TWO
A12 = A11 / XET / XEC
A21 = STRANT(1)+STRANT(2)+STRANT(3)
A22 = (XEC - XET) / XEC / XET * A21
A31 = STRANT(4)**TWO+STRANT(5)**TWO+STRANT(6)**TWO
A32 = A31 / XES**TWO
A1= A12 + A22 + A32
C B11 = STRANT(2)**TWO
C B12 = B11 / XET / XEC
C B21 = STRANT(2)
C B22 = (XEC - XET) / XEC / XET * B21 C B31 = STRANT(5)**TWO
C B32 = B31 / XES**TWO
C B1= B12 + B22 + B32
C C11 = STRANT(3)**TWO
C C12 = C11 / XET / XEC
C C21 = STRANT(3)
C C22 = (XEC - XET) / XEC / XET * C21 C C31 = STRANT(6)**TWO
C C32 = C31 / XES**TWO
C C1= C12 + C22 + C32
STATEV(2)=A1
C STATEV(3)=B1
C STATEV(4)=C1
FMN = ZERO
IF (A1.GT.ZERO) THEN
FMN =SQRT(A1)
C IF (B1.GT.ONE) THEN
C FMN =FMN+SQRT(B1)
C IF(C1.GT.ONE) THEN
C FMN =FMN+SQRT(C1)
C EN
D IF
C EN
D IF
END IF
STATEV(5)=FMN
C write(*,*) FMN
DM = ZERO
DDMDFMN = ZERO
DO I = 1, 6
DFMNDE(I) = ZERO
DDMDE(I) = ZERO
END DO
IF (FMN .GT. ONE) THEN
C CALCULATE DM, DDMDFMN
C WRITE(6,*)FMN
T1 = (C(1,1)-2*V12*C(1,2)) * XET**2 * CELENT / GX
T2 = (ONE - FMN) * T1
DM = ONE - EXP(T2)/FMN
C WRITE(6,*)'T1 ',T1,' T2', T2, ' DM', DM
C write(*,*) DM
C CALCULATE THE DERIVATIVE OF DAMAGE VARIABLE WITH RESPECT TO FAILURE
C RITERION
DDMDFMN = (ONE / FMN + T1) * (ONE - DM)
C CALCULATE DFMNDE
IF (DM .GT. DMOLD) THEN
DFMNDE(1) = HALF/FMN*(TWO*STRANT(1)+XEC-XET)/XET/XEC
DFMNDE(2) = HALF/FMN*(TWO*STRANT(2)+XEC-XET)/XET/XEC
DFMNDE(3) = HALF/FMN*(TWO*STRANT(3)+XEC-XET)/XET/XEC
DFMNDE(4) = ONE/FMN*TWO*STRANT(4)/XES**TWO
DFMNDE(5) = ONE/FMN*TWO*STRANT(5)/XES**TWO
DFMNDE(6) = ONE/FMN*TWO*STRANT(6)/XES**TWO
DO I = 1, 6
DDMDE(I) = DFMNDE(I) * DDMDFMN END DO
END IF
END IF
DM = MAX (DM, DMOLD)
C write(6,*) DM
C SAVE THE OL
D STRESS TO OLD_STRESS
DO I = 1, NTENS
OLD_STRESS(I) = STRESS(I)
END DO
C Effective stiffness
DO I = 1, 6
DO J = 1, 6
CD(I,J)=C(I,J)
END DO
END DO
IF (DM.NE.ZERO) THEN
CD(1,1) = (ONE - DM)*C(1,1)
CD(1,2) = (ONE - DM)*C(1,2)
CD(2,1) = CD(1,2)
CD(2,2) = (ONE - DM)*C(2,2)
CD(1,3) = (ONE - DM)*C(1,3)
CD(2,3) = (ONE- DM)*C(2,3) CD(3,2) = CD(2,3)
CD(4,4) = (ONE - DM)*C(4,4) CD(5,5) = (ONE - DM)*C(5,5) CD(6,6) = (ONE - DM)*C(5,5) END IF
C Elastic derivative
DO I = 1, 6
DO J = 1, 6
DCDDM(I,J) = ZERO
END DO
END DO
C
C CALCULATE DC/DDM
C
DCDDM(1,1) = -C(1,1)
DCDDM(1,2) = -C(1,2)
DCDDM(2,1) = -C(2,1)
DCDDM(2,2) = -C(2,2)
DCDDM(2,3) = -C(2,3)
DCDDM(3,3) = -C(3,3)
DCDDM(3,1) = -C(3,1)
DCDDM(1,3) = -C(1,3)
DCDDM(4,4) = -C(4,4)
DCDDM(5,5) = -C(5,5)
DCDDM(6,6) = -C(6,6)
C UPDATE THE JACOBIAN
DO I = 1, NTENS
ATEMP1(I) = ZERO
DO J = 1, NTENS
ATEMP1(I) = ATEMP1(I) + DCDDM(I,J) * STRANT(J)
END DO
END DO
DDSDDE=0
DO I = 1, NTENS
DO J = 1, NTENS
DDSDDE(I,J)=CD(I,J)+(ATEMP1(I)*DDMDE(J))*DTIME/(DTI ME+ETA)
END DO
END DO
C Update stresses
DO I = 1, NTENS
STRESS(I)=ZERO
DO J = 1, NTENS
C IF(DM.LT.0.5) THEN
STRESS(I)=STRESS(I)+ CD(I,J) * STRANT(J) C ELSE
C STRESS(I)=STRESS(I)+ CD(I,J) * STRANT(J) * (1-DM)
C ENDIF
END DO
END DO
C Energy
DO I = 1, NDI
SSE = SSE + HALF * (STRESS(I) + OLD_STRESS(I)) * DSTRAN(I)
END DO
DO I = NDI+1, NTENS
SSE = SSE + (STRESS(I) + OLD_STRESS(I)) * DSTRAN(I)
END DO
STATEV(1)=DM
RETURN
END