1 : C T1DP SOURCE CHAT 05/01/13 03:29:29 5004 2 : SUBROUTINE T1DP(SIG,DSIG,YOUN,ANU,RT1,RT2,H1,RDP,ADP,HDP, 3 : 1 XLAM1,XLAM2,XLAM3,IDAM,ANG,EPPLDP,ALPHA,KERRE) 4 : IMPLICIT INTEGER(I-N) 5 : IMPLICIT REAL*8(A-H,O-Z) 6 : -INC CCOPTIO 7 : DIMENSION SIG(3),DSIG(3),SGG(3),DSGG(3),SFG(3),DSFG(3), 8 : 1 A(3),B(3),IDAM(3),EPC(3),SSI(3),SI(3),EPC1(3), 9 : 1 EPDP(3),EC(3),EPPLDP(3),SG(3),AT(4),BT(2) 10 : ITER=0 11 : SIREF=1.E-6*YOUN 12 : Y=YOUN/(1.D0-ANU*ANU) 13 : DO 10 ITYP=1,3 14 : EPDP(ITYP)=0.D0 15 : 10 EPPLDP(ITYP)=0.D0 16 : DL3=0.D0 17 : XT=0.D0 18 : IF(IIMPI.EQ.9) WRITE(IOIMP,9999) 19 : 9999 FORMAT(1X,'T1DP COUPLAGE TRACTION 1 DRUCKER',/) 20 : 9998 FORMAT(1X,'T1DP ',I4,'ITERATIONS INTERNES',/) 21 : C 22 : C------------------------------------------------ 23 : C COUPLAGE TRACTION DRUCKER 24 : C------------------------------------------------ 25 : C---------------------------------------------- 26 : C DEFINITION D UNE CONTRAINTE DE REFERENCE 27 : C SI RMAX DU DRUCKER INFERIEURE A CETTE VALEUR 28 : C ALORS SIGMA=0 29 : C---------------------------------------------- 30 : C 31 : RMAX=MAX((RDP/1.73),(RDP/(1.D0-2.D0*ADP))) 32 : IF(SIREF.GT.RMAX) THEN 33 : DO 20 ITYP=1,3 34 : 20 SIG(ITYP)=SIG(ITYP)+DSIG(ITYP) 35 : RDP=0.D0 36 : ADP=0.D0 37 : HDP=0.D0 37 : CALL CPHOIN(SIG,EPPLDP,YOUN,ANU) 37 : CALL NORME(SIG,DL3) 38 : DL3=DL3/YOUN*10.D0 39 : XLAM3=XLAM3+DL3 40 : DO 30 ITYP=1,3 41 : IDAM(ITYP)=0 42 : DSIG(ITYP)=0.D0 43 : 30 SIG(ITYP)=0.D0 44 : IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER 45 : RETURN 46 : ENDIF 47 : C----------------------------------------------- 48 : C ON SEPLACE DANS LE REPERE DE FISSURATION 49 : C----------------------------------------------- 50 : C 50 : CALL CHREP(SIG,SFG,ANG) 50 : CALL CHREP(DSIG,DSFG,ANG) 51 : C 52 : C----------------------------------------------------- 53 : C ESTIMATION DU PAS D'INCREMENT DE CONTRAINTES 54 : C----------------------------------------------------- 55 : C 55 : CALL SIJ(SFG,SI,SIETOI) 55 : CALL EPSDP(SI,SIETOI,ADP,EPC) 55 : CALL CPHOMO(EPC,SSI,YOUN,ANU,ALPHA) 55 : CALL SCAL2(EPC,SSI,VAL) 56 : VAL=VAL-HDP 57 : IF(VAL.LT.0.D0) THEN 58 : KERRE=459 59 : RETURN 60 : ENDIF 60 : CALL SCAL(SSI,DSFG,VAL) 60 : CALL NORME(DSFG,VA1) 60 : CALL NORME(SSI,VA2) 61 : IF(VA1.EQ.0.D0) THEN 62 : IDAM(1)=0 63 : IDAM(2)=0 64 : IDAM(3)=0 65 : RETURN 66 : ENDIF 67 : RMIN=MIN((RDP/1.73),(RDP/(1.D0+2.D0*ADP))) 68 : X=VAL/VA1/VA2 69 : IF(X.GE.1.D0) THEN 70 : X=1.D0 71 : GOTO 41 72 : ENDIF 73 : X=1.D0/SQRT(1.0001D0-X*X)*RMIN/VA1/8.D0 74 : IF(X.GE.1.D0) X=1.D0 75 : C 76 : C---------------------------------------------------------- 77 : C ON ECOULE ET ON REGARDE LA VARIATION DE NORMALE 78 : C---------------------------------------------------------- 79 : C 80 : C 81 : C--------------------------------------- 82 : C RESOLUTION DU SYSTEME EN DL3 83 : C--------------------------------------- 84 : 41 X=X/2.D0 85 : EC(1)=1.D0 86 : EC(2)=0.D0 87 : EC(3)=0.D0 87 : CALL CPHOOB(EC,SG,YOUN,ANU) 88 : 100 X=X*2.D0 88 : CALL SIJ(SFG,SI,SIETOI) 88 : CALL EPSDP(SI,SIETOI,ADP,EPC) 88 : CALL CPHOMO(EPC,SSI,YOUN,ANU,ALPHA) 88 : CALL SCAL2(EPC,SSI,VAL) 89 : VAL=VAL-HDP 90 : IF(VAL.LT.0.D0) THEN 91 : KERRE=459 92 : RETURN 93 : ENDIF 94 : IF((XT+X).GT.1.D0) X=1.D0-XT 95 : 110 RMAX=MAX((RDP/1.73),(RDP/(1.D0-2.D0*ADP))) 96 : C 97 : C -------------------------------------------------- 98 : C CAS OU LE RAYON MAXIMAL DU DRUCKER-PRAGER 99 : C INFERIEUR A LA CONTRAINTE DE REFERENCE 100 : C -------------------------------------------------- 101 : C 102 : IF(SIREF.LT.RMAX) GOTO 131 103 : DO 120 ITYP=1,3 104 : 120 SFG(ITYP)=SFG(ITYP)+(1.D0-XT)*DSFG(ITYP) 105 : RDP=0.D0 106 : ADP=0.D0 107 : HDP=0.D0 107 : CALL NORME(SFG,DL3) 108 : DL3=DL3/YOUN*10.D0 109 : XLAM3=XLAM3+DL3 109 : CALL CPHOIN(SFG,EPPLDP,YOUN,ANU) 110 : DO 130 ITYP=1,3 111 : IDAM(ITYP)=0 112 : EPDP(ITYP)=EPDP(ITYP)+EPPLDP(ITYP) 113 : DSIG(ITYP)=0.D0 114 : 130 SIG(ITYP)=0.D0 114 : CALL CHREP(EPDP,EPPLDP,-ANG) 115 : IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER 116 : RETURN 117 : 131 ITER=ITER+1 118 : IF(ITER.GT.200) THEN 119 : KERRE=460 120 : RETURN 121 : ENDIF 122 : C 123 : C--------------------------------- 124 : C ON REALISE LE COUPLAGE 125 : C--------------------------------- 126 : C 127 : E=-(SFG(1)+X*DSFG(1)-RT1)/(SG(1)-H1) 128 : F=SSI(1)/(SG(1)-H1) 129 : DO 140 ITYP=1,3 130 : DSGG(ITYP)=X*DSFG(ITYP) 131 : A(ITYP)=SFG(ITYP)+DSGG(ITYP)+E*SG(ITYP) 132 : 140 B(ITYP)=-SSI(ITYP)+F*SG(ITYP) 132 : CALL DLAMDP(A,B,DL3,RDP,ADP,HDP,ITEST) 133 : IF(ITEST.EQ.1) THEN 134 : X=X/2.D0 135 : GOTO 131 136 : ENDIF 137 : DL1=-E-F*DL3 138 : DO 145 ITYP=1,3 139 : 145 SGG(ITYP)=SFG(ITYP)+DSGG(ITYP)-DL1*SG(ITYP)-DL3*SSI(ITYP) 140 : C 141 : C---------------------------------------------- 142 : C SI UN DES DLAMDA < 0 FAUX COUPLAGE 143 : C---------------------------------------------- 144 : C 145 : IF(DL1.LT.-1.D-10) IDAM(1)=0 146 : IF(DL3.LT.-1.D-10) IDAM(3)=0 147 : IF(IDAM(1).EQ.0.OR.IDAM(3).EQ.0) THEN 148 : DL1=0.D0 149 : DL3=0.D0 150 : X=0.D0 151 : GOTO 3000 152 : ENDIF 153 : C 154 : C----------------------------------------------------- 155 : C ON VERIFIE QUE LA NORMALE SUIVANT LE DRUCKER 156 : C PRAGER TOURNE PEU 157 : C----------------------------------------------------- 158 : C 158 : CALL SIJ(SGG,SI,SIETOI) 158 : CALL EPSDP(SI,SIETOI,ADP,EPC1) 158 : CALL SCAL(EPC,EPC1,VAL) 158 : CALL NORME(EPC,VA1) 158 : CALL NORME(EPC1,VA2) 159 : CO=VAL/VA1/VA2 160 : IF(CO.LT.0.99) THEN 161 : X=X/2.D0 162 : GOTO 131 163 : ENDIF 164 : IF((RT1-H1*DL1).GE.-1.E-10) GOTO 200 165 : C 166 : C------------------------------------------------ 167 : C CAS OU LA LIMITE EN TRACTION CHUTE EN 168 : C DESSOUS DE 0 169 : C------------------------------------------------ 170 : C 171 : DL1=RT1/H1 172 : IF(DSFG(1).EQ.0.D0) GOTO 160 173 : E=(DL1*SG(1)-SFG(1))/DSFG(1) 174 : F=SSI(1)/DSFG(1) 175 : DO 150 ITYP=1,3 176 : A(ITYP)=SFG(ITYP)-DL1*SG(ITYP)+E*DSFG(ITYP) 177 : 150 B(ITYP)=F*DSFG(ITYP)-SSI(ITYP) 177 : CALL DLAMD(A,B,DL31,DL32,RDP,ADP,HDP,ITEST) 178 : DO 151 ITYP=1,3 179 : SGG(ITYP)=A(ITYP)+DL31*B(ITYP) 180 : 151 DSGG(ITYP)=A(ITYP)+DL32*B(ITYP) 180 : CALL SIJ(SGG,SI,SIETOI) 180 : CALL EPSDP(SI,SIETOI,ADP,EPC1) 180 : CALL SCAL(EPC,EPC1,VAL) 180 : CALL NORME(EPC,VA1) 180 : CALL NORME(EPC1,VA2) 181 : CO1=VAL/VA1/VA2 181 : CALL SIJ(DSGG,SI,SIETOI) 181 : CALL EPSDP(SI,SIETOI,ADP,EPC1) 181 : CALL SCAL(EPC,EPC1,VAL) 181 : CALL NORME(EPC,VA1) 181 : CALL NORME(EPC1,VA2) 182 : CO2=VAL/VA1/VA2 183 : XIN1=E+DL31*F 184 : XIN2=E+DL32*F 185 : IF(XIN1.GT.-1.D-10.AND.XIN1.LE.X.AND.CO1.GT.0.9) THEN 186 : DL3=DL31 187 : X=XIN1 188 : GOTO 200 189 : ENDIF 190 : IF(XIN2.GT.-1.D-10.AND.XIN2.LE.X.AND.CO2.GT.0.9) THEN 191 : DL3=DL32 192 : X=XIN2 193 : GOTO 200 194 : ENDIF 195 : IF(IIMPI.EQ.9) WRITE(IOIMP,10101) 196 : 10101 FORMAT(1X,'ERREUR DANS T1DP DEPA ',/) 197 : DL3=DL31 198 : X=XIN1 199 : GOTO 200 200 : 160 DL3=(SFG(1)-DL1*SG(1))/SSI(1) 201 : RD=RDP-DL3*HDP 202 : IF(RD.LT.0.D0) THEN 203 : RDP=0.D0 204 : GOTO 110 205 : ENDIF 206 : DO 170 ITYP=1,3 207 : A(ITYP)=SFG(ITYP)-DL1*SG(ITYP)-DL3*SSI(ITYP) 208 : 170 B(ITYP)=DSFG(ITYP) 208 : CALL XDP(A,B,RD,ADP,X,ITEST) 209 : C 210 : C----------------------------------------------- 211 : C ON REGARDE SI ON ENDOMMAGE PAS LES AUTRES 212 : C CRITERES 213 : C----------------------------------------------- 214 : C------------------------------------- 215 : C ON ENDOMMAGE LA TRACTION 2 216 : C------------------------------------- 217 : C 218 : 200 CONTINUE 219 : DO 205 ITYP=1,3 220 : DSGG(ITYP)=X*DSFG(ITYP) 221 : 205 SGG(ITYP)=SFG(ITYP)+DSGG(ITYP)-DL1*SG(ITYP)-DL3*SSI(ITYP) 221 : CALL CTRAF(SGG(2),RT2,VCTR) 222 : IF(VCTR.LE.0.D0) GOTO 1000 223 : IDAM(1)=1 224 : IDAM(2)=1 225 : IDAM(3)=1 226 : AT(1)=-DSFG(1) 227 : AT(2)=Y-H1 228 : AT(4)=Y*ANU 229 : AT(3)=-DSFG(2) 230 : BT(1)=SFG(1)-RT1 231 : BT(2)=SFG(2)-RT2 231 : CALL SYLIN2(AT,BT,X1,DL11) 232 : BT(1)=-SSI(1) 233 : BT(2)=-SSI(2) 233 : CALL SYLIN2(AT,BT,X2,DL12) 234 : DO 230 ITYP=1,3 235 : A(ITYP)=SFG(ITYP)-DL11*SG(ITYP)+X1*DSFG(ITYP) 236 : 230 B(ITYP)=-DL12*SG(ITYP)+X2*DSFG(ITYP)-SSI(ITYP) 236 : CALL DLAMD(A,B,DL31,DL32,RDP,ADP,HDP,ITEST) 237 : DO 231 ITYP=1,3 238 : SGG(ITYP)=A(ITYP)+DL31*B(ITYP) 239 : 231 DSGG(ITYP)=A(ITYP)+DL32*B(ITYP) 239 : CALL SIJ(SGG,SI,SIETOI) 239 : CALL EPSDP(SI,SIETOI,ADP,EPC1) 239 : CALL SCAL(EPC,EPC1,VAL) 239 : CALL NORME(EPC,VA1) 239 : CALL NORME(EPC1,VA2) 240 : CO1=VAL/VA1/VA2 240 : CALL SIJ(DSGG,SI,SIETOI) 240 : CALL EPSDP(SI,SIETOI,ADP,EPC1) 240 : CALL SCAL(EPC,EPC1,VAL) 240 : CALL NORME(EPC,VA1) 240 : CALL NORME(EPC1,VA2) 241 : CO2=VAL/VA1/VA2 242 : DLIN1=DL11+DL31*DL12 243 : DLIN2=DL11+DL32*DL12 244 : XIN1=X1+DL31*X2 245 : XIN2=X1+DL32*X2 246 : IF(XIN1.GT.-1.D-10.AND.DLIN1.GT.-1.E-10.AND.XIN1.LE.X 247 : 1 .AND.CO1.GT.0.9) THEN 248 : DL3=DL31 249 : X=XIN1 250 : DL1=DLIN1 251 : GOTO 3000 252 : ENDIF 253 : IF(XIN2.GT.-1.D-10.AND.DLIN2.GT.-1.E-10.AND.XIN2.LE.X 254 : 1 .AND.CO2.GT.0.9) THEN 255 : DL3=DL32 256 : X=XIN2 257 : DL1=DLIN2 258 : GOTO 3000 259 : ENDIF 260 : IF(IIMPI.EQ.9) WRITE(IOIMP,20202) 261 : 20202 FORMAT(1X,'ERREUR DANS T1DP TRAC ',/) 262 : DL3=DL31 263 : DL1=DLIN1 264 : X=XIN1 265 : GOTO 3000 266 : C 267 : C------------------------------------- 268 : C ON ENDOMMAGE LA COMPRESSION 2 269 : C------------------------------------- 270 : C 270 : 1000 CALL CCOAF(SGG(2),XLAM2,VCCO2) 271 : IF(VCCO2.LE.0.D0) GOTO 2000 272 : IDAM(1)=1 273 : IDAM(3)=1 274 : IDAM(2)=-1 275 : AT(1)=-DSFG(1) 276 : AT(2)=Y-H1 277 : AT(4)=Y*ANU 278 : AT(3)=-DSFG(2) 279 : BT(1)=SFG(1)-RT1 280 : BT(2)=SFG(2) 280 : CALL SYLIN2(AT,BT,X1,DL11) 281 : BT(1)=-SSI(1) 282 : BT(2)=-SSI(2) 282 : CALL SYLIN2(AT,BT,X2,DL12) 283 : DO 1010 ITYP=1,3 284 : A(ITYP)=SFG(ITYP)-DL11*SG(ITYP)+X1*DSFG(ITYP) 285 : 1010 B(ITYP)=-DL12*SG(ITYP)+X2*DSFG(ITYP)-SSI(ITYP) 285 : CALL DLAMD(A,B,DL31,DL32,RDP,ADP,HDP,ITEST) 286 : DO 1011 ITYP=1,3 287 : SGG(ITYP)=A(ITYP)+DL31*B(ITYP) 288 : 1011 DSGG(ITYP)=A(ITYP)+DL32*B(ITYP) 288 : CALL SIJ(SGG,SI,SIETOI) 288 : CALL EPSDP(SI,SIETOI,ADP,EPC1) 288 : CALL SCAL(EPC,EPC1,VAL) 288 : CALL NORME(EPC,VA1) 288 : CALL NORME(EPC1,VA2) 289 : CO1=VAL/VA1/VA2 289 : CALL SIJ(DSGG,SI,SIETOI) 289 : CALL EPSDP(SI,SIETOI,ADP,EPC1) 289 : CALL SCAL(EPC,EPC1,VAL) 289 : CALL NORME(EPC,VA1) 289 : CALL NORME(EPC1,VA2) 290 : CO2=VAL/VA1/VA2 291 : DLIN1=DL11+DL31*DL12 292 : DLIN2=DL11+DL32*DL12 293 : XIN1=X1+DL31*X2 294 : XIN2=X1+DL32*X2 295 : IF(XIN1.GT.-1.D-10.AND.DLIN1.GT.-1.E-10.AND.XIN1.LE.X 296 : 1 .AND.CO1.GT.0.9) THEN 297 : DL3=DL31 298 : X=XIN1 299 : DL1=DLIN1 300 : GOTO 3000 301 : ENDIF 302 : IF(XIN2.GT.-1.D-10.AND.DLIN2.GT.-1.E-10.AND.XIN2.LE.X 303 : 1 .AND.CO2.GT.0.9) THEN 304 : DL3=DL32 305 : X=XIN2 306 : DL1=DLIN2 307 : GOTO 3000 308 : ENDIF 309 : IF(IIMPI.EQ.9) WRITE(IOIMP,30303) 310 : 30303 FORMAT(1X,'ERREUR DANS T1DP COMP ',/) 311 : DL3=DL31 312 : DL1=DLIN1 313 : X=XIN1 314 : GOTO 3000 315 : 2000 CONTINUE 316 : C 317 : C------------------------------------------------ 318 : C CAS OU AUCUN AUTRE CRITERE EST ENDOMMAGE 319 : C------------------------------------------------ 320 : C----------------------------------- 321 : C REMISE A JOUR DES VARIABLES 322 : C PUIS NOUVEAU PAS SI TOUT 323 : C L INCREMENT N EST PAS ECOULE 324 : C----------------------------------- 325 : C 326 : RDP=RDP-DL3*HDP 327 : RMAX=MAX((RDP/1.73),(RDP/(1.D0-2.D0*ADP))) 328 : IF(RMAX.LT.SIREF) THEN 329 : RDP=0.D0 330 : GOTO 110 331 : ENDIF 332 : RT1=RT1-H1*DL1 333 : IF(RT1.LE.1.E-10) THEN 334 : H1=0.D0 335 : RT1=0.D0 336 : ENDIF 337 : XT=XT+X 338 : DO 2010 ITYP=1,3 339 : SFG(ITYP)=SFG(ITYP)+X*DSFG(ITYP)-DL1*SG(ITYP)-DL3*SSI(ITYP) 340 : 2010 EPDP(ITYP)=DL3*EPC(ITYP)+EPDP(ITYP) 341 : XLAM1=XLAM1+DL1 342 : XLAM3=XLAM3+DL3 343 : C=1.D0-1.D-10 344 : IF(XT.LT.C) GOTO 100 345 : DO 2020 ITYP=1,3 346 : IDAM(ITYP)=0 347 : 2020 DSIG(ITYP)=0.D0 347 : CALL CHREP(SFG,SIG,-ANG) 347 : CALL CHREP(EPDP,EPPLDP,-ANG) 348 : IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER 349 : RETURN 350 : 3000 RDP=RDP-DL3*HDP 351 : RMAX=MAX((RDP/1.73),(RDP/(1.D0-2.D0*ADP))) 352 : IF(RMAX.LT.SIREF) THEN 353 : RDP=0.D0 354 : GOTO 110 355 : ENDIF 356 : RT1=RT1-H1*DL1 357 : IF(RT1.LE.1.E-10) THEN 358 : H1=0.D0 359 : RT1=0.D0 360 : ENDIF 361 : XT=XT+X 362 : DO 3010 ITYP=1,3 363 : SFG(ITYP)=SFG(ITYP)+X*DSFG(ITYP)-DL1*SG(ITYP)-DL3*SSI(ITYP) 364 : DSFG(ITYP)=DSFG(ITYP)*(1.D0-XT) 365 : 3010 EPDP(ITYP)=DL3*EPC(ITYP)+EPDP(ITYP) 365 : CALL CHREP(SFG,SIG,-ANG) 365 : CALL CHREP(DSFG,DSIG,-ANG) 365 : CALL CHREP(EPDP,EPPLDP,-ANG) 366 : XLAM1=XLAM1+DL1 367 : XLAM3=XLAM3+DL3 368 : IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER 369 : RETURN 370 : END 371 :
© Cast3M 2003 - All rights reserved.
Disclaimer