1 : C DRUCK2 SOURCE LJASON 07/11/12 21:15:13 5965 2 : C DRUCK2 3 : SUBROUTINE DRUCK2(SIG0,NSTRS,DEPST,VAR0,XMAT,IVAL,NCOMAT, 4 : & XCAR,ICARA,NVARI,SIGF,VARF,DEFP,MFR,KERRE, 5 : & IB,IGAU,IFOURB,XLCARA,MELE) 6 : C---------------------------------------------------------------------- 7 : C 8 : C ENTREES: 9 : C ------- 10 : C NSTRS = NBR. DE COMPOSANTES DES CONTR. OU DES DEFORM. 11 : C SIG0(NSTRS) = CONTR. AU DEBUT DU PAS D'INTEGRATION 12 : C DEPST(NSTRS) = INCREMENT DES DEFORM. TOTALES 13 : C NVARI = NBR. DE VARIABLES INTERNES 14 : C VAR0(NVARI) = VARIABLES INTERNES AU DEBUT DU PAS D'INTEGRATION 15 : C 16 : C NCOMAT = NBR. DE CARACTERISTIQUES MECANIQUES DU MATERIAU 17 : C IVAL(NCOMAT) = INDICE DES COMPOSANTES DE MATERIAU 18 : C XMAT(NCOMAT) = CARACTERISTIQUES MECANIQUES DU MATERIAU 19 : C MFR = INDICE DE LA FORMULATION MECANIQUE 20 : C ICARA = NBR. DE CARACT. GEOMETRIQUES DES ELEMENTS FINIS 21 : C XCAR(ICARA) = CARACT. GEOMETRIQUES DES ELEMENTS FINIS 22 : C 23 : C SORTIES: 24 : C ------- 25 : C SIGF(NSTRS)= CONTR. A LA FIN DU PAS D'INTEGRATION 26 : C VARF(NVARI)= VARIABLES INTERNES A LA FIN DU PAS D'INTEGRATION 27 : C DEFP(NSTRS)= INCREMENT DES DEFORM. PLASTIQUES A LA FIN DU PAS 28 : C D'INTEGRATION 29 : C KERRE = INDICE QUI REGIT LES ERREURS 30 : C = 99 SI LA FORMULATION MECANIQUE N'EST PAS DISPONIBLE 31 : C POUR LE MODELE CONSIDERE OU S'IL Y A INCOMPATIBILITE 32 : C ENTRE MFR ET IFOUR 33 : C 34 : C 35 : IMPLICIT REAL*8(A-H,O-Z) 36 : -INC CCOPTIO 37 : -INC CCREEL 38 : C 39 : PARAMETER (XZER=0.D0,UNDEMI=.5D0,UN=1.D0,DEUX=2.D0,TROIS=3.D0) 40 : PARAMETER (QUATRE=4.D0) 41 : C 42 : DIMENSION SIG0(*),DEPST(*),VAR0(*),XMAT(*),XCAR(*),SIGF(*), 43 : & IVAL(*),VARF(*),DEFP(*) 44 : DIMENSION DDAUX(6,6),S(6),SINI(3),DS(6),SD(3),SDINI(3) 45 : DIMENSION SIG33(3,3),SIGPP(3),R(3,3),RT(3,3),XMS(3,3),TRAV(3,3) 46 : 47 : C---------------------------------------------------------------------- 48 : C CARACTERISTIQUES 49 : C---------------------------------------------------------------------- 50 : 51 : KERRE=0 52 : 53 : YOUN=XMAT(1) 54 : XNU=XMAT(2) 55 : FT=XMAT(5) 56 : FC=XMAT(6) 57 : FB=XMAT(7) 58 : FY=XMAT(8) 59 : EPSU=XMAT(9) 60 : EPSU2=XMAT(10) 61 : XLC=XMAT(11) 62 : GF=XMAT(12) 63 : 64 : G=UNDEMI*YOUN/(UN+XNU) 65 : 66 : **-- Parametres pour le critere en compression 67 : 68 : XK=YOUN/TROIS/(UN-DEUX*XNU) 69 : ALPHA=SQRT(UN/TROIS)*(FC-FB)/(FC-DEUX*FB) 70 : XK0=SQRT(UN/TROIS)*FB/(DEUX*FB-FC) 71 : EPSPU=EPSU+FC/YOUN 72 : XU=SQRT(TROIS)*EPSPU/(SQRT(TROIS)*ALPHA-UN) 73 : OMEGA=FY/FC 74 : XU2=SQRT(TROIS)*EPSU2/(SQRT(TROIS)*ALPHA-UN) 75 : OM1=OMEGA-UN 76 : XUXU=XU*XU 77 : 78 : AUX3=G+9.D0*XK*ALPHA*ALPHA 79 : 80 : **-- Longueur caracteristique pour le critere en traction 81 : 82 : **-- Elements lineaires 83 : * IF (MELE.EQ.8) THEN 84 : * XLC=SQRT(DEUX*XLCARA) 85 : * ELSE IF (MELE.EQ.10) THEN 86 : **-- Elements quadratiques 87 : * XLC=SQRT(XLCARA) 88 : * ELSE 89 : * PRINT*,'pas le bon type d element' 90 : * STOP 91 : * ENDIF 92 : 93 : * PRINT*,'XLC1', XLC 94 : * XLC = XLCARA 95 : * PRINT*,'XLC', XLC 96 : * PRINT*,'XLCARA', XLCARA 97 : 98 : IF(IVAL(9).EQ.0) THEN 99 : XKU=XZER 100 : ELSE 101 : ** XKU=DEUX*GF/(XLC*FT) 102 : XKU=GF/(XLC*FT*(UN-EXP(-UN))) 103 : * PRINT*,'GF', GF 104 : * PRINT*,'FT', FT 105 : * PRINT*,'XLC', XLC 106 : * PRINT*,'XKU', XKU 107 : ENDIF 108 : 109 : * print*,'xku : ',xku 110 : 111 : 112 : C---------------------------------------------------------------------- 113 : C CALCUL DU PREDICTEUR ELASTIQUE 114 : C---------------------------------------------------------------------- 115 : 116 : **-- Matrice de Hooke DDAUX 117 : 118 : DO I=1,NSTRS 119 : DO J=1,NSTRS 120 : DDAUX(I,J)=XZER 121 : ENDDO 122 : ENDDO 123 : 124 : AUX0=UN/((UN+XNU)*(UN-DEUX*XNU)) 125 : D11=YOUN*AUX0*(UN-XNU) 126 : D12=YOUN*AUX0*XNU 127 : AUX2=YOUN*UNDEMI/(UN+XNU) 128 : 129 : DDAUX(1,1)=D11 130 : DDAUX(1,2)=D12 131 : DDAUX(1,3)=D12 132 : 133 : DDAUX(2,1)=D12 134 : DDAUX(2,2)=D11 135 : DDAUX(2,3)=D12 136 : 137 : DDAUX(3,1)=D12 138 : DDAUX(3,2)=D12 139 : DDAUX(3,3)=D11 140 : 141 : DO I=4,NSTRS 142 : DO J=4,NSTRS 143 : IF(I.EQ.J) DDAUX(I,J)=AUX2 144 : ENDDO 145 : ENDDO 146 : 147 : DO I=1,NSTRS 148 : SOM=XZER 149 : DO J=1,NSTRS 150 : SOM=SOM+DDAUX(I,J)*DEPST(J) 151 : ENDDO 152 : S(I)=SIG0(I)+SOM 153 : ENDDO 154 : 155 : C--- calcul des contraintes principales 156 : 156 : CALL ENDOCB(S,SIG33,1,IFOURB) 156 : CALL JACOB3 (SIG33,3,SIGPP,R) 157 : 158 : C---------------------------------------------------------------------- 159 : C CALCUL DES CRITERES PHI1 ET PHI2 160 : C PHI1 = DRUCKER-PRAGER , PHI2 = RANKINE 161 : C---------------------------------------------------------------------- 162 : 163 : ICOMPT=0 164 : PREC0=1.D-8 165 : PREC1=1.D-7 166 : XP1=VAR0(1) 167 : XP2=VAR0(2) 168 : DL1=XZER 169 : DL2=XZER 170 : DL2B=XZER 171 : DL2C=XZER 172 : 173 : DIFF=VAR0(1)-XU 174 : IF (VAR0(1).LT.XU) THEN 175 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 176 : ELSE IF ((VAR0(1).GE.XU).AND.(VAR0(1).LT.XU2)) THEN 177 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 178 : ELSE IF (VAR0(1).GE.XU2) THEN 179 : ZETA=XZER 180 : ENDIF 181 : 182 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 183 : 184 : **-- calcul du deviateur des contraintes SD 185 : 186 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 187 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 188 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 189 : 190 : **-- calcul de la contrainte equivalente SEQ =0.5*S*S 191 : 192 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 193 : & +(SIGPP(3)-SIGPP(1))**2 194 : SEQ=SQRT(SEQ/6.D0) 195 : 196 : PHI1=SEQ+ALPHA*XI1-ZETA 197 : 198 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 199 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 200 : PHI2C=SIGPP(3)-FT*TX2(XP2,XKU) 201 : 202 : 203 : PHI10=ABS(PHI1) 204 : PHI20=ABS(PHI2) 205 : 206 : **-- sauvegarde du predicteur 207 : 208 : F1INI=PHI1 209 : F2INI=PHI2 210 : F2BINI=PHI2B 211 : SEQI=SEQ 212 : DO I=1,3 213 : SINI(I)=SIGPP(I) 214 : SDINI(I)=SD(I) 215 : ENDDO 216 : 217 : IF((PHI1.LT.XZER).AND.(PHI2.LT.XZER)) THEN 218 : * print*,'cas elastique' 219 : C---------------------------------------------------------------------- 220 : C ON N A PAS PLASTIFIE 221 : C---------------------------------------------------------------------- 222 : DO I=1,NSTRS 223 : SIGF(I)=S(I) 224 : DEFP(I)=XZER 225 : ENDDO 226 : VARF(1)=VAR0(1) 227 : VARF(2)=VAR0(2) 228 : VARF(3)=VAR0(3) 229 : VARF(4)=VAR0(4) 230 : VARF(5)=VAR0(5) 231 : GO TO 999 232 : 233 : ELSE IF((PHI1.LT.XZER).AND.(PHI2.GE.XZER)) THEN 234 : * print*,'cas Rankine' 235 : IF(PHI2B.LT.XZER) THEN 236 : C---------------------------------------------------------------------- 237 : C RANKINE - UN CRITERE ACTIF (f2) 238 : C---------------------------------------------------------------------- 239 : *** debut bloc 20 240 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 241 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 242 : ICOMPT=ICOMPT+1 243 : 244 : **-- calcul du multiplicateur 245 : 246 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 247 : XP2=XP2+DDL2 248 : DL2=DL2+DDL2 249 : IF(DL2.LT.XZER) THEN 250 : PRINT*,'Multiplicateur negatif (Rankine 1 Critere)' 251 : STOP 252 : ENDIF 253 : 254 : SIGPP(1)=SIGPP(1)-DDL2*D11 255 : SIGPP(2)=SIGPP(2)-DDL2*D12 256 : SIGPP(3)=SIGPP(3)-DDL2*D12 257 : 258 : **-- calcul du nouveau critere phi2 259 : 260 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 261 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 262 : END DO 263 : **-- test de convergence 264 : 265 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 266 : & (ABS(PHI2).LT.PREC1)) THEN 267 : * PRINT*,'boucle rankine', ICOMPT 268 : GO TO 998 269 : ELSE 270 : IF(ICOMPT.GE.20) THEN 271 : PRINT*,'le nombre d iterations internes depasse 20 272 : & (Rankine 1 Critere)' 273 : 274 : STOP 275 : ENDIF 276 : ENDIF 277 : *** fin bloc 20 278 : ELSE IF (PHI2C.LT.XZER) THEN 279 : C---------------------------------------------------------------------- 280 : C RANKINE - 2 CRITERES ACTIFS (f2 et f2b) 281 : C---------------------------------------------------------------------- 282 : * debut bloc 21 283 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 284 : & ((ABS(PHI2).LT.(PREC0*PHI20)).OR.(ABS(PHI2).LT.PREC1))) 285 : ICOMPT=ICOMPT+1 286 : 287 : **-- calcul des 2 multiplicateurs 288 : 289 : A11=D11-DF2DX(XP2,XKU,FT) 290 : A12=D12-DF2DX(XP2,XKU,FT) 291 : XX=A11**2-A12**2 292 : DDL2=(A11*PHI2-A12*PHI2B)/XX 293 : DDL2B=(A11*PHI2B-A12*PHI2)/XX 294 : 295 : XP2=XP2+DDL2+DDL2B 296 : DL2=DL2+DDL2 297 : DL2B=DL2B+DDL2B 298 : 299 : IF(DL2.LT.XZER) THEN 300 : PRINT*,'Multiplicateur negatif (Rankine 2 Critere)' 301 : STOP 302 : ENDIF 303 : 304 : SIGPP(1)=SIGPP(1)-DDL2*D11-DDL2B*D12 305 : SIGPP(2)=SIGPP(2)-DDL2*D12-DDL2B*D11 306 : SIGPP(3)=SIGPP(3)-DDL2*D12-DDL2B*D12 307 : 308 : **-- calcul des nouveaux criteres phi2 et phi2b 309 : 310 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 311 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 312 : END DO 313 : **-- test de convergence 314 : 315 : IF((ABS(PHI2).LT.(PREC0*PHI20)) 316 : & .OR.(ABS(PHI2).LT.PREC1)) THEN 317 : 318 : IF (DL2B.LT.XZER) THEN 319 : ICOMPT=0 320 : DL2=XZER 321 : PHI2=F2INI 322 : XP2=VAR0(2) 323 : DO I=1,3 324 : SIGPP(I)= SINI(I) 325 : END DO 326 : * PRINT*,'** 2e Multiplicateur negatif **' 327 : ** appel bloc 20 328 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 329 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 330 : ICOMPT=ICOMPT+1 331 : 332 : **-- calcul du multiplicateur 333 : 334 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 335 : XP2=XP2+DDL2 336 : DL2=DL2+DDL2 337 : IF(DL2.LT.XZER) THEN 338 : PRINT*,'Multiplicateur negatif (Rankine 1 Critere)' 339 : STOP 340 : ENDIF 341 : 342 : SIGPP(1)=SIGPP(1)-DDL2*D11 343 : SIGPP(2)=SIGPP(2)-DDL2*D12 344 : SIGPP(3)=SIGPP(3)-DDL2*D12 345 : 346 : **-- calcul du nouveau critere phi2 347 : 348 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 349 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 350 : END DO 351 : **-- test de convergence 352 : 353 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 354 : & (ABS(PHI2).LT.PREC1)) THEN 355 : GO TO 998 356 : ELSE 357 : IF(ICOMPT.GE.20) THEN 358 : PRINT*,'le nombre d iterations internes depasse 20 359 : & (Rankine 1 Critere)' 360 : 361 : STOP 362 : ENDIF 363 : ENDIF 364 : *** fin bloc 20 365 : ENDIF 366 : GO TO 998 367 : ELSE IF(ICOMPT.GE.20) THEN 368 : PRINT*,'le nombre d iterations internes depasse 20 369 : & (Rankine 2 Critere)' 370 : STOP 371 : ENDIF 372 : Cfin bloc 21 373 : 374 : *** deux endif en attente 375 : 376 : ELSE 377 : C---------------------------------------------------------------------- 378 : C RANKINE - 3 CRITERES ACTIFS (f2, f2b et f2c) 379 : C---------------------------------------------------------------------- 380 : * début bloc 22 381 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 382 : & ((ABS(PHI2).LT.(PREC0*PHI20)).OR.(ABS(PHI2).LT.PREC1))) 383 : 384 : ICOMPT=ICOMPT+1 385 : 386 : **-- calcul des 3 multiplicateurs 387 : 388 : A11=D11-DF2DX(XP2,XKU,FT) 389 : A12=D12-DF2DX(XP2,XKU,FT) 390 : XX=A11**2-DEUX*A12**2+A11*A12 391 : DDL2=((A11+A12)*PHI2-A12*(PHI2B+PHI2C))/XX 392 : DDL2B=((A11+A12)*PHI2B-A12*(PHI2+PHI2C))/XX 393 : DDL2C=((A11+A12)*PHI2C-A12*(PHI2+PHI2B))/XX 394 : 395 : XP2=XP2+DDL2+DDL2B+DDL2C 396 : DL2=DL2+DDL2 397 : DL2B=DL2B+DDL2B 398 : DL2C=DL2C+DDL2C 399 : 400 : IF(DL2.LT.XZER) THEN 401 : PRINT*,'Multiplicateur negatif (Rankine 3 Critere)' 402 : STOP 403 : ENDIF 404 : 405 : SIGPP(1)=SIGPP(1)-DDL2*D11-DDL2B*D12-DDL2C*D12 406 : SIGPP(2)=SIGPP(2)-DDL2*D12-DDL2B*D11-DDL2C*D12 407 : SIGPP(3)=SIGPP(3)-DDL2*D12-DDL2B*D12-DDL2C*D11 408 : 409 : **-- calcul des nouveaux criteres phi2 et phi2b 410 : 411 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 412 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 413 : PHI2C=SIGPP(3)-FT*TX2(XP2,XKU) 414 : 415 : END DO 416 : **-- test de convergence 417 : 418 : IF((ABS(PHI2).LT.(PREC0*PHI20)) 419 : & .OR.(ABS(PHI2).LT.PREC1)) THEN 420 : 421 : IF ((DL2B.LT.XZER).AND.(DL2C.LT.XZER)) THEN 422 : ICOMPT=0 423 : DL2=XZER 424 : PHI2=F2INI 425 : XP2=VAR0(2) 426 : DO I=1,3 427 : SIGPP(I)=SINI(I) 428 : END DO 429 : * debut appel bloc 20 430 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 431 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 432 : ICOMPT=ICOMPT+1 433 : 434 : **-- calcul du multiplicateur 435 : 436 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 437 : XP2=XP2+DDL2 438 : DL2=DL2+DDL2 439 : IF(DL2.LT.XZER) THEN 440 : PRINT*,'Multiplicateur negatif (Rankine 1 Critere)' 441 : STOP 442 : ENDIF 443 : 444 : SIGPP(1)=SIGPP(1)-DDL2*D11 445 : SIGPP(2)=SIGPP(2)-DDL2*D12 446 : SIGPP(3)=SIGPP(3)-DDL2*D12 447 : 448 : **-- calcul du nouveau critere phi2 449 : 450 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 451 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 452 : END DO 453 : **-- test de convergence 454 : 455 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 456 : & (ABS(PHI2).LT.PREC1)) THEN 457 : GO TO 998 458 : ELSE 459 : IF(ICOMPT.GE.20) THEN 460 : PRINT*,'le nombre d iterations internes depasse 20 461 : & (Rankine 1 Critere)' 462 : 463 : STOP 464 : ENDIF 465 : ENDIF 466 : *** fin bloc 20 467 : 468 : ENDIF 469 : 470 : IF ((DL2B.GE.XZER).AND.(DL2C.LT.XZER)) THEN 471 : ICOMPT=0 472 : DL2=XZER 473 : DL2B=XZER 474 : PHI2=F2INI 475 : PHI2B=F2BINI 476 : XP2=VAR0(2) 477 : DO I=1,3 478 : SIGPP(I)=SINI(I) 479 : END DO 480 : * debut appel bloc 21 481 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 482 : & ((ABS(PHI2).LT.(PREC0*PHI20)).OR.(ABS(PHI2).LT.PREC1))) 483 : ICOMPT=ICOMPT+1 484 : 485 : **-- calcul des 2 multiplicateurs 486 : 487 : A11=D11-DF2DX(XP2,XKU,FT) 488 : A12=D12-DF2DX(XP2,XKU,FT) 489 : XX=A11**2-A12**2 490 : DDL2=(A11*PHI2-A12*PHI2B)/XX 491 : DDL2B=(A11*PHI2B-A12*PHI2)/XX 492 : 493 : XP2=XP2+DDL2+DDL2B 494 : DL2=DL2+DDL2 495 : DL2B=DL2B+DDL2B 496 : 497 : IF(DL2.LT.XZER) THEN 498 : PRINT*,'Multiplicateur negatif (Rankine 2 Critere)' 499 : STOP 500 : ENDIF 501 : 502 : SIGPP(1)=SIGPP(1)-DDL2*D11-DDL2B*D12 503 : SIGPP(2)=SIGPP(2)-DDL2*D12-DDL2B*D11 504 : SIGPP(3)=SIGPP(3)-DDL2*D12-DDL2B*D12 505 : 506 : **-- calcul des nouveaux criteres phi2 et phi2b 507 : 508 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 509 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 510 : END DO 511 : **-- test de convergence 512 : 513 : IF((ABS(PHI2).LT.(PREC0*PHI20)) 514 : & .OR.(ABS(PHI2).LT.PREC1)) THEN 515 : 516 : IF (DL2B.LT.XZER) THEN 517 : ICOMPT=0 518 : DL2=XZER 519 : PHI2=F2INI 520 : XP2=VAR0(2) 521 : DO I=1,3 522 : SIGPP(I)= SINI(I) 523 : END DO 524 : * PRINT*,'** 2e Multiplicateur negatif **' 525 : ** appel bloc 20 526 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 527 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 528 : ICOMPT=ICOMPT+1 529 : 530 : **-- calcul du multiplicateur 531 : 532 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 533 : XP2=XP2+DDL2 534 : DL2=DL2+DDL2 535 : IF(DL2.LT.XZER) THEN 536 : PRINT*,'Multiplicateur negatif (Rankine 1 Critere)' 537 : STOP 538 : ENDIF 539 : 540 : SIGPP(1)=SIGPP(1)-DDL2*D11 541 : SIGPP(2)=SIGPP(2)-DDL2*D12 542 : SIGPP(3)=SIGPP(3)-DDL2*D12 543 : 544 : **-- calcul du nouveau critere phi2 545 : 546 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 547 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 548 : END DO 549 : **-- test de convergence 550 : 551 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 552 : & (ABS(PHI2).LT.PREC1)) THEN 553 : GO TO 998 554 : ELSE 555 : IF(ICOMPT.GE.20) THEN 556 : PRINT*,'le nombre d iterations internes depasse 20 557 : & (Rankine 1 Critere)' 558 : 559 : STOP 560 : ENDIF 561 : ENDIF 562 : *** fin bloc 20 563 : ENDIF 564 : GO TO 998 565 : ELSE IF(ICOMPT.GE.20) THEN 566 : PRINT*,'le nombre d iterations internes depasse 20 567 : & (Rankine 2 Critere)' 568 : STOP 569 : ENDIF 570 : Cfin bloc 21 571 : ENDIF 572 : 573 : GO TO 998 574 : ELSE IF(ICOMPT.GE.20) THEN 575 : PRINT*,'le nombre d iterations internes depasse 20 576 : & (Rankine 3 Critere)' 577 : STOP 578 : ENDIF 579 : * fin bloc 22 580 : 581 : ENDIF 582 : 583 : ELSE IF((PHI1.GE.XZER).AND.(PHI2.LT.XZER)) THEN 584 : 585 : C---------------------------------------------------------------------- 586 : C DRUCKER-PRAGER 587 : C---------------------------------------------------------------------- 588 : ** début bloc 9 589 : IF(XP1.GE.XU2) THEN 590 : DO I=1,NSTRS 591 : SIGF(I)=XZER 592 : DEFP(I)=DEPST(I) 593 : END DO 594 : DDL1=PHI1/AUX3 595 : VARF(1)=XP1+DDL1 596 : VARF(2)=XP2 597 : VARF(3)=VARF(1)+VARF(2) 598 : VARF(4)=VARF(1)/XU 599 : VARF(5)=EXP(-VARF(2)/XKU)/FT 600 : GO TO 999 601 : ENDIF 602 : *** debut bloc 10 603 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 604 : & ((ABS(PHI1).LT.(PREC0*PHI10)).OR.(ABS(PHI1).LT.PREC1))) 605 : ICOMPT=ICOMPT+1 606 : 607 : DIFF=XP1-XU 608 : IF (XP1.LT.XU) THEN 609 : DFDX=-DEUX*XK0*FC*OM1*DIFF/XUXU 610 : ELSE IF((XP1.GE.XU)) THEN 611 : DFDX=DEUX*XK0*FC*DIFF/((XU2-XU)**2) 612 : ENDIF 613 : 614 : **-- calcul du multiplicateur 615 : 616 : DDL1=PHI1/(AUX3-DFDX) 617 : XP1=XP1+DDL1 618 : DL1=DL1+DDL1 619 : 620 : IF(XP1.GE.XU2) THEN 621 : DO I=1,NSTRS 622 : SIGF(I)=XZER 623 : DEFP(I)=DEPST(I) 624 : ENDDO 625 : VARF(1)=XP1-DDL1+(PHI1/AUX3) 626 : VARF(2)=XP2 627 : VARF(3)=VARF(1)+VARF(2) 628 : VARF(4)=VARF(1)/XU 629 : VARF(5)=EXP(-VARF(2)/XKU)/FT 630 : GO TO 999 631 : ENDIF 632 : 633 : SIGPP(1)=SIGPP(1)-DDL1*(G*SD(1)/SEQ+TROIS*XK*ALPHA) 634 : SIGPP(2)=SIGPP(2)-DDL1*(G*SD(2)/SEQ+TROIS*XK*ALPHA) 635 : SIGPP(3)=SIGPP(3)-DDL1*(G*SD(3)/SEQ+TROIS*XK*ALPHA) 636 : 637 : **-- calcul du nouveau critere phi1 638 : 639 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 640 : 641 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 642 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 643 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 644 : 645 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 646 : & +(SIGPP(3)-SIGPP(1))**2 647 : SEQ=SQRT(SEQ/6.D0) 648 : 649 : DIFF=XP1-XU 650 : IF (XP1.LT.XU) THEN 651 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 652 : ELSE IF ((XP1.GE.XU)) THEN 653 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 654 : ENDIF 655 : PHI1=SEQ+ALPHA*XI1-ZETA 656 : END DO 657 : 658 : **-- test de convergence 659 : 660 : IF((ABS(PHI1).LT.(PREC0*PHI10)).OR.(ABS(PHI1).LT.PREC1))THEN 661 : GO TO 998 662 : ELSE IF(ICOMPT.GE.20) THEN 663 : PRINT*,'le nombre d iterations internes depasse 20 664 : & (Drucker-Prager 1 Critere)' 665 : STOP 666 : ENDIF 667 : *** fin bloc 10 668 : *** fin bloc 9 669 : 670 : ELSE IF((PHI1.GE.XZER).AND.(PHI2.GE.XZER)) THEN 671 : 672 : IF(PHI2B.LT.XZER) THEN 673 : 674 : C----------------------------------------------------------------------- 675 : C RANKINE ET DRUCKER-PRAGER - 2 CRITERES ACTIFS (f1 et f2) 676 : C----------------------------------------------------------------------- 677 : * debut bloc 29 678 : IF(XP1.GE.XU2) THEN 679 : DO I=1,NSTRS 680 : SIGF(I)=XZER 681 : DEFP(I)=DEPST(I) 682 : END DO 683 : DDL1=PHI1/AUX3 684 : VARF(1)=XP1+DDL1 685 : VARF(2)=XP2 686 : VARF(3)=VARF(1)+VARF(2) 687 : VARF(4)=VARF(1)/XU 688 : VARF(5)=EXP(-VARF(2)/XKU)/FT 689 : GO TO 999 690 : ENDIF 691 : 692 : * debut bloc 30 693 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 694 : & ((ABS(PHI1).LT.PREC0*PHI10).AND.(ABS(PHI2).LT.PREC0*PHI20) 695 : & .OR.(ABS(PHI1).LT.PREC1).AND.(ABS(PHI2).LT.PREC1))) 696 : ICOMPT=ICOMPT+1 697 : DIFF=XP1-XU 698 : IF (XP1.LT.XU) THEN 699 : DFDX=-DEUX*XK0*FC*OM1*DIFF/XUXU 700 : ELSE IF((XP1.GE.XU)) THEN 701 : DFDX=DEUX*XK0*FC*DIFF/((XU2-XU)**2) 702 : ENDIF 703 : 704 : **-- calcul des 2 multiplicateurs dl1 et dl2 705 : 706 : A11=AUX3-DFDX 707 : A12=G*SD(1)/SEQ+TROIS*XK*ALPHA 708 : A22=D11-DF2DX(XP2,XKU,FT) 709 : XX=A11*A22-A12**2 710 : DDL1=(A22*PHI1-A12*PHI2)/XX 711 : DDL2=(A11*PHI2-A12*PHI1)/XX 712 : 713 : XP1=XP1+DDL1 714 : DL1=DL1+DDL1 715 : XP2=XP2+DDL2 716 : DL2=DL2+DDL2 717 : 718 : IF(XP1.GE.XU2) THEN 719 : DO I=1,NSTRS 720 : SIGF(I)=XZER 721 : DEFP(I)=DEPST(I) 722 : END DO 723 : VARF(1)=XP1-DDL1+(PHI1/AUX3) 724 : VARF(2)=XP2 725 : VARF(3)=VARF(1)+VARF(2) 726 : VARF(4)=VARF(1)/XU 727 : VARF(5)=EXP(-VARF(2)/XKU)/FT 728 : GO TO 999 729 : ENDIF 730 : 731 : SIGPP(1)=SIGPP(1)-DDL1*(G*SD(1)/SEQ+TROIS*XK*ALPHA) 732 : & -DDL2*D11 733 : SIGPP(2)=SIGPP(2)-DDL1*(G*SD(2)/SEQ+TROIS*XK*ALPHA) 734 : & -DDL2*D12 735 : SIGPP(3)=SIGPP(3)-DDL1*(G*SD(3)/SEQ+TROIS*XK*ALPHA) 736 : & -DDL2*D12 737 : 738 : **-- calcul des nouveaux criteres phi1 et phi2 739 : 740 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 741 : 742 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 743 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 744 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 745 : 746 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 747 : & +(SIGPP(3)-SIGPP(1))**2 748 : SEQ=SQRT(SEQ/6.D0) 749 : 750 : DIFF=XP1-XU 751 : IF (XP1.LT.XU) THEN 752 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 753 : ELSE IF ((XP1.GE.XU)) THEN 754 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 755 : ENDIF 756 : PHI1=SEQ+ALPHA*XI1-ZETA 757 : 758 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 759 : 760 : END DO 761 : **-- test de convergence 762 : 763 : IF((ABS(PHI1).LT.PREC0*PHI10).AND.(ABS(PHI2).LT.PREC0*PHI20) 764 : & .OR.(ABS(PHI1).LT.PREC1).AND.(ABS(PHI2).LT.PREC1)) THEN 765 : 766 : IF(DL1.LT.XZER) THEN 767 : * PRINT*,'1er multiplicateur negatif (D.P.-Rankine)' 768 : ICOMPT=0 769 : DL2=XZER 770 : PHI2=F2INI 771 : XP1=VAR0(1) 772 : XP2=VAR0(2) 773 : DO I=1,3 774 : SIGPP(I)=SINI(I) 775 : END DO 776 : * debut appel bloc 20 777 : *** debut bloc 20 778 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 779 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 780 : ICOMPT=ICOMPT+1 781 : 782 : **-- calcul du multiplicateur 783 : 784 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 785 : XP2=XP2+DDL2 786 : DL2=DL2+DDL2 787 : IF(DL2.LT.XZER) THEN 788 : PRINT*,'Multiplicateur negatif (Rankine 1 Critere)' 789 : STOP 790 : ENDIF 791 : 792 : SIGPP(1)=SIGPP(1)-DDL2*D11 793 : SIGPP(2)=SIGPP(2)-DDL2*D12 794 : SIGPP(3)=SIGPP(3)-DDL2*D12 795 : 796 : **-- calcul du nouveau critere phi2 797 : 798 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 799 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 800 : END DO 801 : **-- test de convergence 802 : 803 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 804 : & (ABS(PHI2).LT.PREC1)) THEN 805 : GO TO 998 806 : ELSE 807 : IF(ICOMPT.GE.20) THEN 808 : PRINT*,'le nombre d iterations internes depasse 20 809 : & (Rankine 1 Critere)' 810 : 811 : STOP 812 : ENDIF 813 : ENDIF 814 : *** fin bloc 20 815 : * fin appel bloc 20 816 : ENDIF 817 : 818 : IF (DL2.LT.XZER) THEN 819 : * PRINT*,'2eme multiplicateur negatif (D.P.-Rankine)' 820 : ICOMPT=0 821 : DL1=XZER 822 : PHI1=F1INI 823 : XP1=VAR0(1) 824 : XP2=VAR0(2) 825 : SEQ=SEQI 826 : DO I=1,3 827 : SIGPP(I)=SINI(I) 828 : SD(I)=SDINI(I) 829 : END DO 830 : * debut appel bloc 9 831 : ** début bloc 9 832 : IF(XP1.GE.XU2) THEN 833 : DO I=1,NSTRS 834 : SIGF(I)=XZER 835 : DEFP(I)=DEPST(I) 836 : END DO 837 : DDL1=PHI1/AUX3 838 : VARF(1)=XP1+DDL1 839 : VARF(2)=XP2 840 : VARF(3)=VARF(1)+VARF(2) 841 : VARF(4)=VARF(1)/XU 842 : VARF(5)=EXP(-VARF(2)/XKU)/FT 843 : GO TO 999 844 : ENDIF 845 : *** debut bloc 10 846 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 847 : & ((ABS(PHI1).LT.(PREC0*PHI10)).OR.(ABS(PHI1).LT.PREC1))) 848 : ICOMPT=ICOMPT+1 849 : 850 : DIFF=XP1-XU 851 : IF (XP1.LT.XU) THEN 852 : DFDX=-DEUX*XK0*FC*OM1*DIFF/XUXU 853 : ELSE IF((XP1.GE.XU)) THEN 854 : DFDX=DEUX*XK0*FC*DIFF/((XU2-XU)**2) 855 : ENDIF 856 : 857 : **-- calcul du multiplicateur 858 : 859 : DDL1=PHI1/(AUX3-DFDX) 860 : XP1=XP1+DDL1 861 : DL1=DL1+DDL1 862 : 863 : IF(XP1.GE.XU2) THEN 864 : DO I=1,NSTRS 865 : SIGF(I)=XZER 866 : DEFP(I)=DEPST(I) 867 : END DO 868 : VARF(1)=XP1-DDL1+(PHI1/AUX3) 869 : VARF(2)=XP2 870 : VARF(3)=VARF(1)+VARF(2) 871 : VARF(4)=VARF(1)/XU 872 : VARF(5)=EXP(-VARF(2)/XKU)/FT 873 : GO TO 999 874 : ENDIF 875 : 876 : SIGPP(1)=SIGPP(1)-DDL1*(G*SD(1)/SEQ+TROIS*XK*ALPHA) 877 : SIGPP(2)=SIGPP(2)-DDL1*(G*SD(2)/SEQ+TROIS*XK*ALPHA) 878 : SIGPP(3)=SIGPP(3)-DDL1*(G*SD(3)/SEQ+TROIS*XK*ALPHA) 879 : 880 : **-- calcul du nouveau critere phi1 881 : 882 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 883 : 884 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 885 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 886 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 887 : 888 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 889 : & +(SIGPP(3)-SIGPP(1))**2 890 : SEQ=SQRT(SEQ/6.D0) 891 : 892 : DIFF=XP1-XU 893 : IF (XP1.LT.XU) THEN 894 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 895 : ELSE IF ((XP1.GE.XU)) THEN 896 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 897 : ENDIF 898 : PHI1=SEQ+ALPHA*XI1-ZETA 899 : END DO 900 : 901 : **-- test de convergence 902 : 903 : IF((ABS(PHI1).LT.(PREC0*PHI10)).OR. 904 : & (ABS(PHI1).LT.PREC1))THEN 905 : GO TO 998 906 : ELSE IF(ICOMPT.GE.20) THEN 907 : PRINT*,'le nombre d iterations internes depasse 20 908 : & (Drucker-Prager 1 Critere)' 909 : STOP 910 : ENDIF 911 : *** fin bloc 10 912 : *** fin bloc 9 913 : * fin appel bloc 9 914 : ENDIF 915 : GO TO 998 916 : ELSE IF(ICOMPT.GE.20) THEN 917 : PRINT*,'le nombre d iterations internes depasse 20 918 : & (D.P.-Rankine)' 919 : STOP 920 : ENDIF 921 : * fin bloc 30 922 : * fin bloc 29 923 : ENDIF 924 : 925 : ELSE 926 : 927 : C----------------------------------------------------------------------- 928 : C RANKINE ET DRUCKER-PRAGER - 3 CRITERES ACTIFS (f1, f2 et f2b) 929 : C----------------------------------------------------------------------- 930 : 931 : IF(XP1.GE.XU2) THEN 932 : * print*,'On est au dessus de xu2' 933 : DO I=1,NSTRS 934 : SIGF(I)=XZER 935 : DEFP(I)=DEPST(I) 936 : END DO 937 : DDL1=PHI1/AUX3 938 : VARF(1)=XP1+DDL1 939 : VARF(2)=XP2 940 : VARF(3)=VARF(1)+VARF(2) 941 : VARF(4)=VARF(1)/XU 942 : VARF(5)=EXP(-VARF(2)/XKU)/FT 943 : GO TO 999 944 : ENDIF 945 : ** debut bloc 31 946 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 947 : & ((ABS(PHI1).LT.PREC0*PHI10).AND.(ABS(PHI2).LT.PREC0*PHI20) 948 : & .OR.(ABS(PHI1).LT.PREC1).AND.(ABS(PHI2).LT.PREC1))) 949 : 950 : ICOMPT=ICOMPT+1 951 : DIFF=XP1-XU 952 : IF (XP1.LT.XU) THEN 953 : DFDX=-DEUX*XK0*FC*OM1*DIFF/XUXU 954 : ELSE IF((XP1.GE.XU)) THEN 955 : DFDX=DEUX*XK0*FC*DIFF/((XU2-XU)**2) 956 : ENDIF 957 : 958 : **-- calcul des 3 multiplicateurs dl1, dl2 et dl2b 959 : 960 : A11=AUX3-DFDX 961 : A12=G*SD(1)/SEQ+TROIS*XK*ALPHA 962 : A13=G*SD(2)/SEQ+TROIS*XK*ALPHA 963 : A22=D11-DF2DX(XP2,XKU,FT) 964 : A23=D12-DF2DX(XP2,XKU,FT) 965 : XX=A11*A22**2-A22*A13**2-A22*A12**2-A11*A23**2 966 : & +DEUX*A23*A12*A13 967 : DDL1=((A22**2-A23**2)*PHI1+(A13*A23-A12*A22)*PHI2 968 : & +(A12*A23-A13*A22)*PHI2B)/XX 969 : DDL2=((A23*A13-A12*A22)*PHI1+(A11*A22-A13**2)*PHI2 970 : & +(A13*A12-A11*A23)*PHI2B)/XX 971 : DDL2B=((A23*A12-A22*A13)*PHI1+(A13*A12-A11*A23)*PHI2 972 : & +(A11*A22-A12**2)*PHI2B)/XX 973 : XP1=XP1+DDL1 974 : DL1=DL1+DDL1 975 : XP2=XP2+DDL2+DDL2B 976 : DL2=DL2+DDL2 977 : DL2B=DL2B+DDL2B 978 : 979 : 980 : IF(XP1.GE.XU2) THEN 981 : * print*,'On est au dessus de xu2' 982 : DO I=1,NSTRS 983 : SIGF(I)=XZER 984 : DEFP(I)=DEPST(I) 985 : END DO 986 : VARF(1)=XP1-DDL1+(PHI1/AUX3) 987 : VARF(2)=XP2 988 : VARF(3)=VARF(1)+VARF(2) 989 : VARF(4)=VARF(1)/XU 990 : VARF(5)=EXP(-VARF(2)/XKU)/FT 991 : GO TO 999 992 : ENDIF 993 : 994 : SIGPP(1)=SIGPP(1)-DDL1*(G*SD(1)/SEQ+TROIS*XK*ALPHA)-DDL2*D11 995 : & -DDL2B*D12 996 : SIGPP(2)=SIGPP(2)-DDL1*(G*SD(2)/SEQ+TROIS*XK*ALPHA)-DDL2*D12 997 : & -DDL2B*D11 998 : SIGPP(3)=SIGPP(3)-DDL1*(G*SD(3)/SEQ+TROIS*XK*ALPHA)-DDL2*D12 999 : & -DDL2B*D12 1000 : 1001 : **-- calcul des nouveaux criteres phi1 et phi2 1002 : 1003 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 1004 : 1005 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 1006 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 1007 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 1008 : 1009 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 1010 : & +(SIGPP(3)-SIGPP(1))**2 1011 : SEQ=SQRT(SEQ/6.D0) 1012 : 1013 : DIFF=XP1-XU 1014 : IF (XP1.LT.XU) THEN 1015 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 1016 : ELSE IF ((XP1.GE.XU)) THEN 1017 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 1018 : ENDIF 1019 : PHI1=SEQ+ALPHA*XI1-ZETA 1020 : 1021 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 1022 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 1023 : END DO 1024 : 1025 : **-- test de convergence 1026 : 1027 : IF((ABS(PHI1).LT.PREC0*PHI10).AND.(ABS(PHI2).LT.PREC0*PHI20) 1028 : & .OR.(ABS(PHI1).LT.PREC1).AND.(ABS(PHI2).LT.PREC1)) THEN 1029 : 1030 : IF(DL1.LT.XZER) THEN 1031 : IF(DL2.LT.XZER)THEN 1032 : print*,'mutiplicateurs negatifs' 1033 : STOP 1034 : ELSE 1035 : IF(DL2B.LT.XZER)THEN 1036 : * PRINT*,'multiplicateur negatif (D.P.-Rankine) 20' 1037 : ICOMPT=0 1038 : DL2=XZER 1039 : PHI2=F2INI 1040 : XP1=VAR0(1) 1041 : XP2=VAR0(2) 1042 : DO I=1,3 1043 : SIGPP(I)=SINI(I) 1044 : END DO 1045 : *** debut appel bloc 20 1046 : *** debut bloc 20 1047 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 1048 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 1049 : ICOMPT=ICOMPT+1 1050 : 1051 : **-- calcul du multiplicateur 1052 : 1053 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 1054 : XP2=XP2+DDL2 1055 : DL2=DL2+DDL2 1056 : IF(DL2.LT.XZER) THEN 1057 : PRINT*,'Multiplicateur negatif (Rankine 1 Critere)' 1058 : STOP 1059 : ENDIF 1060 : 1061 : SIGPP(1)=SIGPP(1)-DDL2*D11 1062 : SIGPP(2)=SIGPP(2)-DDL2*D12 1063 : SIGPP(3)=SIGPP(3)-DDL2*D12 1064 : 1065 : **-- calcul du nouveau critere phi2 1066 : 1067 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 1068 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 1069 : END DO 1070 : **-- test de convergence 1071 : 1072 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 1073 : & (ABS(PHI2).LT.PREC1)) THEN 1074 : GO TO 998 1075 : ELSE 1076 : IF(ICOMPT.GE.20) THEN 1077 : PRINT*,'le nombre d iterations internes depasse 20 1078 : & (Rankine 1 Critere)' 1079 : 1080 : STOP 1081 : ENDIF 1082 : ENDIF 1083 : *** fin bloc 20 1084 : *** fin appel bloc 20 1085 : ELSE 1086 : * PRINT*,'multiplicateur negatif (D.P.-Rankine) 21 ' 1087 : ICOMPT=0 1088 : DL2=XZER 1089 : DL2B=XZER 1090 : PHI2=F2INI 1091 : PHI2B=F2BINI 1092 : XP1=VAR0(1) 1093 : XP2=VAR0(2) 1094 : DO I=1,3 1095 : SIGPP(I)=SINI(I) 1096 : END DO 1097 : * debut appel bloc 21 1098 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 1099 : & ((ABS(PHI2).LT.(PREC0*PHI20)).OR.(ABS(PHI2).LT.PREC1))) 1100 : ICOMPT=ICOMPT+1 1101 : 1102 : **-- calcul des 2 multiplicateurs 1103 : 1104 : A11=D11-DF2DX(XP2,XKU,FT) 1105 : A12=D12-DF2DX(XP2,XKU,FT) 1106 : XX=A11**2-A12**2 1107 : DDL2=(A11*PHI2-A12*PHI2B)/XX 1108 : DDL2B=(A11*PHI2B-A12*PHI2)/XX 1109 : 1110 : XP2=XP2+DDL2+DDL2B 1111 : DL2=DL2+DDL2 1112 : DL2B=DL2B+DDL2B 1113 : 1114 : IF(DL2.LT.XZER) THEN 1115 : PRINT*,'Multiplicateur negatif (Rankine 2 Critere)' 1116 : STOP 1117 : ENDIF 1118 : 1119 : SIGPP(1)=SIGPP(1)-DDL2*D11-DDL2B*D12 1120 : SIGPP(2)=SIGPP(2)-DDL2*D12-DDL2B*D11 1121 : SIGPP(3)=SIGPP(3)-DDL2*D12-DDL2B*D12 1122 : 1123 : **-- calcul des nouveaux criteres phi2 et phi2b 1124 : 1125 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 1126 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 1127 : END DO 1128 : **-- test de convergence 1129 : 1130 : IF((ABS(PHI2).LT.(PREC0*PHI20)) 1131 : & .OR.(ABS(PHI2).LT.PREC1)) THEN 1132 : 1133 : IF (DL2B.LT.XZER) THEN 1134 : ICOMPT=0 1135 : DL2=XZER 1136 : PHI2=F2INI 1137 : XP2=VAR0(2) 1138 : DO I=1,3 1139 : SIGPP(I)= SINI(I) 1140 : END DO 1141 : * PRINT*,'** 2e Multiplicateur negatif **' 1142 : ** appel bloc 20 1143 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 1144 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 1145 : ICOMPT=ICOMPT+1 1146 : 1147 : **-- calcul du multiplicateur 1148 : 1149 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 1150 : XP2=XP2+DDL2 1151 : DL2=DL2+DDL2 1152 : IF(DL2.LT.XZER) THEN 1153 : PRINT*,'Multiplicateur negatif 1154 : & (Rankine 1 Critere)' 1155 : STOP 1156 : ENDIF 1157 : 1158 : SIGPP(1)=SIGPP(1)-DDL2*D11 1159 : SIGPP(2)=SIGPP(2)-DDL2*D12 1160 : SIGPP(3)=SIGPP(3)-DDL2*D12 1161 : 1162 : **-- calcul du nouveau critere phi2 1163 : 1164 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 1165 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 1166 : END DO 1167 : **-- test de convergence 1168 : 1169 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 1170 : & (ABS(PHI2).LT.PREC1)) THEN 1171 : GO TO 998 1172 : ELSE 1173 : IF(ICOMPT.GE.20) THEN 1174 : PRINT*,'le nombre d iterations 1175 : & internes depasse 20 1176 : & (Rankine 1 Critere)' 1177 : 1178 : STOP 1179 : ENDIF 1180 : ENDIF 1181 : *** fin bloc 20 1182 : ENDIF 1183 : GO TO 998 1184 : ELSE IF(ICOMPT.GE.20) THEN 1185 : PRINT*,'le nombre d iterations internes depasse 20 1186 : & (Rankine 2 Critere)' 1187 : STOP 1188 : ENDIF 1189 : Cfin bloc 21 1190 : * fin appel bloc 21 1191 : ENDIF 1192 : ENDIF 1193 : ELSE 1194 : 1195 : IF (DL2.LT.XZER) THEN 1196 : * PRINT*,'multiplicateur negatif (D.P.-Rankine) 9' 1197 : ICOMPT=0 1198 : DL1=XZER 1199 : PHI1=F1INI 1200 : XP1=VAR0(1) 1201 : XP2=VAR0(2) 1202 : SEQ=SEQI 1203 : DO I=1,3 1204 : SIGPP(I)=SINI(I) 1205 : SD(I)=SDINI(I) 1206 : END DO 1207 : * debut appel bloc 9 1208 : ** début bloc 9 1209 : IF(XP1.GE.XU2) THEN 1210 : DO I=1,NSTRS 1211 : SIGF(I)=XZER 1212 : DEFP(I)=DEPST(I) 1213 : END DO 1214 : DDL1=PHI1/AUX3 1215 : VARF(1)=XP1+DDL1 1216 : VARF(2)=XP2 1217 : VARF(3)=VARF(1)+VARF(2) 1218 : VARF(4)=VARF(1)/XU 1219 : VARF(5)=EXP(-VARF(2)/XKU)/FT 1220 : GO TO 999 1221 : ENDIF 1222 : *** debut bloc 10 1223 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 1224 : & ((ABS(PHI1).LT.(PREC0*PHI10)).OR.(ABS(PHI1).LT.PREC1))) 1225 : ICOMPT=ICOMPT+1 1226 : 1227 : DIFF=XP1-XU 1228 : IF (XP1.LT.XU) THEN 1229 : DFDX=-DEUX*XK0*FC*OM1*DIFF/XUXU 1230 : ELSE IF((XP1.GE.XU)) THEN 1231 : DFDX=DEUX*XK0*FC*DIFF/((XU2-XU)**2) 1232 : ENDIF 1233 : 1234 : **-- calcul du multiplicateur 1235 : 1236 : DDL1=PHI1/(AUX3-DFDX) 1237 : XP1=XP1+DDL1 1238 : DL1=DL1+DDL1 1239 : 1240 : IF(XP1.GE.XU2) THEN 1241 : DO I=1,NSTRS 1242 : SIGF(I)=XZER 1243 : DEFP(I)=DEPST(I) 1244 : ENDDO 1245 : VARF(1)=XP1-DDL1+(PHI1/AUX3) 1246 : VARF(2)=XP2 1247 : VARF(3)=VARF(1)+VARF(2) 1248 : VARF(4)=VARF(1)/XU 1249 : VARF(5)=EXP(-VARF(2)/XKU)/FT 1250 : GO TO 999 1251 : ENDIF 1252 : 1253 : SIGPP(1)=SIGPP(1)-DDL1*(G*SD(1)/SEQ+TROIS*XK*ALPHA) 1254 : SIGPP(2)=SIGPP(2)-DDL1*(G*SD(2)/SEQ+TROIS*XK*ALPHA) 1255 : SIGPP(3)=SIGPP(3)-DDL1*(G*SD(3)/SEQ+TROIS*XK*ALPHA) 1256 : 1257 : **-- calcul du nouveau critere phi1 1258 : 1259 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 1260 : 1261 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 1262 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 1263 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 1264 : 1265 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 1266 : & +(SIGPP(3)-SIGPP(1))**2 1267 : SEQ=SQRT(SEQ/6.D0) 1268 : 1269 : DIFF=XP1-XU 1270 : IF (XP1.LT.XU) THEN 1271 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 1272 : ELSE IF ((XP1.GE.XU)) THEN 1273 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 1274 : ENDIF 1275 : PHI1=SEQ+ALPHA*XI1-ZETA 1276 : END DO 1277 : 1278 : **-- test de convergence 1279 : 1280 : IF((ABS(PHI1).LT.(PREC0*PHI10)).OR.(ABS(PHI1).LT.PREC1))THEN 1281 : GO TO 998 1282 : ELSE IF(ICOMPT.GE.20) THEN 1283 : PRINT*,'le nombre d iterations internes depasse 20 1284 : & (Drucker-Prager 1 Critere)' 1285 : STOP 1286 : ENDIF 1287 : *** fin bloc 10 1288 : *** fin bloc 9 1289 : * fin appel bloc 9 1290 : ELSE 1291 : IF(DL2B.LT.XZER)THEN 1292 : * PRINT*,'multiplicateur negatif (D.P.-Rankine) 29 ' 1293 : ICOMPT=0 1294 : DL1=XZER 1295 : DL2=XZER 1296 : PHI1=F1INI 1297 : PHI2=F2INI 1298 : XP1=VAR0(1) 1299 : XP2=VAR0(2) 1300 : SEQ=SEQI 1301 : DO I=1,3 1302 : SIGPP(I)=SINI(I) 1303 : SD(I)=SDINI(I) 1304 : END DO 1305 : * debut appel bloc 29 1306 : IF(XP1.GE.XU2) THEN 1307 : DO I=1,NSTRS 1308 : SIGF(I)=XZER 1309 : DEFP(I)=DEPST(I) 1310 : END DO 1311 : DDL1=PHI1/AUX3 1312 : VARF(1)=XP1+DDL1 1313 : VARF(2)=XP2 1314 : VARF(3)=VARF(1)+VARF(2) 1315 : VARF(4)=VARF(1)/XU 1316 : VARF(5)=EXP(-VARF(2)/XKU)/FT 1317 : GO TO 999 1318 : ENDIF 1319 : 1320 : * debut bloc 30 1321 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 1322 : & ((ABS(PHI1).LT.PREC0*PHI10).AND.(ABS(PHI2).LT.PREC0*PHI20) 1323 : & .OR.(ABS(PHI1).LT.PREC1).AND.(ABS(PHI2).LT.PREC1))) 1324 : ICOMPT=ICOMPT+1 1325 : DIFF=XP1-XU 1326 : IF (XP1.LT.XU) THEN 1327 : DFDX=-DEUX*XK0*FC*OM1*DIFF/XUXU 1328 : ELSE IF((XP1.GE.XU)) THEN 1329 : DFDX=DEUX*XK0*FC*DIFF/((XU2-XU)**2) 1330 : ENDIF 1331 : 1332 : **-- calcul des 2 multiplicateurs dl1 et dl2 1333 : 1334 : A11=AUX3-DFDX 1335 : A12=G*SD(1)/SEQ+TROIS*XK*ALPHA 1336 : A22=D11-DF2DX(XP2,XKU,FT) 1337 : XX=A11*A22-A12**2 1338 : DDL1=(A22*PHI1-A12*PHI2)/XX 1339 : DDL2=(A11*PHI2-A12*PHI1)/XX 1340 : 1341 : XP1=XP1+DDL1 1342 : DL1=DL1+DDL1 1343 : XP2=XP2+DDL2 1344 : DL2=DL2+DDL2 1345 : 1346 : IF(XP1.GE.XU2) THEN 1347 : DO I=1,NSTRS 1348 : SIGF(I)=XZER 1349 : DEFP(I)=DEPST(I) 1350 : END DO 1351 : VARF(1)=XP1-DDL1+(PHI1/AUX3) 1352 : VARF(2)=XP2 1353 : VARF(3)=VARF(1)+VARF(2) 1354 : VARF(4)=VARF(1)/XU 1355 : VARF(5)=EXP(-VARF(2)/XKU)/FT 1356 : GO TO 999 1357 : ENDIF 1358 : 1359 : SIGPP(1)=SIGPP(1)-DDL1*(G*SD(1)/SEQ+TROIS*XK*ALPHA) 1360 : & -DDL2*D11 1361 : SIGPP(2)=SIGPP(2)-DDL1*(G*SD(2)/SEQ+TROIS*XK*ALPHA) 1362 : & -DDL2*D12 1363 : SIGPP(3)=SIGPP(3)-DDL1*(G*SD(3)/SEQ+TROIS*XK*ALPHA) 1364 : & -DDL2*D12 1365 : 1366 : **-- calcul des nouveaux criteres phi1 et phi2 1367 : 1368 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 1369 : 1370 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 1371 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 1372 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 1373 : 1374 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 1375 : & +(SIGPP(3)-SIGPP(1))**2 1376 : SEQ=SQRT(SEQ/6.D0) 1377 : 1378 : DIFF=XP1-XU 1379 : IF (XP1.LT.XU) THEN 1380 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 1381 : ELSE IF ((XP1.GE.XU)) THEN 1382 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 1383 : ENDIF 1384 : PHI1=SEQ+ALPHA*XI1-ZETA 1385 : 1386 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 1387 : 1388 : END DO 1389 : **-- test de convergence 1390 : 1391 : IF((ABS(PHI1).LT.PREC0*PHI10).AND.(ABS(PHI2).LT.PREC0*PHI20) 1392 : & .OR.(ABS(PHI1).LT.PREC1).AND.(ABS(PHI2).LT.PREC1)) THEN 1393 : 1394 : IF(DL1.LT.XZER) THEN 1395 : * PRINT*,'1er multiplicateur negatif (D.P.-Rankine)' 1396 : ICOMPT=0 1397 : DL2=XZER 1398 : PHI2=F2INI 1399 : XP1=VAR0(1) 1400 : XP2=VAR0(2) 1401 : DO I=1,3 1402 : SIGPP(I)=SINI(I) 1403 : END DO 1404 : * debut appel bloc 20 1405 : *** debut bloc 20 1406 : DO WHILE ((ICOMPT .LT. 20) . AND. .NOT. 1407 : & ((ABS(PHI2).LT.PREC0*PHI20).OR.(ABS(PHI2).LT.PREC1))) 1408 : ICOMPT=ICOMPT+1 1409 : 1410 : **-- calcul du multiplicateur 1411 : 1412 : DDL2=PHI2/(D11-DF2DX(XP2,XKU,FT)) 1413 : XP2=XP2+DDL2 1414 : DL2=DL2+DDL2 1415 : IF(DL2.LT.XZER) THEN 1416 : PRINT*,'Multiplicateur negatif (Rankine 1 Critere)' 1417 : STOP 1418 : ENDIF 1419 : 1420 : SIGPP(1)=SIGPP(1)-DDL2*D11 1421 : SIGPP(2)=SIGPP(2)-DDL2*D12 1422 : SIGPP(3)=SIGPP(3)-DDL2*D12 1423 : 1424 : **-- calcul du nouveau critere phi2 1425 : 1426 : PHI2=SIGPP(1)-FT*TX2(XP2,XKU) 1427 : PHI2B=SIGPP(2)-FT*TX2(XP2,XKU) 1428 : END DO 1429 : **-- test de convergence 1430 : 1431 : IF((ABS(PHI2).LT.PREC0*PHI20).OR. 1432 : & (ABS(PHI2).LT.PREC1)) THEN 1433 : GO TO 998 1434 : ELSE 1435 : IF(ICOMPT.GE.20) THEN 1436 : PRINT*,'le nombre d iterations internes depasse 20 1437 : & (Rankine 1 Critere)' 1438 : 1439 : STOP 1440 : ENDIF 1441 : ENDIF 1442 : *** fin bloc 20 1443 : * fin appel bloc 20 1444 : ENDIF 1445 : 1446 : IF (DL2.LT.XZER) THEN 1447 : * PRINT*,'2eme multiplicateur negatif (D.P.-Rankine)' 1448 : ICOMPT=0 1449 : DL1=XZER 1450 : PHI1=F1INI 1451 : XP1=VAR0(1) 1452 : XP2=VAR0(2) 1453 : SEQ=SEQI 1454 : DO I=1,3 1455 : SIGPP(I)=SINI(I) 1456 : SD(I)=SDINI(I) 1457 : END DO 1458 : * debut appel bloc 9 1459 : ** début bloc 9 1460 : IF(XP1.GE.XU2) THEN 1461 : DO I=1,NSTRS 1462 : SIGF(I)=XZER 1463 : DEFP(I)=DEPST(I) 1464 : END DO 1465 : DDL1=PHI1/AUX3 1466 : VARF(1)=XP1+DDL1 1467 : VARF(2)=XP2 1468 : VARF(3)=VARF(1)+VARF(2) 1469 : VARF(4)=VARF(1)/XU 1470 : VARF(5)=EXP(-VARF(2)/XKU)/FT 1471 : GO TO 999 1472 : ENDIF 1473 : *** debut bloc 10 1474 : DO WHILE ((ICOMPT .LT. 20) .AND. .NOT. 1475 : & ((ABS(PHI1).LT.(PREC0*PHI10)).OR.(ABS(PHI1).LT.PREC1))) 1476 : ICOMPT=ICOMPT+1 1477 : 1478 : DIFF=XP1-XU 1479 : IF (XP1.LT.XU) THEN 1480 : DFDX=-DEUX*XK0*FC*OM1*DIFF/XUXU 1481 : ELSE IF((XP1.GE.XU)) THEN 1482 : DFDX=DEUX*XK0*FC*DIFF/((XU2-XU)**2) 1483 : ENDIF 1484 : 1485 : **-- calcul du multiplicateur 1486 : 1487 : DDL1=PHI1/(AUX3-DFDX) 1488 : XP1=XP1+DDL1 1489 : DL1=DL1+DDL1 1490 : 1491 : IF(XP1.GE.XU2) THEN 1492 : DO I=1,NSTRS 1493 : SIGF(I)=XZER 1494 : DEFP(I)=DEPST(I) 1495 : END DO 1496 : VARF(1)=XP1-DDL1+(PHI1/AUX3) 1497 : VARF(2)=XP2 1498 : VARF(3)=VARF(1)+VARF(2) 1499 : VARF(4)=VARF(1)/XU 1500 : VARF(5)=EXP(-VARF(2)/XKU)/FT 1501 : GO TO 999 1502 : ENDIF 1503 : 1504 : SIGPP(1)=SIGPP(1)-DDL1*(G*SD(1)/SEQ+TROIS*XK*ALPHA) 1505 : SIGPP(2)=SIGPP(2)-DDL1*(G*SD(2)/SEQ+TROIS*XK*ALPHA) 1506 : SIGPP(3)=SIGPP(3)-DDL1*(G*SD(3)/SEQ+TROIS*XK*ALPHA) 1507 : 1508 : **-- calcul du nouveau critere phi1 1509 : 1510 : XI1=SIGPP(1)+SIGPP(2)+SIGPP(3) 1511 : 1512 : SD(1)=(DEUX*SIGPP(1)-SIGPP(2)-SIGPP(3))/TROIS 1513 : SD(2)=(DEUX*SIGPP(2)-SIGPP(1)-SIGPP(3))/TROIS 1514 : SD(3)=(DEUX*SIGPP(3)-SIGPP(1)-SIGPP(2))/TROIS 1515 : 1516 : SEQ=(SIGPP(1)-SIGPP(2))**2+(SIGPP(2)-SIGPP(3))**2 1517 : & +(SIGPP(3)-SIGPP(1))**2 1518 : SEQ=SQRT(SEQ/6.D0) 1519 : 1520 : DIFF=XP1-XU 1521 : IF (XP1.LT.XU) THEN 1522 : ZETA=XK0*FC*(UN+OM1*(DIFF**2)/XUXU) 1523 : ELSE IF ((XP1.GE.XU)) THEN 1524 : ZETA=XK0*FC*( UN -((DIFF**2)/((XU2-XU)**2))) 1525 : ENDIF 1526 : PHI1=SEQ+ALPHA*XI1-ZETA 1527 : END DO 1528 : 1529 : **-- test de convergence 1530 : 1531 : IF((ABS(PHI1).LT.(PREC0*PHI10)).OR. 1532 : & (ABS(PHI1).LT.PREC1))THEN 1533 : GO TO 998 1534 : ELSE IF(ICOMPT.GE.20) THEN 1535 : PRINT*,'le nombre d iterations internes depasse 20 1536 : & (Drucker-Prager 1 Critere)' 1537 : STOP 1538 : ENDIF 1539 : *** fin bloc 10 1540 : *** fin bloc 9 1541 : * fin appel bloc 9 1542 : ENDIF 1543 : GO TO 998 1544 : ELSE IF(ICOMPT.GE.20) THEN 1545 : PRINT*,'le nombre d iterations internes depasse 20 1546 : & (D.P.-Rankine)' 1547 : STOP 1548 : ENDIF 1549 : * fin bloc 30 1550 : * fin bloc 29 1551 : * fin appel bloc 29 1552 : ENDIF 1553 : ENDIF 1554 : ENDIF 1555 : GO TO 998 1556 : ELSE 1557 : IF(ICOMPT.GE.20) THEN 1558 : PRINT*,'le nombre d iterations internes depasse 20 1559 : & (D.P.-Rankine)' 1560 : STOP 1561 : ENDIF 1562 : ENDIF 1563 : ** fin bloc 31 1564 : ENDIF 1565 : 1566 : * PRINT*,'pas dedans' 1567 : 998 DO I=1,3 1568 : DO J=1,3 1569 : XMS(I,J)=XZER 1570 : IF(I.EQ.J) XMS(I,J)=SIGPP(I) 1571 : ENDDO 1572 : ENDDO 1573 : * PRINT*,'dedans' 1574 : **-- calcul de sigf = R * sigpp * RT 1575 : 1575 : CALL TRSPOD (R,3,3,RT) 1575 : CALL MULMAT(TRAV,XMS,RT,3,3,3) 1575 : CALL MULMAT(XMS,R,TRAV,3,3,3) 1576 : 1577 : SIGF(1)=XMS(1,1) 1578 : SIGF(2)=XMS(2,2) 1579 : SIGF(3)=XMS(3,3) 1580 : SIGF(4)=XMS(1,2) 1581 : IF (IFOURB.EQ.2) THEN 1582 : SIGF(5)=XMS(1,3) 1583 : SIGF(6)=XMS(2,3) 1584 : ENDIF 1585 : DO I=1,NSTRS 1586 : DS(I)=S(I)-SIGF(I) 1587 : ENDDO 1587 : CALL EPSIG0(DS,DEFP,MFR,IFOURB,YOUN,XNU,XCAR,NSTRS) 1588 : VARF(1)=XP1 1589 : VARF(2)=XP2 1590 : VARF(3)=VARF(1)+VARF(2) 1591 : VARF(4)=VARF(1)/XU 1592 : VARF(5)=EXP(-VARF(2)/XKU)/FT 1593 : 999 RETURN 1594 : END 1595 : 1596 : FUNCTION DF2DX(XP2,XKU,FT) 1597 : IMPLICIT REAL*8(A-H,O-Z) 1598 : IF (XKU.EQ.0.D0) THEN 1599 : DF2DX=0.D0 1600 : ELSE 1601 : DF2DX=FT*EXP(-XP2/XKU)/XKU 1602 : ENDIF 1603 : END 1604 : 1605 : FUNCTION TX2(XP2,XKU) 1606 : IMPLICIT REAL*8(A-H,O-Z) 1607 : IF (XKU.EQ.0.D0) THEN 1608 : TX2=1.D0 1609 : ELSE 1610 : TX2=EXP(-XP2/XKU) 1611 : ENDIF 1612 : END 1613 : 1614 : 1615 : 1616 : 1617 : 1618 : 1619 :
© Cast3M 2003 - All rights reserved.
Disclaimer