1 : C ECUBI2 SOURCE CHAT 05/01/12 23:28:05 5004 2 : SUBROUTINE ECUBI2 (SOG1,DSOG1,GAMMA1,PHI1,PSI1,NCAS,TRA1, 3 : 1 XNU,YOUNG,COHE1,HACHE1,CO11,SI11, 4 : 2 PHI2,PSI2,COHE2,CO22,SI22,HACHE2,TRA2, 5 : 3 CO21,SI21,TSUG,DEFPL,DLAM,ICHR,SOGMA,DSOGMA,KERRE) 6 : C 7 : IMPLICIT INTEGER(I-N) 8 : IMPLICIT REAL*8(A-H,O-Z) 9 : -INC CCOPTIO 10 : C 11 : DIMENSION DSOG1(4),SOG1(4),DEFPL(4),SOGMA(4),DSOGMA(4), 12 : 1 SUG1(4),DSUG1(4),SIP1(4),SAP1(4),SYG1(4),DSAG1(4),SAG1(4), 13 : 2 TSUG2(4),SUG2(4),SIP2(4),SAP2(4),EPPLA1(4),TSUG(4),DSUG(4) 14 : 3 ,TSUG1(4),DSEG1(4),TSEG1(4),TSEG2(4),DSOG2(4),DSIG1(4) 15 : ICHR=0 16 : ZER=0.D0 17 : ITER=0 18 : DLA=0.D0 19 : C------------------------------------------------- 20 : C ON EFFECTUE L'ECOULEMENT1 SANS COUPLAGE 21 : C------------------------------------------------- 22 : C -------------------------------- 23 : DO 90 ITYP=1,4 24 : C 25 : C ON SE PLACE SUR LE CRITERE 1 26 : C 27 : SAG1(ITYP)=SOG1(ITYP)+GAMMA1*DSOG1(ITYP) 28 : DSAG1(ITYP)=(1-GAMMA1)*DSOG1(ITYP) 29 : 90 CONTINUE 30 : C -------------------------------- 31 : C ----------------------------------------------- 32 : C 33 : C CALCUL DE DLAMDA 34 : C 35 : CIS1=SIGN(1.D0,SAG1(4)) 36 : 457 VCRIT1=VCRITE(SAG1(2),SAG1(4),PHI1,COHE1) 37 : DLA=(PHI1*DSAG1(2)+CIS1 38 : 1 *DSAG1(4)+VCRIT1)/HACHE1 39 : C---------------------------------------------- 40 : C--------------------------- 41 : C CAS DEFO PLANES 42 : C 42 : IF(NCAS.NE.1) CALL DPHOOK(ZER,PSI1*DLA,ZER,CIS1 43 : 1 *DLA,DSUG1(1),DSUG1(2),DSUG1(3),DSUG1(4),XNU,YOUNG) 44 : C--------------------------- 45 : C--------------------------- 46 : C CONT PLANES 47 : C 47 : IF(NCAS.EQ.1) CALL CPHOOK(ZER,PSI1*DLA,ZER,CIS1 48 : 1 *DLA,DSUG1(1),DSUG1(2),DSUG1(3),DSUG1(4),XNU,YOUNG) 49 : C--------------------------- 50 : C--------------------------------------------------- 51 : C ON EFFECTUE L ECOULEMENT 52 : C 53 : DO 91 ITYP=1,4 54 : TSUG1(ITYP)=DSAG1(ITYP)+SAG1(ITYP)-DSUG1(ITYP) 55 : 91 CONTINUE 56 : C--------------------------------------------------- 57 : C ON REGARDE SI LE CRITERE 2 EST ENDOMMAGE 58 : C OU SI ON DEPASSE LA LIMITE EN TRACTION 59 : C 59 : CALL CHREPE(CO21,SI21,TSUG1,TSUG2) 60 : VCRI1=VCRITE(TSUG1(2),TSUG1(4),PHI1,COHE1) 61 : IF(IIMPI.EQ.28) WRITE(IOIMP,1978) VCRI1 62 : VCRIT2=VCRITE(TSUG2(2),TSUG2(4),PHI2,COHE2) 63 : C 64 : IF(TSUG1(2).LE.TRA1.AND.VCRIT2.LE.0.D0) DLAM=DLAM+DLA 65 : IF(TSUG1(2).LE.TRA1.AND.VCRIT2.LE.0.D0) GO TO 500 66 : C----------------------- 67 : C CAS DEFO PLANES 68 : C 69 : DLB=DLA 69 : IF(NCAS.NE.1) CALL DPHOOK(ZER,PSI1,ZER,CIS1 70 : 1 ,SYG1(1),SYG1(2),SYG1(3),SYG1(4),XNU,YOUNG) 71 : C----------------------- 72 : C CAS CONT PLANES 73 : C 73 : IF(NCAS.EQ.1) CALL CPHOOK(ZER,PSI1,ZER,CIS1 74 : 1 ,SYG1(1),SYG1(2),SYG1(3),SYG1(4),XNU,YOUNG) 75 : C----------------------- 76 : C--------------------------------------------------- 77 : C------------------------------------------------- 78 : C CAS OU ON DEPASSE LA LIMITE EN TRACTION 79 : C ON CALCULE UN NOUVEAU DLAMDA 80 : C 81 : IF(TSUG1(2).LE.TRA1) GO TO 700 82 : XAT=(TRA1-SAG1(2)+VCRIT1*SYG1(2)/HACHE1)/(DSAG1(2)-SYG1(2)*(DLA- 83 : 1 VCRIT1/HACHE1)) 84 : DLA=XAT*(DLA-VCRIT1/HACHE1)+VCRIT1/HACHE1 85 : DO 102 ITYP=1,4 86 : TSUG1(ITYP)=SAG1(ITYP)+XAT*DSAG1(ITYP)-DLA*SYG1(ITYP) 87 : 102 CONTINUE 88 : IF(IIMPI.EQ.28) WRITE(IOIMP,3333) 89 : 3333 FORMAT('0 ON A DEPASSE LA LIMITE EN TRACTION') 90 : C------------------------------------------------- 91 : C------------------------------------------------------------- 92 : C LE CRITERE 2 EST ENDOMMAGE QUAND ON EST A LA 93 : C POINTE DU CRITERE 1 94 : C ON CALCULE UN NOUVEAU DLAMDA 95 : C 95 : CALL CHREPE(CO21,SI21,TSUG1,TSUG2) 96 : 1978 FORMAT( '0 ECUBI4 VCRIT1 =',1PE12.5/) 97 : VCR=VCRITE(TSUG2(2),TSUG2(4),PHI2,COHE2) 98 : IF(VCR.LE.0) GO TO 750 99 : 700 VCRIT1=VCRITE(SAG1(2),SAG1(4),PHI1,COHE1) 100 : IF(IIMPI.EQ.28) WRITE(IOIMP,1978) VCRIT1 101 : DLA=(PHI1*DSAG1(2)+CIS1 102 : 1 *DSAG1(4)+VCRIT1)/HACHE1 103 : DO 101 ITYP=1,4 104 : SIP1(ITYP)=SAG1(ITYP)-VCRIT1/HACHE1*SYG1(ITYP) 105 : SAP1(ITYP)=DSAG1(ITYP)-(DLA-VCRIT1/HACHE1)*SYG1(ITYP) 106 : 101 CONTINUE 107 : IF(IIMPI.EQ.28) WRITE(IOIMP,4333) 108 : 4333 FORMAT('0 ON A ENDOMMAGE LE 2 CRITERE') 108 : CALL CHREPE(CO21,SI21,SIP1,SIP2) 108 : CALL CHREPE(CO21,SI21,SAP1,SAP2) 109 : XBT=-(PHI2*SIP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/( 110 : 1 PHI2*SAP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) 111 : XCT=-(PHI2*SIP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/( 112 : 1 PHI2*SAP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) 113 : XAT=MIN(XBT,XCT) 114 : IF(XBT.LE.0.D0) XAT=XCT 115 : IF(XCT.LE.0.D0) XAT=XBT 116 : IF(IIMPI.EQ.28) WRITE(IOIMP,3888) DLA,XAT,XBT,XCT 117 : 3888 FORMAT( '0 ECUBI4 DLA XAT XBT XCT ',4(1X,1PE12.5)/) 118 : DLA=XAT*(DLA-VCRIT1/HACHE1)+VCRIT1/HACHE1 119 : IF(IIMPI.EQ.28) WRITE(IOIMP,3889) DLA 120 : 3889 FORMAT('0 ECUBI4 DLA=',1PE12.5/) 121 : DO 105 ITYP=1,4 122 : SUG1(ITYP)=SAG1(ITYP)+XAT*DSAG1(ITYP)-DLA*SYG1(ITYP) 123 : SAG1(ITYP)=SAG1(ITYP)+XAT*DSAG1(ITYP)-SUG1(ITYP) 124 : DSUG1(ITYP)=DSAG1(ITYP)*(1.D0-XAT) 125 : 105 CONTINUE 126 : C------------------------- 127 : C IL FAUT APPELLE LE COUPLAGE 128 : C------------------------- 128 : CALL CHREPE(CO11,-SI11,SUG1,SOGMA) 129 : VCRIT1=VCRITE(SUG1(2),SUG1(4),PHI1,COHE1) 130 : IF(IIMPI.EQ.28) WRITE(IOIMP,1978) VCRIT1 130 : CALL CHREPE(CO11,-SI11,DSUG1,DSOGMA) 130 : CALL CHREPE(CO11,-SI11,SAG1,DSAG1) 130 : IF (NCAS.NE.1) CALL DPCONT(DSAG1,DEFPL,XNU,YOUNG) 130 : IF (NCAS.EQ.1) CALL CPCONT(DSAG1,DEFPL,XNU,YOUNG) 131 : ICHR=1 132 : DLAM=DLAM+DLA 133 : IF(IIMPI.EQ.28) WRITE(IOIMP,1111) ITER 134 : 1111 FORMAT('0 ITER= ',E12.5/) 135 : RETURN 136 : C------------------------------------------------------------- 137 : 750 DO 103 ITYP=1,4 138 : DSIG1(ITYP)=DSAG1(ITYP)*(1.D0-XAT) 139 : 103 CONTINUE 140 : VCRT=-DSIG1(2)+PSI1*ABS(DSIG1(4)) 141 : IF(VCRT.LE.0.D0.OR.ITER.LE.1) GO TO 456 142 : ITER=ITER+1 143 : CIS1=SIGN(1.D0,DSIG1(4)) 144 : DO 97 ITYP=1,4 145 : SAG1(ITYP)=TSUG1(ITYP) 146 : DSAG1(ITYP)=DSIG1(ITYP) 147 : 97 CONTINUE 148 : DLAM=DLAM+DLA 149 : GO TO 457 150 : 456 IF (NCAS.NE.1) HOOK=YOUNG/(1.D0+XNU)/(1.D0-2.D0*XNU) 151 : IF (NCAS.NE.1) ALPHA=(DSIG1(2)+TSUG1(2)-TRA1)/HOOK /(1.D0-XNU) 152 : IF (NCAS.NE.1) BETA=(DSIG1(4)+TSUG1(4))/HOOK/(1.D0-2.D0*XNU) 153 : C-------------------------- 154 : C CAS CONT PLANES 155 : C 156 : IF (NCAS.EQ.1) HOOK=YOUNG/(1.D0-XNU*XNU) 157 : IF (NCAS.EQ.1) ALPHA=(DSIG1(2)+TSUG1(2)-TRA1)/HOOK 158 : IF (NCAS.EQ.1) BETA=(DSIG1(4)+TSUG1(4))/HOOK/(1.D0-XNU) 159 : C----------------------------------------------------- 159 : IF(NCAS.EQ.1) CALL CPHOOK(ZER,ALPHA,ZER, 160 : # BETA,DSEG1(1),DSEG1(2) 161 : # ,DSEG1(3),DSEG1(4),XNU,YOUNG) 161 : IF(NCAS.NE.1) CALL DPHOOK(ZER,ALPHA,ZER, 162 : # BETA,DSEG1(1),DSEG1(2) 163 : # ,DSEG1(3),DSEG1(4),XNU,YOUNG) 164 : DO 453 ITYP=1,4 165 : TSEG1(ITYP)=TSUG1(ITYP)+DSIG1(ITYP)-DSEG1(ITYP) 166 : DSUG1(ITYP)=DSOG1(ITYP)+SOG1(ITYP)-TSEG1(ITYP) 167 : 453 CONTINUE 168 : C------------------------------------- 169 : C ON REGARDE SI ON ENDOMMAGE LE 2 EME CRITERE 170 : C------------------------------------- 171 : C 171 : CALL CHREPE(CO21,SI21,TSEG1,TSUG2) 172 : VCRTA=VCRITE(TSUG2(2),TSUG2(4),PHI2,COHE2) 173 : IF(IIMPI.EQ.28) WRITE(IOIMP,2978) VCRTA 174 : 2978 FORMAT( '0 ECUBI VCRTA =',1PE12.5/) 175 : IF(VCRTA.LE.0.) GO TO 522 175 : C CALL CHREPE(CO21,SI21,DSUG1,DSUG2) 176 : C 177 : IF(ABS(DSEG1(4)).LE.1.D-7) ATS=0.D0 178 : IF(ABS(DSEG1(4)).LE.1.D-7) BTS=1.D0 179 : IF(ABS(DSEG1(4)).LE.1.D-7) GOTO 543 180 : C ATS=(TSUG1(4)-TRA1)/DSEG1(4) 181 : ATS=TSUG1(4)/DSEG1(4) 182 : C BTS=DSIG1(2)/DSEG1(2) 183 : BTS=DSIG1(4)/DSEG1(4) 184 : 543 DO 987 ITYP=1,4 185 : SIP1(ITYP)=TSUG1(ITYP)-ATS*DSEG1(ITYP) 186 : SAP1(ITYP)=DSIG1(ITYP)-BTS*DSEG1(ITYP) 187 : 987 CONTINUE 187 : CALL CHREPE(CO21,SI21,SIP1,SIP2) 187 : CALL CHREPE(CO21,SI21,SAP1,SAP2) 188 : XCT=-(PHI2*SIP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/ 189 : 1 (PHI2*SAP2(2)+SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) 190 : XDT=-(PHI2*SIP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SIP2(4)-COHE2)/ 191 : 1 (PHI2*SAP2(2)-SIGN(1.D0,SIP2(4)+SAP2(4))*SAP2(4)) 192 : IF(XCT.LE.0.D0) XCT=100.D0 193 : IF(XDT.LE.0.D0) XDT=100.D0 194 : XAT=MIN(XCT,XDT) 195 : IF(XAT.GT.1.) XAT=0.D0 196 : XBT=ATS+XAT*BTS 197 : IF(IIMPI.EQ.28) WRITE(IOIMP,3978) XAT,XBT,XCT,XDT 198 : 3978 FORMAT( '0 ECUBI XAT XBT XCT XDT ',4(1X,1PE12.5)/) 199 : DO 195 ITYP=1,4 200 : TSUG1(ITYP)=TSUG1(ITYP)+XAT*DSIG1(ITYP)-DSEG1(ITYP)*XBT 201 : DSUG1(ITYP)=DSOG1(ITYP)+SOG1(ITYP)-TSUG1(ITYP) 202 : DSOG1(ITYP)=(1.D0-XAT)*DSIG1(ITYP) 203 : 195 CONTINUE 203 : CALL CHREPE(CO11,-SI11,TSUG1,SOGMA) 203 : CALL CHREPE(CO11,-SI11,DSOG1,DSOGMA) 203 : CALL CHREPE(CO11,-SI11,DSUG1,DSUG) 204 : ICHR=1 204 : IF (NCAS.NE.1) CALL DPCONT(DSUG,DEFPL,XNU,YOUNG) 204 : IF (NCAS.EQ.1) CALL CPCONT(DSUG,DEFPL,XNU,YOUNG) 205 : DLAM=DLAM+DLA*XBT 206 : IF(IIMPI.EQ.28) WRITE(IOIMP,1111) ITER 207 : RETURN 207 : 522 CALL CHREPE(CO11,-SI11,TSEG1,TSUG) 208 : DO 32 ITYP=1,4 209 : TSUG1(ITYP)=TSEG1(ITYP) 210 : 32 CONTINUE 211 : DLAM=DLA+SQRT((ALPHA*ALPHA+BETA*BETA)/(PSI1*PSI1+1.D0))+DLAM 212 : GOTO 523 212 : 500 CALL CHREPE(CO11,-SI11,TSUG1,TSUG) 213 : 523 CONTINUE 214 : DO 31 ITYP=1,4 215 : DSUG1(ITYP)=-TSUG1(ITYP)+DSOG1(ITYP)+SOG1(ITYP) 216 : 31 CONTINUE 216 : CALL CHREPE(CO11,-SI11,DSUG1,DSUG) 216 : IF (NCAS.NE.1) CALL DPCONT(DSUG,DEFPL,XNU,YOUNG) 216 : IF (NCAS.EQ.1) CALL CPCONT(DSUG,DEFPL,XNU,YOUNG) 217 : IF(IIMPI.EQ.28) WRITE(IOIMP,1111) ITER 218 : RETURN 219 : END 220 :
© Cast3M 2003 - All rights reserved.
Disclaimer