SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD, 1 RPL,DDSDDT,DRPLDE,DRPLDT, 2 STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME, 3 NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT, 4 CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC) C INCLUDE 'ABA_PARAM.INC' C CHARACTER*80 CMNAME DIMENSION STRESS(NTENS),STATEV(NSTATV), 1 DDSDDE(NTENS,NTENS),DDSDDT(NTENS),DRPLDE(NTENS), 2 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1), 3 PROPS(NPROPS),COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3) c ****Ende des Headers**** c-----DEFINITION DER VARIABLEN--------------------------------------------------------------------- c DOUBLE PRECISION BI1, BI2, SELA(6), SCOM(6) DIMENSION B(3,3),DISTGR(3,3), A(3,3) DIMENSION UI1(4), UI2(4) DOUBLE PRECISION H, CS(6,6), BZ DOUBLE COMPLEX, PARAMETER :: i =(0.0D0,1.0D0) DOUBLE COMPLEX :: BI1UNIC, BI1BIAXC, WURZEL, BI2BIAXC, WURZEL2 DOUBLE COMPLEX :: BI2UNIC, BI1UNIC2 DOUBLE PRECISION BI1PS, BI1UNI(3),BI1BIAX(3) DOUBLE PRECISION BI2UNI(3), BI2BIAX(3), D1, NU, H1 INTEGER O, J, K, L, M, N INTEGER KOEFF(6,2) DOUBLE PRECISION D(3,3), DSDEEL(6,6), DSDECO(6,6) DATA KOEFF /1,2,3,1,1,2,1,2,3,2,3,3/ DATA D /1.D0,0.D0,0.D0,0.D0,1.D0,0.D0,0.D0,0.D0,1.D0/ c ****Erweiterung für YEOH-Modell**** DOUBLE PRECISION C10int, C20int, C30int c c-----MATERIALPARAMETER--------------------------------------------------------------------------- c c ****Materialparameter der reinen Beanspruchungszustände für Yeoh-Modell**** c c ****Uniaxial**** C10int = 2.8 C20int = 0.D0 C30int = 0.0123 c ****Querkontraktionszahl**** NU=0.5 c-----AUFSTELLEN DES DEFORMATIONSTENSORS----------------------------------------------------------- c c ****Berechnung des Deformationsgradienten DISTGR(i,j)**** c c ****Skalierung zur Elimination der Kompressibilität**** c DET=DFGRD1(1,1)*DFGRD1(2,2)*DFGRD1(3,3) 1 +DFGRD1(1,2)*DFGRD1(2,3)*DFGRD1(3,1) 2 +DFGRD1(1,3)*DFGRD1(2,1)*DFGRD1(3,2) 3 -DFGRD1(3,1)*DFGRD1(2,2)*DFGRD1(1,3) 4 -DFGRD1(3,2)*DFGRD1(2,3)*DFGRD1(1,1) 5 -DFGRD1(3,3)*DFGRD1(2,1)*DFGRD1(1,2) c SCALE=DET**(-1.D0/3.D0) DISTGR(1, 1)=SCALE*DFGRD1(1,1) DISTGR(2, 2)=SCALE*DFGRD1(2,2) DISTGR(3, 3)=SCALE*DFGRD1(3,3) DISTGR(1, 2)=SCALE*DFGRD1(1,2) DISTGR(1, 3)=SCALE*DFGRD1(1,3) DISTGR(2, 1)=SCALE*DFGRD1(2,1) DISTGR(2, 3)=SCALE*DFGRD1(2,3) DISTGR(3, 1)=SCALE*DFGRD1(3,1) DISTGR(3, 2)=SCALE*DFGRD1(3,2) c c c c ****Berechnung der ersten und zweiten Invarianten aus dem linken Cauchy-Green Tensor**** BI1 = B(1,1)+B(2,2)+B(3,3) BI2 = B(1,1)*B(2,2)+B(2,2)*B(3,3) 1 +B(1,1)*B(3,3) 2 -B(1,2)**2-B(2,3)**2-B(1,3)**2 c c ****Berechnung des linken Cauchy-Green Tensors**** c B(1,1)=DISTGR(1, 1)**2+DISTGR(1, 2)**2+DISTGR(1, 3)**2 B(2,2)=DISTGR(2, 1)**2+DISTGR(2, 2)**2+DISTGR(2, 3)**2 B(3,3)=DISTGR(3, 3)**2+DISTGR(3, 1)**2+DISTGR(3, 2)**2 B(1,2)=DISTGR(1, 1)*DISTGR(2, 1)+DISTGR(1, 2)*DISTGR(2, 2) 1 +DISTGR(1, 3)*DISTGR(2, 3) B(1,3)=DISTGR(1, 1)*DISTGR(3, 1)+DISTGR(1, 2)*DISTGR(3, 2) 1 +DISTGR(1, 3)*DISTGR(3, 3) B(2,3)=DISTGR(2, 1)*DISTGR(3, 1)+DISTGR(2, 2)*DISTGR(3, 2) 1 +DISTGR(2, 3)*DISTGR(3, 3) B(3,2)=B(2,3) B(3,1)=B(1,3) B(2,1)=B(1,2) c ****Matrixmultiplikationen**** A(1,1) = B(1,1)*B(1,1)+B(1,2)*B(2,1)+B(1,3)*B(3,1) A(2,2) = B(2,1)*B(1,2)+B(2,2)*B(2,2)+B(2,3)*B(3,2) A(3,3) = B(3,1)*B(1,3)+B(3,2)*B(2,3)+B(3,3)*B(3,3) A(1,2) = B(1,1)*B(1,2)+B(1,2)*B(2,2)+B(1,3)*B(3,2) A(2,1) = B(2,1)*B(1,1)+B(2,2)*B(2,1)+B(2,3)*B(3,1) A(1,3) = B(1,1)*B(1,3)+B(1,2)*B(2,3)+B(1,3)*B(3,3) A(3,1) = B(3,1)*B(1,1)+B(3,2)*B(2,1)+B(3,3)*B(3,1) A(2,3) = B(2,1)*B(1,3)+B(2,2)*B(2,3)+B(2,3)*B(3,3) A(3,2) = B(3,1)*B(1,2)+B(3,2)*B(2,2)+B(3,3)*B(3,2) c-----BILDUNG DER ABLEITUNGEN DER FORMÄNDERUNGSENERGIEDICHTE U------------------------------------- c c ****Ableitung von U nach BI1**** c UI1(1) = C10(1)+2.D0*C20*BI1+3.D0*C30*BI1**2.D0 UI1(1) = C10int+2*C20int*(BI1-3) F +3*C30int*(BI1-3)**2 c ****Ableitung von U nach BI2**** UI1(2) = 0 c c ****Zweite Ableitung von U nach BI1**** c UI2(1) = 2.D0*C20+6.D0*C30*BI1 UI2(1) = 2*C20int+6*C30int*(BI1-3) c c ****Zweite Ableitung von U nach BI2**** UI2(2) = 0 c ****Zweite Ableitung von U nach BI1 und BI2**** UI2(4) = 0 c c c-----BERECHNUNG DES SPANNUNGSTENSORS-------------------------------------------------------------- c c c ****Hyperelastischer Anteil des Spannungstensors**** SELA(1) = 2.D0/DET*(B(1,1)*(UI1(1) + BI1*UI1(2)) 1 - UI1(2)*A(1,1)) SELA(2) = 2.D0/DET*(B(2,2)*(UI1(1) + BI1*UI1(2)) 1 - UI1(2)*A(2,2)) SELA(3) = 2.D0/DET*(B(3,3)*(UI1(1) + BI1*UI1(2)) 1 - UI1(2)*A(3,3)) SELA(4) = 2.D0/DET*(B(1,2)*(UI1(1) + BI1*UI1(2)) 1 - UI1(2)*A(1,2)) SELA(5) = 2.D0/DET*(B(1,3)*(UI1(1) + BI1*UI1(2)) 1 - UI1(2)*A(1,3)) SELA(6) = 2.D0/DET*(B(2,3)*(UI1(1) + BI1*UI1(2)) 1 - UI1(2)*A(2,3)) c c ****Kugelanteil abziehen**** H = (SELA(1) + SELA(2) + SELA(3)) / 3.D0 SELA(1) = SELA(1) - H SELA(2) = SELA(2) - H SELA(3) = SELA(3) - H c ****Spannungsanteil aus Kompressibilität (linearer Ansatz)**** c c ****Berechnung des Kompressionsmoduls aus der Querkontraktionszahl NU**** c H1 = (2*NU-2)/(6*NU-3+0.00000000001) D1 = 1/(H1*C10int) c H = 2.D0/D1 * (DET - 1.D0) SCOM(1) = H SCOM(2) = H SCOM(3) = H SCOM(4) = 0.D0 SCOM(5) = 0.D0 SCOM(6) = 0.D0 c c c-----BERECHNUNG DER TANGENTENSTEIFIGKEITSMATRIX--------------------------------------------------- c c c ****Hyperelastischer Anteil der Tangentensteifigkeitsmatrix**** c DO 107 M = 1, 6 DO 107 N = 1, 6 O = KOEFF(M,1) J = KOEFF(M,2) K = KOEFF(N,1) L = KOEFF(N,2) CS (M,N) = 2.D0/DET * (UI1(1) + BI1*UI1(2)) * F (D(O,K)*B(J,L) + B(O,K)*D(J,L)) - F 2.D0/DET * UI1(2) * F (D(O,K)*A(J,L)+A(O,K)*D(J,L)+2.D0*B(O,K)*B(J,L)) + F 4.D0/DET * (UI2(1) + UI1(2)+ 2.D0*BI1*UI2(4) + F BI1*BI1*UI2(2)) * B(O,J) * B(K,L) - F 4.D0/DET * (UI2(4) + BI1*UI2(2)) * F (A(O,J)*B(K,L) + B(O,J)*A(K,L)) + F 4.D0/DET * UI2(2) * A(O,J) * A(K,L) 107 CONTINUE DO 108 M = 1, 6 DO 108 N = 1, 6 O = KOEFF(M,1) J = KOEFF(M,2) K = KOEFF(N,1) L = KOEFF(N,2) DSDEEL(M,N) = CS(M,N) - F D(O,J)/3.D0 * (CS(1,N) + CS(2,N) + CS(3,N)) - F D(K,L)/3.D0 * (CS(M,1) + CS(M,2) + CS(M,3)) + F D(O,J)*D(K,L)/9.D0 * (CS(1,1) + CS(2,1) + CS(3,1) + F CS(1,2) + CS(2,2) + CS(3,2) + F CS(1,3) + CS(2,3) + CS(3,3) ) 108 CONTINUE c ****Anteil der Tangentensteifigkeitsmatrix aus Kompressibilität (linearer Ansatz)**** H = 2.D0/D1 * (2.D0*DET - 1.D0) DO 2010 M = 1, 3 DO 2015 N = 1, 3 DSDECO(M,N) = H 2015 CONTINUE DO 2010 N = 4, 6 DSDECO(M,N) = 0.D0 2010 CONTINUE DO 2020 M = 4, 6 DO 2020 N = 1, 6 DSDECO(M,N) = 0.D0 2020 CONTINUE c ****Hyperelastische Anteile und Kompressionsanteil zusammensetzten**** DO 2030 O = 1, 6 STRESS(O) = SELA(O) + SCOM(O) DO 2030 J = 1, 6 DDSDDE(O,J) = DSDEEL(O,J) + DSDECO(O,J) c IF (O.GT.3 .OR. J.GT.3) DDSDDE(O,J) = 0.5D0*DDSDDE(O,J) 2030 CONTINUE RETURN END