1 : C INCREP SOURCE KICH 05/11/09 21:17:40 5232 2 : C INCREP SOURCE BV 03/11/04 21:29:05 1995 3 : SUBROUTINE INCREP(SIG,VAR,EPSVPT,VARPT,XMAT,NSTRS0, 4 : & MFR,NVARI,NCOMAT) 5 : IMPLICIT REAL*8(A-H,O-Z) 6 : DIMENSION SIG(*),VAR(*),EPSVPT(*),VARPT(*),XMAT(*) 7 : DIMENSION XX(6),EPS0(6),XP(3,3) 8 : DIMENSION ST(3,3),D(3),ZMU(3,3) 9 : DIMENSION TEMP(3),EPSD(3),VPT(6) 10 : PARAMETER (AMAX = 1.0D20 , AMIN = 1.D-10) 11 : DETIER = 2.0D0/3.0D0 12 : ROOT = SQRT(DETIER) 13 : C-------------------------------------------------------------------| 14 : C******* EVALUATION OF J2 SIGMA | 15 : C-------------------------------------------------------------------I 16 : TRACE =(SIG(1)+SIG(2)+SIG(3))/3.0D0 17 : DO 70 I=1,NSTRS0 18 : A = 0.0D0 19 : IF (I.LE.3) A=1.0D0 20 : XX(I) = SIG(I)-A*TRACE 21 : 70 CONTINUE 22 : AJ2 = PROCON (XX,XX,NSTRS0) 23 : AJ2 = SQRT(1.5D0*AJ2) 24 : C------------------------------------------------------ 25 : C******* EVALUATION OF J0 CONTRAINTE PRINCIPALE MAX 26 : C------------------------------------------------------ 26 : CALL ZERO(ST,3,3) 27 : DO 50 I=1,3 28 : ST(I,I) = SIG(I) 29 : 50 CONTINUE 30 : ST(1,2) = SIG(4) 31 : IF(NSTRS0.NE.6) GOTO 54 32 : ST(1,3) = SIG(5) 33 : ST(2,3) = SIG(6) 34 : 54 CONTINUE 35 : ST(2,1) = ST(1,2) 36 : ST(3,1) = ST(1,3) 37 : ST(3,2) = ST(2,3) 37 : CALL ZERO(XP,3,3) 38 : * 38 : CALL JACOB3(ST,3,D,XP) 39 : C 40 : AJ0 = D(1) 41 : * AJ0 = DMAX1(D(1),D(2),D(3)) 42 : C-------------------------------------------------------- 43 : C******* CALCUL DE LA CONTRAINTE EQUIVALENTE AU SENS 44 : C******* DE L'ENERGIE DE CISAILLEMENT (VON MISES) 45 : C-------------------------------------------------------- 46 : DD = VAR(3) 47 : DC = VAR(4) 48 : IF (DD.GT.0.98D0) THEN 49 : DD = 0.98D0 50 : ENDIF 51 : IF (DC.GT.0.98D0) THEN 52 : DC = 0.98D0 53 : ENDIF 54 : ZMU(1,1) = 2.*(1.+(XMAT(2)*(1.-DC))) 55 : ZMU(1,2) = (DC / 2.0D0)-(1.0D0+XMAT(2)*(1.0D0-DC)) 56 : ZMU(1,3) = ZMU(1,2) 57 : ZMU(2,2) = 2.0D0*(1.0D0+XMAT(2))*(1.0D0-DC) 58 : ZMU(2,3) = -1.0D0*(1.0D0+XMAT(2))*(1.0D0-DC) 59 : ZMU(3,3) = ZMU(2,2) 60 : ZMU(2,1) = ZMU(1,2) 61 : ZMU(3,1) = ZMU(1,3) 62 : ZMU(3,2) = ZMU(2,3) 63 : DO 18 J=1,3 64 : DO 17 I=1,3 65 : ZMU(I,J) = ZMU(I,J)/(2.0D0*(1.0D0+XMAT(2))) 66 : 17 CONTINUE 67 : 18 CONTINUE 68 : TEMP(1) = 0.0D0 69 : TEMP(2) = 0.0D0 70 : TEMP(3) = 0.0D0 71 : TEMP(1) = ZMU(1,1)*D(1)+ZMU(1,2)*D(2)+ZMU(1,3)*D(3) 72 : TEMP(2) = ZMU(2,1)*D(1)+ZMU(2,2)*D(2)+ZMU(2,3)*D(3) 73 : TEMP(3) = ZMU(3,1)*D(1)+ZMU(3,2)*D(2)+ZMU(3,3)*D(3) 74 : SEQE = TEMP(1)*D(1)+TEMP(2)*D(2)+TEMP(3)*D(3) 75 : SEQE = SEQE/((1.0D0-DD)*(1.0D0-DC)) 76 : SEQE = SQRT(SEQE) 77 : C--------------------------------------------------------------------| 78 : C******* CALCUL OF SIGV/K | 79 : C -------------------------------------------------------------------I 80 : RR = VAR(1) 81 : SK = (SEQE - RR-XMAT(5))/( XMAT(8) ) 82 : C--------------------------------------------------------------------| 83 : C******* CALCULATION OF EFFECTIVE INELASTIC STRAIN INCREMENT (P) | 84 : C -------------------------------------------------------------------I 85 : IF (SK.GT.0.0D0) THEN 86 : RPT= SK**XMAT(9) 87 : ELSE 88 : RPT= 0.0D0 89 : ENDIF 90 : C-----------------------------------------------------------------| 91 : C******* CALCULATION OF DR | 92 : C-----------------------------------------------------------------I 93 : DR= XMAT(6)*(XMAT(7)-VAR(1))*RPT 94 : VARPT(1)=DR 95 : C---------------------------------------------------------------------| 96 : C******* CALCULATION OF INELASTIC STRAIN INCREMENTS (EPSVPT) | 97 : C---------------------------------------------------------------------I 98 : C Calcul des deformations principales 99 : EPSD(1) = ZMU(1,1)*D(1)+ZMU(1,2)*D(2)+ZMU(1,3)*D(3) 100 : EPSD(1) = EPSD(1)*RPT/(SEQE*(1.0D0-DD)*(1.0D0-DC)) 101 : EPSD(2) = ZMU(2,1)*D(1)+ZMU(2,2)*D(2)+ZMU(2,3)*D(3) 102 : EPSD(2) = EPSD(2)*RPT/(SEQE*(1.0D0-DD)*(1.0D0-DC)) 103 : EPSD(3) = ZMU(3,1)*D(1)+ZMU(3,2)*D(2)+ZMU(3,3)*D(3) 104 : EPSD(3) = EPSD(3)*RPT/(SEQE*(1.0D0-DD)*(1.0D0-DC)) 105 : C 106 : IF (SEQE.EQ.0.0D0) THEN 107 : PPT = 0.0D0 108 : ELSE 109 : PPT = DETIER*(EPSD(1)*EPSD(1)+EPSD(2)*EPSD(2)+EPSD(3)*EPSD(3)) 110 : PPT = SQRT(PPT) 111 : ENDIF 112 : C 113 : VARPT(2) = PPT 114 : DO 11 I=1,NSTRS0 115 : EPS0(I)=0.0D0 116 : 11 CONTINUE 117 : C 118 : EPS0(1)=XP(1,1)*XP(1,1)*EPSD(1) 119 : EPS0(1)=EPS0(1)+XP(1,2)*XP(1,2)*EPSD(2) 120 : EPS0(1)=EPS0(1)+XP(1,3)*XP(1,3)*EPSD(3) 121 : EPS0(2)=XP(2,1)*XP(2,1)*EPSD(1) 122 : EPS0(2)=EPS0(2)+XP(2,2)*XP(2,2)*EPSD(2) 123 : EPS0(2)=EPS0(2)+XP(2,3)*XP(2,3)*EPSD(3) 124 : EPS0(3)=XP(3,1)*XP(3,1)*EPSD(1) 125 : EPS0(3)=EPS0(3)+XP(3,2)*XP(3,2)*EPSD(2) 126 : EPS0(3)=EPS0(3)+XP(3,3)*XP(3,3)*EPSD(3) 127 : EPS0(4)=XP(1,1)*XP(2,1)*EPSD(1) 128 : EPS0(4)=EPS0(4)+XP(1,2)*XP(2,2)*EPSD(2) 129 : EPS0(4)=EPS0(4)+XP(1,3)*XP(2,3)*EPSD(3) 130 : IF(NSTRS0.NE.6) GOTO 56 131 : EPS0(5)=XP(1,1)*XP(3,1)*EPSD(1) 132 : EPS0(5)=EPS0(5)+XP(1,2)*XP(3,2)*EPSD(2) 133 : EPS0(5)=EPS0(5)+XP(1,3)*XP(3,3)*EPSD(3) 134 : EPS0(6)=XP(2,1)*XP(3,1)*EPSD(1) 135 : EPS0(6)=EPS0(6)+XP(2,2)*XP(3,2)*EPSD(2) 136 : EPS0(6)=EPS0(6)+XP(2,3)*XP(3,3)*EPSD(3) 137 : 56 CONTINUE 138 : C 139 : DO 47 I=1,NSTRS0,1 140 : IF (RPT.EQ.0.0) THEN 141 : EPSVPT(I)=0.0D0 142 : ELSE 143 : EPSVPT(I)=EPS0(I) 144 : ENDIF 145 : 47 CONTINUE 146 : C-------------------------------------------------------- 147 : C******* CALCUL de D ductile AVEC SEUIL D ENDOMMAGEMENT 148 : C-------------------------------------------------------- 149 : YD = (D(1)*D(1)/(1.0D0-DC))+(D(2)*D(2))+(D(3)*D(3)) 150 : YD = YD-2.0D0*XMAT(2)*(D(1)*D(2))+(D(2)*D(3)+D(1)*D(3)) 151 : YD = YD/(2.0D0*XMAT(1)*(1.0D0-DD)*(1.0D0-DD)) 152 : C 153 : P = VAR(2) 154 : SE1 = XMAT(12) 155 : IF ((P.GT.SE1).AND.(DC.LT.0.98)) THEN 156 : DDD = ((YD/XMAT(10))**XMAT(11))*RPT 157 : ELSE 158 : DDD = 0.0D0 159 : ENDIF 160 : VARPT(3)= DDD 161 : C write(6,*)'DDD',DDD 162 : C----------------------------------------------------- 163 : C******* CALCUL de D creep AVEC SEUIL D ENDOMMAGEMENT 164 : C----------------------------------------------------- 165 : YC = (D(1)*D(1)/((1.0D0-DC)**2)) 166 : YC = YC/(2.0D0*XMAT(1)*(1.0D0-DD)) 167 : C 168 : SE2 = XMAT(15) 169 : IF ((P.GT.SE2).AND.(DC.LT.0.98)) THEN 170 : DDC = (YC/XMAT(13))**XMAT(14) 171 : ELSE 172 : DDC = 0.0D0 173 : ENDIF 174 : VARPT(4)= DDC 175 : C write(6,*)'DDC',DDC 176 : C 177 : RETURN 178 : END 179 :
© Cast3M 2003 - All rights reserved.
Disclaimer