C ECUBI2 SOURCE CHAT 05/01/12 23:28:05 5004 SUBROUTINE ECUBI2 (SOG1,DSOG1,GAMMA1,PHI1,PSI1,NCAS,TRA1, 1 XNU,YOUNG,COHE1,HACHE1,CO11,SI11, 2 PHI2,PSI2,COHE2,CO22,SI22,HACHE2,TRA2, 3 CO21,SI21,TSUG,DEFPL,DLAM,ICHR,SOGMA,DSOGMA,KERRE) C IMPLICIT INTEGER(I-N) IMPLICIT REAL*8(A-H,O-Z) -INC CCOPTIO C DIMENSION DSOG1(4),SOG1(4),DEFPL(4),SOGMA(4),DSOGMA(4), 1 SUG1(4),DSUG1(4),SIP1(4),SAP1(4),SYG1(4),DSAG1(4),SAG1(4), 2 TSUG2(4),SUG2(4),SIP2(4),SAP2(4),EPPLA1(4),TSUG(4),DSUG(4) 3 ,TSUG1(4),DSEG1(4),TSEG1(4),TSEG2(4),DSOG2(4),DSIG1(4) ICHR=0 ZER=0.D0 ITER=0 DLA=0.D0 C------------------------------------------------- C ON EFFECTUE L'ECOULEMENT1 SANS COUPLAGE C------------------------------------------------- C -------------------------------- DO 90 ITYP=1,4 C C ON SE PLACE SUR LE CRITERE 1 C SAG1(ITYP)=SOG1(ITYP)+GAMMA1*DSOG1(ITYP) DSAG1(ITYP)=(1-GAMMA1)*DSOG1(ITYP) 90 CONTINUE C -------------------------------- C ----------------------------------------------- C C CALCUL DE DLAMDA C CIS1=SIGN(1.D0,SAG1(4)) 457 VCRIT1=VCRITE(SAG1(2),SAG1(4),PHI1,COHE1) DLA=(PHI1*DSAG1(2)+CIS1 1 *DSAG1(4)+VCRIT1)/HACHE1 C---------------------------------------------- C--------------------------- C CAS DEFO PLANES C IF(NCAS.NE.1) CALL DPHOOK(ZER,PSI1*DLA,ZER,CIS1 1 *DLA,DSUG1(1),DSUG1(2),DSUG1(3),DSUG1(4),XNU,YOUNG) C--------------------------- C--------------------------- C CONT PLANES C IF(NCAS.EQ.1) CALL CPHOOK(ZER,PSI1*DLA,ZER,CIS1 1 *DLA,DSUG1(1),DSUG1(2),DSUG1(3),DSUG1(4),XNU,YOUNG) C--------------------------- C--------------------------------------------------- C ON EFFECTUE L ECOULEMENT C DO 91 ITYP=1,4 TSUG1(ITYP)=DSAG1(ITYP)+SAG1(ITYP)-DSUG1(ITYP) 91 CONTINUE C--------------------------------------------------- C ON REGARDE SI LE CRITERE 2 EST ENDOMMAGE C OU SI ON DEPASSE LA LIMITE EN TRACTION C CALL CHREPE(CO21,SI21,TSUG1,TSUG2) VCRI1=VCRITE(TSUG1(2),TSUG1(4),PHI1,COHE1) IF(IIMPI.EQ.28) WRITE(IOIMP,1978) VCRI1 VCRIT2=VCRITE(TSUG2(2),TSUG2(4),PHI2,COHE2) C IF(TSUG1(2).LE.TRA1.AND.VCRIT2.LE.0.D0) DLAM=DLAM+DLA IF(TSUG1(2).LE.TRA1.AND.VCRIT2.LE.0.D0) GO TO 500 C----------------------- C CAS DEFO PLANES C DLB=DLA IF(NCAS.NE.1) CALL DPHOOK(ZER,PSI1,ZER,CIS1 1 ,SYG1(1),SYG1(2),SYG1(3),SYG1(4),XNU,YOUNG) C----------------------- C CAS CONT PLANES C IF(NCAS.EQ.1) CALL CPHOOK(ZER,PSI1,ZER,CIS1 1 ,SYG1(1),SYG1(2),SYG1(3),SYG1(4),XNU,YOUNG) C----------------------- C--------------------------------------------------- C------------------------------------------------- C CAS OU ON DEPASSE LA LIMITE EN TRACTION C ON CALCULE UN NOUVEAU DLAMDA C IF(TSUG1(2).LE.TRA1) GO TO 700 XAT=(TRA1-SAG1(2)+VCRIT1*SYG1(2)/HACHE1)/(DSAG1(2)-SYG1(2)*(DLA- 1 VCRIT1/HACHE1)) DLA=XAT*(DLA-VCRIT1/HACHE1)+VCRIT1/HACHE1 DO 102 ITYP=1,4 TSUG1(ITYP)=SAG1(ITYP)+XAT*DSAG1(ITYP)-DLA*SYG1(ITYP) 102 CONTINUE IF(IIMPI.EQ.28) WRITE(IOIMP,3333) 3333 FORMAT('0 ON A DEPASSE LA LIMITE EN TRACTION') C------------------------------------------------- C------------------------------------------------------------- C LE CRITERE 2 EST ENDOMMAGE QUAND ON EST A LA C POINTE DU CRITERE 1 C ON CALCULE UN NOUVEAU DLAMDA C CALL CHREPE(CO21,SI21,TSUG1,TSUG2) 1978 FORMAT( '0 ECUBI4 VCRIT1 =',1PE12.5/) VCR=VCRITE(TSUG2(2),TSUG2(4),PHI2,COHE2) IF(VCR.LE.0) GO TO 750 700 VCRIT1=VCRITE(SAG1(2),SAG1(4),PHI1,COHE1) IF(IIMPI.EQ.28) WRITE(IOIMP,1978) VCRIT1 DLA=(PHI1*DSAG1(2)+CIS1 1 *DSAG1(4)+VCRIT1)/HACHE1 DO 101 ITYP=1,4 SIP1(ITYP)=SAG1(ITYP)-VCRIT1/HACHE1*SYG1(ITYP) SAP1(ITYP)=DSAG1(ITYP)-(DLA-VCRIT1/HACHE1)*SYG1(ITYP) 101 CONTINUE IF(IIMPI.EQ.28) WRITE(IOIMP,4333) 4333 FORMAT('0 ON A ENDOMMAGE LE 2 CRITERE') CALL CHREPE(CO21,SI21,SIP1,SIP2) CALL CHREPE(CO21,SI21,SAP1,SAP2) XBT=-(PHI2*SIP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/( 1 PHI2*SAP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) XCT=-(PHI2*SIP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/( 1 PHI2*SAP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) XAT=MIN(XBT,XCT) IF(XBT.LE.0.D0) XAT=XCT IF(XCT.LE.0.D0) XAT=XBT IF(IIMPI.EQ.28) WRITE(IOIMP,3888) DLA,XAT,XBT,XCT 3888 FORMAT( '0 ECUBI4 DLA XAT XBT XCT ',4(1X,1PE12.5)/) DLA=XAT*(DLA-VCRIT1/HACHE1)+VCRIT1/HACHE1 IF(IIMPI.EQ.28) WRITE(IOIMP,3889) DLA 3889 FORMAT('0 ECUBI4 DLA=',1PE12.5/) DO 105 ITYP=1,4 SUG1(ITYP)=SAG1(ITYP)+XAT*DSAG1(ITYP)-DLA*SYG1(ITYP) SAG1(ITYP)=SAG1(ITYP)+XAT*DSAG1(ITYP)-SUG1(ITYP) DSUG1(ITYP)=DSAG1(ITYP)*(1.D0-XAT) 105 CONTINUE C------------------------- C IL FAUT APPELLE LE COUPLAGE C------------------------- CALL CHREPE(CO11,-SI11,SUG1,SOGMA) VCRIT1=VCRITE(SUG1(2),SUG1(4),PHI1,COHE1) IF(IIMPI.EQ.28) WRITE(IOIMP,1978) VCRIT1 CALL CHREPE(CO11,-SI11,DSUG1,DSOGMA) CALL CHREPE(CO11,-SI11,SAG1,DSAG1) IF (NCAS.NE.1) CALL DPCONT(DSAG1,DEFPL,XNU,YOUNG) IF (NCAS.EQ.1) CALL CPCONT(DSAG1,DEFPL,XNU,YOUNG) ICHR=1 DLAM=DLAM+DLA IF(IIMPI.EQ.28) WRITE(IOIMP,1111) ITER 1111 FORMAT('0 ITER= ',E12.5/) RETURN C------------------------------------------------------------- 750 DO 103 ITYP=1,4 DSIG1(ITYP)=DSAG1(ITYP)*(1.D0-XAT) 103 CONTINUE VCRT=-DSIG1(2)+PSI1*ABS(DSIG1(4)) IF(VCRT.LE.0.D0.OR.ITER.LE.1) GO TO 456 ITER=ITER+1 CIS1=SIGN(1.D0,DSIG1(4)) DO 97 ITYP=1,4 SAG1(ITYP)=TSUG1(ITYP) DSAG1(ITYP)=DSIG1(ITYP) 97 CONTINUE DLAM=DLAM+DLA GO TO 457 456 IF (NCAS.NE.1) HOOK=YOUNG/(1.D0+XNU)/(1.D0-2.D0*XNU) IF (NCAS.NE.1) ALPHA=(DSIG1(2)+TSUG1(2)-TRA1)/HOOK /(1.D0-XNU) IF (NCAS.NE.1) BETA=(DSIG1(4)+TSUG1(4))/HOOK/(1.D0-2.D0*XNU) C-------------------------- C CAS CONT PLANES C IF (NCAS.EQ.1) HOOK=YOUNG/(1.D0-XNU*XNU) IF (NCAS.EQ.1) ALPHA=(DSIG1(2)+TSUG1(2)-TRA1)/HOOK IF (NCAS.EQ.1) BETA=(DSIG1(4)+TSUG1(4))/HOOK/(1.D0-XNU) C----------------------------------------------------- IF(NCAS.EQ.1) CALL CPHOOK(ZER,ALPHA,ZER, # BETA,DSEG1(1),DSEG1(2) # ,DSEG1(3),DSEG1(4),XNU,YOUNG) IF(NCAS.NE.1) CALL DPHOOK(ZER,ALPHA,ZER, # BETA,DSEG1(1),DSEG1(2) # ,DSEG1(3),DSEG1(4),XNU,YOUNG) DO 453 ITYP=1,4 TSEG1(ITYP)=TSUG1(ITYP)+DSIG1(ITYP)-DSEG1(ITYP) DSUG1(ITYP)=DSOG1(ITYP)+SOG1(ITYP)-TSEG1(ITYP) 453 CONTINUE C------------------------------------- C ON REGARDE SI ON ENDOMMAGE LE 2 EME CRITERE C------------------------------------- C CALL CHREPE(CO21,SI21,TSEG1,TSUG2) VCRTA=VCRITE(TSUG2(2),TSUG2(4),PHI2,COHE2) IF(IIMPI.EQ.28) WRITE(IOIMP,2978) VCRTA 2978 FORMAT( '0 ECUBI VCRTA =',1PE12.5/) IF(VCRTA.LE.0.) GO TO 522 C CALL CHREPE(CO21,SI21,DSUG1,DSUG2) C IF(ABS(DSEG1(4)).LE.1.D-7) ATS=0.D0 IF(ABS(DSEG1(4)).LE.1.D-7) BTS=1.D0 IF(ABS(DSEG1(4)).LE.1.D-7) GOTO 543 C ATS=(TSUG1(4)-TRA1)/DSEG1(4) ATS=TSUG1(4)/DSEG1(4) C BTS=DSIG1(2)/DSEG1(2) BTS=DSIG1(4)/DSEG1(4) 543 DO 987 ITYP=1,4 SIP1(ITYP)=TSUG1(ITYP)-ATS*DSEG1(ITYP) SAP1(ITYP)=DSIG1(ITYP)-BTS*DSEG1(ITYP) 987 CONTINUE CALL CHREPE(CO21,SI21,SIP1,SIP2) CALL CHREPE(CO21,SI21,SAP1,SAP2) XCT=-(PHI2*SIP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/ 1 (PHI2*SAP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) XDT=-(PHI2*SIP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/ 1 (PHI2*SAP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) IF(XCT.LE.0.D0) XCT=100.D0 IF(XDT.LE.0.D0) XDT=100.D0 XAT=MIN(XCT,XDT) IF(XAT.GT.1.) XAT=0.D0 XBT=ATS+XAT*BTS IF(IIMPI.EQ.28) WRITE(IOIMP,3978) XAT,XBT,XCT,XDT 3978 FORMAT( '0 ECUBI XAT XBT XCT XDT ',4(1X,1PE12.5)/) DO 195 ITYP=1,4 TSUG1(ITYP)=TSUG1(ITYP)+XAT*DSIG1(ITYP)-DSEG1(ITYP)*XBT DSUG1(ITYP)=DSOG1(ITYP)+SOG1(ITYP)-TSUG1(ITYP) DSOG1(ITYP)=(1.D0-XAT)*DSIG1(ITYP) 195 CONTINUE CALL CHREPE(CO11,-SI11,TSUG1,SOGMA) CALL CHREPE(CO11,-SI11,DSOG1,DSOGMA) CALL CHREPE(CO11,-SI11,DSUG1,DSUG) ICHR=1 IF (NCAS.NE.1) CALL DPCONT(DSUG,DEFPL,XNU,YOUNG) IF (NCAS.EQ.1) CALL CPCONT(DSUG,DEFPL,XNU,YOUNG) DLAM=DLAM+DLA*XBT IF(IIMPI.EQ.28) WRITE(IOIMP,1111) ITER RETURN 522 CALL CHREPE(CO11,-SI11,TSEG1,TSUG) DO 32 ITYP=1,4 TSUG1(ITYP)=TSEG1(ITYP) 32 CONTINUE DLAM=DLA+SQRT((ALPHA*ALPHA+BETA*BETA)/(PSI1*PSI1+1.D0))+DLAM GOTO 523 500 CALL CHREPE(CO11,-SI11,TSUG1,TSUG) 523 CONTINUE DO 31 ITYP=1,4 DSUG1(ITYP)=-TSUG1(ITYP)+DSOG1(ITYP)+SOG1(ITYP) 31 CONTINUE CALL CHREPE(CO11,-SI11,DSUG1,DSUG) IF (NCAS.NE.1) CALL DPCONT(DSUG,DEFPL,XNU,YOUNG) IF (NCAS.EQ.1) CALL CPCONT(DSUG,DEFPL,XNU,YOUNG) IF(IIMPI.EQ.28) WRITE(IOIMP,1111) ITER RETURN END