C T1C2DP SOURCE CHAT 05/01/13 03:29:16 5004 SUBROUTINE T1C2DP(SIG,DSIG,YOUN,ANU,RT1,RT2,H1,RDP,ADP,HDP, 1 XLAM1,XLAM2,XLAM3,IDAM,ANG,EPPLDP,ALPHA,KERRE) IMPLICIT INTEGER(I-N) IMPLICIT REAL*8(A-H,O-Z) -INC CCOPTIO DIMENSION SIG(3),DSIG(3),SGG(3),DSGG(3),SFG(3),DSFG(3), 1 A(3),B(3),IDAM(3),EPC(3),SSI(3),SI(3),EPC1(3),SG1(3), 1 SG2(3),EC1(3),EC2(3),AT(4),BT(2),EPDP(3),EPPLDP(3),TENS(3) ITER=0 Y=YOUN/(1.D0-ANU*ANU) DO 5 ITYP=1,3 EPPLDP(ITYP)=0.D0 5 EPDP(ITYP)=0.D0 XT=0.D0 IF(IIMPI.EQ.9) WRITE(IOIMP,9999) 9999 FORMAT(1X,'T1C2DP TRIPLAGE TRACTION 1 COMP 2 DRUCKER',/) 9998 FORMAT(1X,'T1C2DP ',I4,'ITERATIONS INTERNES',/) C C------------------------------------------------ C TRIPLAGE TRACTION 1 COMPRESSION 2 DRUCKER C------------------------------------------------ C---------------------------------------------- C DEFINITION D UNE CONTRAINTE DE REFERENCE C SI RMAX DU DRUCKER INFERIEURE A CETTE VALEUR C ALORS SIGMA=0 C---------------------------------------------- C RMAX=MAX((RDP/1.73),(RDP/(1.D0-2.D0*ADP))) SIREF=1.E-6*YOUN C C-------------------------------------------------- C CAS OU LE RAYON DU DRUCKER EST INFERIEUR C A UNE VALEUR DE REFERENCE C-------------------------------------------------- C IF(SIREF.GT.RMAX) THEN DO 10 ITYP=1,3 10 SIG(ITYP)=SIG(ITYP)+DSIG(ITYP) RDP=0.D0 HDP=0.D0 ADP=0.D0 CALL NORME(SIG,DL3) DL3=DL3/YOUN*10.D0 XLAM3=XLAM3+DL3 CALL CPHOIN(SIG,EPPLDP,YOUN,ANU) DO 20 ITYP=1,3 IDAM(ITYP)=0 DSIG(ITYP)=0.D0 20 SIG(ITYP)=0.D0 IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER RETURN ENDIF C----------------------------------------------- C ON SE PLACE DANS LE REPERE DE FISSURATION C----------------------------------------------- C CALL CHREP(SIG,SFG,ANG) CALL CHREP(DSIG,DSFG,ANG) C C----------------------------------------------------- C ESTIMATION DU PAS D'INCREMENT DE CONTRAINTES C----------------------------------------------------- C CALL SIJ(SFG,SI,SIETOI) CALL EPSDP(SI,SIETOI,ADP,EPC) CALL CPHOMO(EPC,SSI,YOUN,ANU,ALPHA) CALL SCAL2(EPC,SSI,VAL) VAL=VAL-HDP IF(VAL.LT.0.D0) THEN KERRE=459 RETURN ENDIF CALL SCAL(SSI,DSFG,VAL) CALL NORME(DSFG,VA1) CALL NORME(SSI,VA2) IF(VA1.EQ.0.D0) THEN IDAM(1)=0 IDAM(2)=0 IDAM(3)=0 RETURN ENDIF RMIN=MIN((RDP/1.73),(RDP/(1.D0+2.D0*ADP))) X=VAL/VA1/VA2 IF(X.GE.1.D0) THEN X=1.D0 GOTO 25 ENDIF X=1.D0/SQRT(1.0001D0-X*X)*RMIN/VA1/8.D0 IF(X.GT.1.D0) X=1.D0 C C---------------------------------------------------------- C ON ECOULE ET ON REGARDE LA VARIATION DE NORMALE C---------------------------------------------------------- C C C--------------------------------------- C RESOLUTION DU SYSTEME EN DL3 C--------------------------------------- 25 X=X/2.D0 EC1(1)=1.D0 EC1(2)=0.D0 EC1(3)=0.D0 EC2(1)=0.D0 EC2(2)=1.D0 EC2(3)=0.D0 CALL CPHOOB(EC1,SG1,YOUN,ANU) CALL CPHOOB(EC2,SG2,YOUN,ANU) 51 X=X*2.D0 52 RMAX=MAX((RDP/1.73),(RDP/(1.D0-2.D0*ADP))) ITER=ITER+1 653 IF(ITER.GT.200) THEN KERRE=460 RETURN ENDIF IF(SIREF.GT.RMAX) THEN RDP=0.D0 ADP=0.D0 HDP=0.D0 DO 60 ITYP=1,3 60 SFG(ITYP)=SFG(ITYP)+(1.D0-XT)*DSFG(ITYP) CALL NORME(SFG,DL3) DL3=DL3/YOUN*10.D0 XLAM3=XLAM3+DL3 CALL CPHOIN(SFG,EPPLDP,YOUN,ANU) DO 70 ITYP=1,3 70 EPDP(ITYP)=EPDP(ITYP)+EPPLDP(ITYP) CALL CHREP(EPDP,EPPLDP,-ANG) DO 80 ITYP=1,3 IDAM(ITYP)=0 DSIG(ITYP)=0.D0 80 SIG(ITYP)=0.D0 IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER RETURN ENDIF IF((XT+X).GT.1.D0) X=1.D0-XT CALL SIJ(SFG,SI,SIETOI) CALL EPSDP(SI,SIETOI,ADP,EPC) CALL CPHOMO(EPC,SSI,YOUN,ANU,ALPHA) CALL SCAL2(EPC,SSI,VAL) VAL=VAL-HDP IF(VAL.LT.0.D0) THEN KERRE=459 RETURN ENDIF C 90 ITER=ITER+1 IF(ITER.GT.201) GOTO 653 IDAM(2)=-1 AT(1)=Y-H1 AT(4)=Y AT(2)=Y*ANU AT(3)=AT(2) DO 95 ITYP=1,3 95 DSGG(ITYP)=X*DSFG(ITYP) BT(1)=SFG(1)+DSGG(1)-RT1 BT(2)=SFG(2)+DSGG(2) CALL SYLIN2(AT,BT,DL11,DL21) BT(1)=-SSI(1) BT(2)=-SSI(2) CALL SYLIN2(AT,BT,DL12,DL22) DO 100 ITYP=1,3 A(ITYP)=SFG(ITYP)+DSGG(ITYP)-DL11*SG1(ITYP)-DL21*SG2(ITYP) 100 B(ITYP)=-DL12*SG1(ITYP)-DL22*SG2(ITYP)-SSI(ITYP) CALL DLAMDP(A,B,DL3,RDP,ADP,HDP,ITEST) IF(ITEST.EQ.1) THEN X=X/2.D0 GOTO 90 ENDIF DL1=DL11+DL3*DL12 DL2=DL21+DL3*DL22 DO 105 ITYP=1,3 SGG(ITYP)=SFG(ITYP)-DL1*SG1(ITYP)-DL2*SG2(ITYP) 105 SGG(ITYP)=SGG(ITYP)-DL3*SSI(ITYP)+X*DSFG(ITYP) IF(DL1.LT.-1.D-10) IDAM(1)=0 IF(DL2.GT.1.D-10) IDAM(2)=0 IF(DL3.LT.-1.D-10) IDAM(3)=0 IF(IDAM(1).EQ.0.OR.IDAM(2).EQ.0.OR.IDAM(3).EQ.0) THEN DL1=0.D0 DL2=0.D0 DL3=0.D0 X=0.D0 GOTO 2000 ENDIF C C------------------------------------------------------- C ON VERIFIE SI LA NORMALE DU DRUCKER VARIE PEU C------------------------------------------------------- C CALL SIJ(SGG,SI,SIETOI) CALL EPSDP(SI,SIETOI,ADP,EPC1) CALL SCAL(EPC,EPC1,VAL) CALL NORME(EPC,VA1) CALL NORME(EPC1,VA2) CO=VAL/VA1/VA2 IF(CO.LT.0.99) THEN X=X/2.D0 GOTO 90 ENDIF RD=RDP-DL3*HDP RMAX=MAX((RD/1.73),(RD/(1.D0-2.D0*ADP))) IF(RMAX.LT.SIREF) THEN RDP=0.D0 GOTO 52 ENDIF C C------------------------------------------------- C CAS OU LA TRACTION PASSE EN DESSOUS DE 0 C------------------------------------------------- C IF((RT1-H1*DL1).GE.-1.E-10) GOTO 200 C IF(H1.EQ.0.D0) H1=1.D0 DL1=RT1/H1 AT(1)=-DSFG(1) AT(3)=-DSFG(2) AT(4)=Y AT(2)=Y*ANU BT(1)=SFG(1)-DL1*SG1(1) BT(2)=SFG(2)-DL1*SG1(2) CALL SYLIN2(AT,BT,X1,DL21) BT(1)=-SSI(1) BT(2)=-SSI(2) CALL SYLIN2(AT,BT,X2,DL22) DO 110 ITYP=1,3 A(ITYP)=SFG(ITYP)+X1*DSFG(ITYP)-DL1*SG1(ITYP)-DL21*SG2(ITYP) 110 B(ITYP)=X2*DSFG(ITYP)-DL22*SG2(ITYP)-SSI(ITYP) CALL DLAMD(A,B,DL31,DL32,RDP,ADP,HDP,ITEST) DO 111 ITYP=1,3 SGG(ITYP)=A(ITYP)+DL31*B(ITYP) 111 DSGG(ITYP)=A(ITYP)+DL32*B(ITYP) CALL SIJ(SGG,SI,SIETOI) CALL EPSDP(SI,SIETOI,ADP,EPC1) CALL SCAL(EPC,EPC1,VAL) CALL NORME(EPC,VA1) CALL NORME(EPC1,VA2) CO1=VAL/VA1/VA2 CALL SIJ(DSGG,SI,SIETOI) CALL EPSDP(SI,SIETOI,ADP,EPC1) CALL SCAL(EPC,EPC1,VAL) CALL NORME(EPC,VA1) CALL NORME(EPC1,VA2) CO2=VAL/VA1/VA2 XIN1=X1+X2*DL31 XIN2=X1+X2*DL32 DLIN1=DL21+DL31*DL22 DLIN2=DL21+DL32*DL22 IF(XIN1.GT.-1.E-10.AND.DLIN1.LT.1.E-10.AND.XIN1.LE.X 1 .AND.CO1.GT.0.9) THEN DL3=DL31 X=XIN1 DL2=DLIN1 GOTO 200 ENDIF IF(XIN2.GT.-1.E-10.AND.DLIN2.LT.1.E-10.AND.XIN2.LE.X 1 .AND.CO2.GT.0.9) THEN DL3=DL32 X=XIN2 DL2=DLIN2 GOTO 200 ENDIF IF(IIMPI.EQ.9) WRITE(IOIMP,10101) 10101 FORMAT(1X,'ERREUR DANS T1C2DP TRAC',/) DL3=DL31 X=XIN1 DL2=DLIN1 200 IF((XLAM2+DL2).GE.0.D0) GOTO 300 C C------------------------------------------------- C CAS OU LA FISSURE 2 EST TOTALEMENT FERME C------------------------------------------------- C C IDAM(2)=0 DL2=-XLAM2 AT(1)=-DSFG(1) AT(3)=-DSFG(2) AT(2)=Y-H1 AT(4)=Y*ANU BT(1)=SFG(1)-DL2*SG2(1)-RT1 BT(2)=SFG(2)-DL2*SG2(2) CALL SYLIN2(AT,BT,X1,DL11) BT(1)=-SSI(1) BT(2)=-SSI(2) CALL SYLIN2(AT,BT,X2,DL12) DO 210 ITYP=1,3 A(ITYP)=SFG(ITYP)+X1*DSFG(ITYP)-DL11*SG1(ITYP)-DL2*SG2(ITYP) 210 B(ITYP)=X2*DSFG(ITYP)-DL12*SG1(ITYP)-SSI(ITYP) CALL DLAMD(A,B,DL31,DL32,RDP,ADP,HDP,ITEST) DO 211 ITYP=1,3 SGG(ITYP)=A(ITYP)+DL31*B(ITYP) 211 DSGG(ITYP)=A(ITYP)+DL32*B(ITYP) CALL SIJ(SGG,SI,SIETOI) CALL EPSDP(SI,SIETOI,ADP,EPC1) CALL SCAL(EPC,EPC1,VAL) CALL NORME(EPC,VA1) CALL NORME(EPC1,VA2) CO1=VAL/VA1/VA2 CALL SIJ(DSGG,SI,SIETOI) CALL EPSDP(SI,SIETOI,ADP,EPC1) CALL SCAL(EPC,EPC1,VAL) CALL NORME(EPC,VA1) CALL NORME(EPC1,VA2) CO2=VAL/VA1/VA2 XIN1=X1+X2*DL31 XIN2=X1+X2*DL32 DLIN1=DL11+DL31*DL12 DLIN2=DL11+DL32*DL12 IF(XIN1.GT.-1.E-10.AND.DLIN1.GT.-1.E-10.AND.XIN1.LE.X 1 .AND.CO1.GT.0.9) THEN DL3=DL31 X=XIN1 DL1=DLIN1 GOTO 300 ENDIF IF(XIN2.GT.-1.E-10.AND.DLIN2.GT.-1.E-10.AND.XIN2.LE.X 1 .AND.CO2.GT.0.9) THEN DL3=DL32 X=XIN2 DL1=DLIN2 GOTO 300 ENDIF IF(IIMPI.EQ.9) WRITE(IOIMP,20202) 20202 FORMAT(1X,'ERREUR DANS T1C2DP FERM',/) DL3=DL31 X=XIN1 DL1=DLIN1 300 DO 310 ITYP=1,3 SGG(ITYP)=SFG(ITYP)-DL1*SG1(ITYP)-DL2*SG2(ITYP) 310 SGG(ITYP)=SGG(ITYP)-DL3*SSI(ITYP)+X*DSFG(ITYP) RDP=RDP-DL3*HDP RMAX=MAX((RDP/1.73),(RDP/(1.D0-2.D0*ADP))) IF(RMAX.LT.SIREF) THEN RDP=0.D0 GOTO 52 ENDIF XT=XT+X DO 330 ITYP=1,3 EPDP(ITYP)=EPDP(ITYP)+DL3*EPC(ITYP) 330 SFG(ITYP)=SGG(ITYP) XLAM1=XLAM1+DL1 XLAM2=XLAM2+DL2 XLAM3=XLAM3+DL3 RT1=RT1-H1*DL1 IF(RT1.LT.1.E-10) THEN H1=0.D0 RT1=0.D0 ENDIF IF(XLAM2.LT.1.E-10) THEN IDAM(2)=0 XLAM2=0.D0 ENDIF IF(IDAM(2).EQ.0) THEN XLAM2=0.D0 C DO 340 ITYP=1,3 DSFG(ITYP)=(1.D0-XT)*DSFG(ITYP) 340 TENS(ITYP)=SFG(ITYP)+DSFG(ITYP) GAMTR1=10.D0 GAMTR2=10.D0 GAMCO1=10.D0 C GAMDP=10.D0 CALL CTRAF(TENS(1),RT1,VCTR1) CALL CTRAF(TENS(2),RT2,VCTR2) CALL CCOAF(TENS(1),XLAM1,VCCO1) C CALL CDP(TENS,ADP,RDP,VCDP) IF(VCTR1.GT.0.D0) CALL GAMTAF(SFG(1),DSFG(1),RT1,GAMTR1) IF(VCTR2.GT.0.D0) CALL GAMTAF(SFG(2),DSFG(2),RT2,GAMTR2) IF(VCCO1.GT.0.D0) CALL GAMCAF(SFG(1),DSFG(1),GAMCO1) C IF(VCDP.GT.0.D0) CALL GDP(SFG,DSFG,RDP,ADP,GAMDP) IDAM(1)=0 IDAM(2)=0 IDAM(3)=0 GAM=MIN(GAMTR1,GAMTR2,GAMCO1,GAMDP) IF(GAM.GE.1.D0) THEN DO 341 ITYP=1,3 SFG(ITYP)=SFG(ITYP)+DSFG(ITYP) 341 DSIG(ITYP)=0.D0 CALL CHREP(SFG,SIG,-ANG) CALL CHREP(EPDP,EPPLDP,-ANG) IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER RETURN ENDIF IF(ABS(GAM-GAMTR1).LE.1.E-10) IDAM(1)=1 IF(ABS(GAM-GAMCO1).LE.1.E-10) IDAM(1)=-1 IF(ABS(GAM-GAMTR2).LE.1.E-10) IDAM(2)=1 C IF(ABS(GAM-GAMDP).LE.1.E-10) IDAM(3)=1 DO 342 ITYP=1,3 SFG(ITYP)=SFG(ITYP)-GAM*DSFG(ITYP) 342 DSFG(ITYP)=(1.D0-GAM)*DSFG(ITYP) CALL CHREP(SFG,SIG,-ANG) CALL CHREP(DSFG,DSIG,-ANG) CALL CHREP(EPDP,EPPLDP,-ANG) IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER RETURN ENDIF C=1.D0-1.D-10 IF(XT.GE.C) THEN DO 350 ITYP=1,3 IDAM(ITYP)=0 350 DSIG(ITYP)=0.D0 CALL CHREP(SFG,SIG,-ANG) CALL CHREP(EPDP,EPPLDP,-ANG) IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER RETURN ENDIF GOTO 51 2000 CONTINUE CALL CHREP(SFG,SIG,-ANG) DO 2010 ITYP=1,3 2010 DSFG(ITYP)=(1.D0-XT)*DSFG(ITYP) CALL CHREP(DSFG,DSIG,-ANG) CALL CHREP(EPDP,EPPLDP,-ANG) IF(IIMPI.EQ.9) WRITE(IOIMP,9998) ITER RETURN END