1 : C INCRE6 SOURCE CHAT 05/01/13 00:35:27 5004 2 : SUBROUTINE INCRE6(SIG,VAR,EPSVPT,VARPT,XMAT,ALFA,NSTRS0, 3 : & NVARI,NCOMAT) 4 : IMPLICIT INTEGER(I-N) 5 : IMPLICIT REAL*8(A-H,O-Z) 6 : DIMENSION SIG(*),VAR(*),EPSVPT(*),VARPT(*),XMAT(*) 7 : DIMENSION SM(4),SB(4),XXM1(4),XXB1(4),XXM2(4),XXB2(4) 8 : DIMENSION GAMM(4),GABB(4),EPMM(4),EPBB(4),XXM(4),XXB(4) 9 : DIMENSION ANSM(4),ANSB(4),ANM(4),ANB(4) 10 : DIMENSION XK1M(4),XK1B(4),XK2M(4),XK2B(4) 11 : DIMENSION XEM(4),XEB(4),EPS0(8) 12 : PARAMETER (AMAX = 1.0E20) 13 : DETIER=2.0D0/3.0D0 14 : ROOT=SQRT(1.5D0) 15 : C-------------------------------------------------------------------| 16 : C STEP-1 SEPARATING MEMBRANE & BENDING PARTS | 17 : C-------------------------------------------------------------------I 18 : NSTRS = NSTRS0+2 19 : NS = NSTRS/2 20 : SM(1)=SIG(1) 21 : SM(2)=SIG(2) 22 : SM(3)=0.D0 23 : IF (NS.EQ.4) SM(4)=SIG(3) 24 : SB(1)=SIG(NS) 25 : SB(2)=SIG(NS+1) 26 : SB(3)=0.D0 27 : IF (NS.EQ.4) SB(4)=SIG(NS+2) 28 : DO 20 I=1,NS 29 : XXM1(I) = VAR( I) 30 : XXM2(I) = VAR(2*NS+I) 31 : XXB1(I) = VAR( NS+I) 32 : XXB2(I) = VAR(3*NS+I) 33 : GAMM(I) = VAR(4*NS+I) 34 : GABB(I) = VAR(5*NS+I) 35 : EPMM(I) = VAR(6*NS+I) 36 : EPBB(I) = VAR(7*NS+I) 37 : 20 CONTINUE 38 : C-------------------------------------------------------------------| 39 : C STEP-2: EVALUATION OF J2 SIGMA - X | 40 : C-------------------------------------------------------------------I 41 : XXM(1)= (2.0D0*SM(1)-SM(2))/3.0D0 -XXM1(1) -XXM2(1) 42 : XXM(2)= (2.0D0*SM(2)-SM(1))/3.0D0 -XXM1(2) -XXM2(2) 43 : XXM(3)= -1.D0*(SM(1)+SM(2))/3.0D0 -XXM1(3) -XXM2(3) 44 : XXM(4)= 0.D0 45 : IF(NS.EQ.4) XXM(4)= SM(4) -XXM1(4) -XXM2(4) 46 : XXB(1)= (2.0D0*SB(1)-SB(2))/3.0D0 -XXB1(1) -XXB2(1) 47 : XXB(2)= (2.0D0*SB(2)-SB(1))/3.0D0 -XXB1(2) -XXB2(2) 48 : XXB(3)= -1.D0*(SB(1)+SB(2))/3.0D0 -XXB1(3) -XXB2(3) 49 : XXB(4)= 0.D0 50 : IF(NS.EQ.4) XXB(4)= SB(4) -XXB1(4) -XXB2(4) 51 : * T1 = XXM(1)*XXM(1) + XXM(2)*XXM(2) + XXM(1)*XXM(2) 52 : * IF(NS.EQ.3) T1 = T1 + XXM(3)*XXM(3) 53 : T1 = PROCON(XXM,XXM,NS) 54 : AJM2 = 3.0D0*T1/2.D0 55 : * T1 = XXB(1)*XXB(1) + XXB(2)*XXB(2) + XXB(1)*XXB(2) 56 : * IF(NS.EQ.3) T1 = T1 + XXB(3)*XXB(3) 57 : T1 = PROCON(XXB,XXB,NS) 58 : AJB2 = 3.0D0*T1/2.D0 59 : AJ2 = SQRT( AJM2 + ALFA*ALFA*AJB2 ) 60 : C--------------------------------------------------------------------| 61 : C STEP-3: CALCUL OF SIGV/K | 62 : C -------------------------------------------------------------------I 63 : RR = VAR (4*NSTRS+2) 64 : RS = XMAT(10)*RR 65 : SK = (AJ2 - RS-XMAT(7))/( XMAT(8) + XMAT(9)*RR ) 66 : C--------------------------------------------------------------------| 67 : C STEP-4: CALCULATION OF EFFECTIVE INELASTIC STRAIN INCREMENT (P) | 68 : C -------------------------------------------------------------------I 69 : IF (SK.GT.0.0D0) THEN 70 : PPT1= XPUISS(SK,XMAT(6) ,AMAX) 71 : PPT2= XPUISS(SK,XMAT(6)+1,AMAX) 72 : ELSE 73 : PPT1= 0.0D0 74 : PPT2= 0.0D0 75 : ENDIF 76 : T1 = XMAT(5)*PPT2 77 : PPT2=1.0D20 78 : IF(DABS(T1).LT.40.0D0) PPT2= DEXP(T1) 79 : PPT = PPT1 * PPT2 80 : VARPT(4*NSTRS+1) = PPT 81 : C---------------------------------------------------------------------| 82 : C******* CALCULATION OF INELASTIC STRAIN INCREMENTS (EPSVPT) | 83 : C---------------------------------------------------------------------I 84 : IF(AJ2.EQ.0.0D0) AJ2=1.0D-20 85 : T1 = 1.5D0*PPT/AJ2 86 : EPS0(1) = T1*XXM(1) 87 : EPS0(2) = T1*XXM(2) 88 : EPS0(3) = T1*XXM(3) 89 : IF(NS.EQ.4) EPS0(4)= T1*XXM(4) 90 : EPS0(NS+1) = T1*XXB(1) *ALFA *ALFA 91 : EPS0(NS+2) = T1*XXB(2) *ALFA *ALFA 92 : EPS0(NS+3) = T1*XXB(3) *ALFA *ALFA 93 : IF(NS.EQ.4) EPS0(NS+4)= T1*XXB(4) *ALFA*ALFA 94 : DO 30 I=1,NSTRS 95 : 30 VARPT(3*NSTRS+I)=EPS0(I) 96 : C 97 : C----------------------------------------------------------------| 98 : C ******* CALCULATION OF PI(P) | 99 : C----------------------------------------------------------------I 100 : P = VAR(4*NSTRS+1) 101 : T1 = XMAT(16)*P 102 : T2 = DEXP(-T1) 103 : PIP= XMAT(13)+(1.0D0-XMAT(13)) * T2 104 : C----------------------------------------------------------------| 105 : C******* CALCULATION OF XII | 106 : C----------------------------------------------------------------I 107 : C CORRECTIONS MLR 26/2/93 COMME DANS INCRE2 108 : C 109 : * T1 = 3.0D0*(XXM1(1)*XXM1(1)+XXM1(2)*XXM1(2)+XXM1(1)*XXM1(2)) 110 : * IF(NS.EQ.3) T1 = T1 + 3.0D0 * XXM1(3) * XXM1(3) 111 : T1M = PROCON(XXM1,XXM1,NS) 112 : T1M = 3.D0*T1M/2.D0 113 : * T2 = 3.0D0*(XXB1(1)*XXB1(1)+XXB1(2)*XXB1(2)+XXB1(1)*XXB1(2)) 114 : * IF(NS.EQ.3) T2 = T2 + 3.0D0 * XXB1(3) * XXB1(3) 115 : T1B = PROCON(XXB1,XXB1,NS) 116 : T1B = 3.D0*T1B/2.D0 117 : X1II = 0.D0 118 : IF(XMAT(17).NE.0.D0) 119 : . X1II = SQRT(T1M+ALFA*ALFA*T1B)/XMAT(17) 120 : * T1 = 3.0D0*(XXM2(1)*XXM2(1)+XXM2(2)*XXM2(2)+XXM2(1)*XXM2(2)) 121 : * IF(NS.EQ.3) T1 = T1 + 3.0D0 * XXM2(3) * XXM2(3) 122 : T2M = PROCON(XXM2,XXM2,NS) 123 : T2M = 3.D0*T2M/2.D0 124 : * T2 = 3.0D0*(XXB2(1)*XXB2(1)+XXB2(2)*XXB2(2)+XXB2(1)*XXB2(2)) 125 : * IF(NS.EQ.3) T2 = T2 + 3.0D0 * XXB2(3) * XXB2(3) 126 : T2B = PROCON(XXB2,XXB2,NS) 127 : T2B = 3.D0*T2B/2.D0 128 : X2II = 0.D0 129 : IF(XMAT(18).NE.0.D0) 130 : . X2II = SQRT(T2M+ALFA*ALFA*T2B)/XMAT(18) 131 : COX1 = XPUISS(X1II,XMAT(19)-1.0D0,AMAX) 132 : COX2 = XPUISS(X2II,XMAT(20)-1.0D0,AMAX) 133 : C-----------------------------------------------------------------| 134 : C***** CALCULATION OF THE SPECIFIC TERMS OF THE OHNO'S MODEL *****| 135 : C-----------------------------------------------------------------| 136 : XLIM1 = XMAT(12)/PIP 137 : XLIM2 = XMAT(15)/PIP 138 : * T1 = 3.0D0*(XXM1(1)*XXM1(1)+XXM1(2)*XXM1(2)+XXM1(1)*XXM1(2)) 139 : * IF(NS.EQ.3) T1 = T1 + 3.0D0 * XXM1(3) * XXM1(3) 140 : * T2 = 3.0D0*(XXB1(1)*XXB1(1)+XXB1(2)*XXB1(2)+XXB1(1)*XXB1(2)) 141 : * IF(NS.EQ.3) T2 = T2 + 3.0D0 * XXB1(3) * XXB1(3) 142 : X1II = SQRT(T1M+ALFA*ALFA*T1B) 143 : * T1 = 3.0D0*(XXM2(1)*XXM2(1)+XXM2(2)*XXM2(2)+XXM2(1)*XXM2(2)) 144 : * IF(NS.EQ.3) T1 = T1 + 3.0D0 * XXM2(3) * XXM2(3) 145 : * T2 = 3.0D0*(XXB2(1)*XXB2(1)+XXB2(2)*XXB2(2)+XXB2(1)*XXB2(2)) 146 : * IF(NS.EQ.3) T2 = T2 + 3.0D0 * XXB2(3) * XXB2(3) 147 : X2II = SQRT(T2M+ALFA*ALFA*T2B) 148 : RAPP1 = (X1II/XLIM1)**XMAT(27) 149 : RAPP2 = (X2II/XLIM2)**XMAT(28) 150 : DO 11 I=1,NS 151 : IF (X1II.GT.1.D-20) THEN 152 : XK1M(I)=XXM1(I)/X1II 153 : ELSE 154 : XK1M(I)=0.D0 155 : ENDIF 156 : IF (X2II.GT.1.D-20) THEN 157 : XK2M(I)=XXM2(I)/X2II 158 : ELSE 159 : XK2M(I)=0.D0 160 : ENDIF 161 : XEM(I)=EPS0(I) 162 : 11 CONTINUE 163 : DO 21 I=1,NS 164 : IF (X1II.GT.1.D-20) THEN 165 : XK1B(I)=XXB1(I)/X1II 166 : ELSE 167 : XK1B(I)=0.D0 168 : ENDIF 169 : IF (X2II.GT.1.D-20) THEN 170 : XK2B(I)=XXB2(I)/X2II 171 : ELSE 172 : XK2B(I)=0.D0 173 : ENDIF 174 : XEB(I)=EPS0(NS+I) 175 : 21 CONTINUE 176 : PRODM1 = PROCON(XEM,XK1M,NS) 177 : * PRODM1 = XEM(1)*XK1M(1)+XEM(2)*XK1M(2) 178 : * IF (NS.EQ.3) PRODM1 = PRODM1+2.D0*XEM(3)*XK1M(3) 179 : PRODB1 = PROCON(XEB,XK1B,NS) 180 : * PRODB1 = XEB(1)*XK1B(1)+XEB(2)*XK1B(2) 181 : * IF (NS.EQ.3) PRODB1 = PRODB1+2.D0*XEB(3)*XK1B(3) 182 : PRODM2 = PROCON(XEM,XK2M,NS) 183 : * PRODM2 = XEM(1)*XK2M(1)+XEM(2)*XK2M(2) 184 : * IF (NS.EQ.3) PRODM2 = PRODM2+2.D0*XEM(3)*XK2M(3) 185 : PRODB2 = PROCON(XEB,XK2B,NS) 186 : * PRODB2 = XEB(1)*XK2B(1)+XEB(2)*XK2B(2) 187 : * IF (NS.EQ.3) PRODB2 = PRODB2+2.D0*XEB(3)*XK2B(3) 188 : PROD1 = PRODM1+PRODB1 189 : IF (PROD1.LT.0.D0) PROD1=0.D0 190 : PROD2 = PRODM2+PRODB2 191 : IF (PROD2.LT.0.D0) PROD2=0.D0 192 : PPT10 = RAPP1*PROD1 193 : PPT20 = RAPP2*PROD2 194 : C ----------------------------------------------------------------| 195 : C******* CALCULATION OF DX1 AND DX2 | 196 : C-----------------------------------------------------------------I 197 : DO 73 I=1,NSTRS 198 : VARPT( I)=0.0D0 199 : 73 VARPT(NSTRS+I)=0.0D0 200 : DO 74 I=1,NSTRS 201 : T1 = DETIER * XMAT(11) * XMAT(12) * EPS0(I) 202 : T2 = XMAT(11) * PIP * VAR(I) * PPT10 203 : T3 = COX1 * VAR(I) 204 : VARPT(I)= T1 - T2 -T3 205 : T1 = DETIER * XMAT(14) * XMAT(15) * EPS0(I) 206 : T2 = XMAT(14) * PIP * VAR(NSTRS+I) * PPT20 207 : T3 = COX2 * VAR(NSTRS+I) 208 : 74 VARPT(NSTRS+I)= T1 - T2 - T3 209 : C-----------------------------------------------------------------| 210 : C******* CALCULATION OF DR | 211 : C-----------------------------------------------------------------I 212 : T1 = 1.0D0-(1.0D0-VAR(4*NSTRS+3)/XMAT(25))**2 213 : QR = VAR(4*NSTRS+3)-XMAT(26)*T1 214 : CO1= QR-VAR(4*NSTRS+2) 215 : CO=DABS(CO1) 216 : CO = XPUISS(CO,XMAT(22)-1.0D0,AMAX) 217 : T1= XMAT(16)*(VAR(4*NSTRS+3)-VAR(4*NSTRS+2))*PPT 218 : T2= XMAT(21)*CO*CO1 219 : DR= T1 + T2 220 : VARPT(4*NSTRS+2)=DR 221 : C----------------------------------------------------------------| 222 : C******* CALCULATION OF PROD (N X N*) | 223 : C----------------------------------------------------------------I 224 : TM1 = EPMM(1)-GAMM(1) 225 : TM2 = EPMM(2)-GAMM(2) 226 : TM3 = EPMM(3)-GAMM(3) 227 : TM4 = 0.0D0 228 : IF(NS.EQ.4) TM4 = EPMM(4)-GAMM(4) 229 : AIM2 = TM1*TM1+TM2*TM2+TM3*TM3+2.D0*TM4*TM4 230 : AIM2 = AIM2*3.D0/2.D0 231 : TB1 = EPBB(1)-GABB(1) 232 : TB2 = EPBB(2)-GABB(2) 233 : TB3 = EPBB(3)-GABB(3) 234 : TB4 = 0.0D0 235 : IF(NS.EQ.4) TB4 = EPBB(4)-GABB(4) 236 : AIB2 = TB1*TB1+TB2*TB2+TB3*TB3+2.D0*TB4*TB4 237 : AIB2 = AIB2*3.D0/2.D0 238 : AI2 = SQRT(AIM2+AIB2/ALFA/ALFA) 239 : IF (AI2.LT.1.D-10) THEN 240 : IF (AI2.EQ.0.D0) AI2=1.D-10 241 : AI20=AI2*1.D20 242 : TM1=TM1*1.D20 243 : TM2=TM2*1.D20 244 : TM3=TM3*1.D20 245 : TM4=TM4*1.D20 246 : TB1=TB1*1.D20 247 : TB2=TB2*1.D20 248 : TB3=TB3*1.D20 249 : TB4=TB4*1.D20 250 : ELSE 251 : AI20=AI2 252 : ENDIF 253 : ANSM(1) = ROOT*TM1/AI20 254 : ANSM(2) = ROOT*TM2/AI20 255 : ANSM(3) = ROOT*TM3/AI20 256 : ANSM(4) = 0.0D0 257 : IF(NS.EQ.4) ANSM(4) = ROOT*TM4/AI20 258 : ANSB(1) = ROOT*TB1/AI20 259 : ANSB(2) = ROOT*TB2/AI20 260 : ANSB(3) = ROOT*TB3/AI20 261 : ANSB(4) = 0.0D0 262 : IF(NS.EQ.4) ANSB(4) = ROOT*TB4/AI20 263 : ANM (1) = ROOT*XXM(1)/AJ2 264 : ANM (2) = ROOT*XXM(2)/AJ2 265 : ANM (3) = ROOT*XXM(3)/AJ2 266 : ANM (4) = 0.0D0 267 : IF(NS.EQ.4) ANM (4) = ROOT*XXM(4)/AJ2 268 : ANB (1) = ROOT*XXB(1)/AJ2 269 : ANB (2) = ROOT*XXB(2)/AJ2 270 : ANB (3) = ROOT*XXB(3)/AJ2 271 : ANB (4) = 0.0D0 272 : IF(NS.EQ.4) ANB (4) = ROOT*XXB(4)/AJ2 273 : 76 CONTINUE 274 : FF = DETIER*AI2-VAR(4*NSTRS+4) 275 : HF=0.0D0 276 : IF(FF.GT.0.0D0) HF=1.0D0 277 : ANAN1=PROCON(ANM,ANSM,NS) 278 : * ANAN1 = ANM(1)*ANSM(1)+ANM(2)*ANSM(2) 279 : * ANAN1 = ANAN1+2.D0*ANM(3)*ANSM(3) 280 : ANANSM=0.0D0 281 : IF(ANAN1.GT.0.0D0) ANANSM=ANAN1 282 : ANAN1=PROCON(ANB,ANSB,NS) 283 : * ANAN1 = ANB(1)*ANSB(1)+ANB(2)*ANSB(2) 284 : * ANAN1 = ANAN1+2.D0*ANB(3)*ANSB(3) 285 : ANANSB=0.0D0 286 : IF(ANAN1.GT.0.0D0) ANANSB=ANAN1 287 : C----------------------------------------------------------------| 288 : C******* CALCULATION OF DQ | 289 : C----------------------------------------------------------------I 290 : ANANS=ANANSM+ANANSB/ALFA/ALFA 291 : DQ = XMAT(23)*HF*ANANS*PPT 292 : VARPT(4*NSTRS+4)=DQ 293 : C----------------------------------------------------------------| 294 : C******* CALCULATION OF DGETA | 295 : C----------------------------------------------------------------I 296 : T1 = ROOT*(1.0D0-XMAT(23))*HF*ANANS*PPT 297 : DO 78 I=1,NS 298 : VARPT(2*NSTRS+I) = T1*ANSM(I) 299 : VARPT(2*NSTRS+NS+I)= T1*ANSB(I) 300 : 78 CONTINUE 301 : C----------------------------------------------------------------| 302 : C******* CALCULATION OF DQQ | 303 : C----------------------------------------------------------------I 304 : VARPT(4*NSTRS+3)=2.0D0*XMAT(24)*(XMAT(25)-VAR(4*NSTRS+3))*DQ 305 : C 306 : EPSVPT(1)=EPS0(1) 307 : EPSVPT(2)=EPS0(2) 308 : IF (NS.EQ.4) EPSVPT(3)=EPS0(4) 309 : EPSVPT(NS)=EPS0(NS+1) 310 : EPSVPT(NS+1)=EPS0(NS+2) 311 : IF (NS.EQ.4) EPSVPT(NS+2)=EPS0(NS+4) 312 : C 313 : RETURN 314 : END 315 : 316 : 317 : 318 :
© Cast3M 2003 - All rights reserved.
Disclaimer