1 : C NMAZAR SOURCE LJASON 07/10/24 21:15:03 5930 2 : 3 : SUBROUTINE NMAZAR(XMAT,NMATT,VAR0,VARF,NVARI,SIG0, 4 : & SIGF,DEPST,NSTRS,ISTEP) 5 : 6 : IMPLICIT INTEGER(I-N) 7 : IMPLICIT REAL*8(A-H,O-Z) 8 : 9 : * Entrees 10 : INTEGER NMATT, NVARI, NSTRS, ISTEP 11 : REAL*8 XMAT(NMATT), VAR0(NVARI), SIG0(NSTRS), DEPST(NSTRS) 12 : REAL*8 SIGF(NSTRS), VARF(NVARI) 13 : 14 : * Variables intermediaires 15 : REAL*8 EPS33(3,3),EPSIPP(3),EPSILT(3),VALP33(3,3) 16 : REAL*8 SIGP(3),SIGPT(3),SIGPC(3),TRSIGT,TRSIGC 17 : REAL*8 SIG00(6), EPSTIM 18 : REAL*8 YOUN,XNU,EPSD0,ACOM,BCOM,ATRA,BTRA,BETA 19 : REAL*8 XZERO,UN,DEUX,XPETIT, DDHOOK(6,6) 20 : REAL*8 DINI,D,DT,DC,EPSTIL,ALFAT,ALFAC,GAMMA 21 : REAL*8 EPSILO(6), PREC 22 : INTEGER ISTRS,JSTRS,KCAS,IRTD 23 : CHARACTER CMATE 24 : PARAMETER (XZERO=0.D0 , UN=1.D0 , DEUX=2.D0, XPETIT=1.D-12) 25 : 26 : CMATE = 'ISOTROPE' 27 : YOUN = XMAT(1) 28 : XNU = XMAT(2) 29 : EPSD0= XMAT(9) 30 : ACOM = XMAT(5) 31 : BCOM = XMAT(7) 32 : ATRA = XMAT(6) 33 : BTRA = XMAT(8) 34 : BETA = XMAT(20) 35 : DINI = VAR0(2) 36 : 37 : C 38 : C calcul des seuils 39 : C 40 : C 41 : C calcul de la deformation élastique 42 : C 43 : * Calcul de la matrice elastique 44 : KCAS=2 44 : CALL DOHMAS(XMAT,'ISOTROPE',2,NSTRS,KCAS,DDHOOK,IRTD) 45 : 46 : PREC=1.D-20 46 : CALL INVALM(DDHOOK,NSTRS,NSTRS,KERRE,PREC) 47 : 48 : * Calcul de La partie élastique des déformations 48 : CALL MATVE1 (DDHOOK,SIG0,NSTRS,NSTRS,EPSILO,1) 49 : 50 : * DO ISTRS = 1,NSTRS 51 : * EPSILO(ISTRS) = 1/(1-DINI)*EPSILO(ISTRS) + DEPST(ISTRS) 52 : * END DO 53 : C 54 : C calcul des deformations principales 55 : C 56 : C 57 : C on reecrit les deformations sous forme matricielle 58 : C 58 : CALL ENDOCB (EPSILO,EPS33,2,2) 59 : C 60 : C et on diagonalise 61 : C 61 : CALL JACOB3 (EPS33,3,EPSIPP,VALP33) 62 : 63 : 64 : IF (ISTEP .EQ. 0 .OR. ISTEP.EQ.2) THEN 65 : C 66 : C on calcule la matrice de hooke et les contraintes ppales 67 : C 68 : CMATE = 'ISOTROPE' 69 : KCAS=2 69 : CALL DOHMAS(XMAT,'ISOTROPE',2,NSTRS,KCAS,DDHOOK,IRTD) 70 : DO 200 ISTRS=1,3 71 : SIGP(ISTRS)= XZERO 72 : DO 210 JSTRS=1,3 73 : SIGP(ISTRS)=SIGP(ISTRS)+DDHOOK(ISTRS,JSTRS)*EPSIPP(JSTRS) 74 : 210 CONTINUE 75 : 200 CONTINUE 76 : END IF 77 : 78 : C on calcule le epsilontild 79 : C 80 : EPSTIL=MAX( XZERO , EPSIPP(1) )**2 + 81 : 1 MAX( XZERO , EPSIPP(2) )**2 + 82 : 2 MAX( XZERO , EPSIPP(3) )**2 83 : EPSTIL=SQRT (EPSTIL) 84 : 85 : IF (ISTEP .EQ. 0) THEN 86 : EPSTIM=EPSTIL 87 : VARF(1)=EPSTIL 88 : ELSE IF (ISTEP .EQ. 1) THEN 89 : VARF(9)=VAR0(9) 90 : VARF(3)=VAR0(3) 91 : VARF(4)=VAR0(4) 92 : VARF(5)=VAR0(5) 93 : VARF(6)=VAR0(6) 94 : VARF(7)=VAR0(7) 95 : VARF(8)=VAR0(8) 96 : VARF(2)=DINI 97 : VARF(1)=EPSTIL 98 : ELSE IF (ISTEP .EQ. 2) THEN 99 : EPSTIM=VAR0(1) 100 : VARF(1)=EPSTIM 101 : ELSE 102 : PRINT*,'DANS MAZARS ISTEP = 0,1,2 ET PAS ',ISTEP 103 : END IF 104 : 105 : IF ( (ISTEP .EQ. 0) .OR. (ISTEP .EQ. 2)) THEN 106 : C 107 : C on calcule l'endommagement et les contraintes 108 : C 109 : C 110 : C on calcule ALFAT ALFAC DT et DC puis D 111 : C dans le cas ou le seuil initial est depasse 112 : C 113 : IF ( EPSTIM .GT. EPSD0) THEN 114 : C 115 : C calcul de l'endommagement 116 : C 117 : C 118 : C on calcule le signe des contraintes elastiques 119 : C 120 : DO 300 ISTRS=1,3 121 : IF (SIGP(ISTRS).LT. XZERO) THEN 122 : SIGPC(ISTRS)=SIGP(ISTRS) 123 : SIGPT(ISTRS)=XZERO 124 : ELSE 125 : SIGPT(ISTRS)=SIGP(ISTRS) 126 : SIGPC(ISTRS)=XZERO 127 : END IF 128 : 300 CONTINUE 129 : TRSIGT=SIGPT(1)+SIGPT(2)+SIGPT(3) 130 : TRSIGC=SIGPC(1)+SIGPC(2)+SIGPC(3) 131 : C 132 : C on calcule les deformations dues aux contraintes positives 133 : C 134 : DO 400 ISTRS=1,3 135 : EPSILT(ISTRS)=(SIGPT(ISTRS)*(UN+XNU)-TRSIGT*XNU)/YOUN 136 : 400 CONTINUE 137 : C 138 : C on en deduit ALFAT et ALFAC 139 : C 140 : ALFAT = MAX(XZERO,EPSIPP(1))*EPSILT(1) + 141 : 1 MAX(XZERO,EPSIPP(2))*EPSILT(2) + 142 : 2 MAX(XZERO,EPSIPP(3))*EPSILT(3) 143 : ALFAT = ALFAT/(EPSTIL*EPSTIL) 144 : ALFAC = UN - ALFAT 145 : C 146 : C modification pour la bi ou tricompression 147 : C 148 : * IF (TRSIGC.LT. -XPETIT .AND. TRSIGT.LT.XPETIT) THEN 149 : * GAMMA=SIGPC(1)*SIGPC(1)+SIGPC(2)*SIGPC(2)+ 150 : * 1 SIGPC(3)*SIGPC(3) 151 : * GAMMA=-SQRT(GAMMA)/TRSIGC 152 : * EPSTIM=EPSTIM*GAMMA 153 : * END IF 154 : C 155 : C amelioration de la reponse en cisaillement pour beta > 1. 156 : C 157 : IF (BETA .GT. UN) THEN 158 : IF ( ALFAT .GT. XPETIT ) THEN 159 : ALFAT=ALFAT**BETA 160 : END IF 161 : IF ( ALFAC .GT. XPETIT ) THEN 162 : ALFAC=ALFAC**BETA 163 : END IF 164 : END IF 165 : 166 : C 167 : C on calcule DT et DC puis D 168 : C dans le cas ou le seuil initial est depasse 169 : C 170 : C on est oblige de verifier car on a pu multiplier par gamma 171 : C 172 : IF (EPSTIM .GT. EPSD0) THEN 173 : DT=UN - EPSD0*(UN-ATRA)/EPSTIM - 174 : 1 ATRA*EXP(-BTRA*(EPSTIM-EPSD0)) 175 : DC=UN - EPSD0*(UN-ACOM)/EPSTIM - 176 : 1 ACOM*EXP(-BCOM*(EPSTIM-EPSD0)) 177 : ELSE 178 : DT=XZERO 179 : DC=XZERO 180 : END IF 181 : D = ALFAT*DT + ALFAC*DC 182 : C 183 : C on borne la valeur de D a 0.99999 184 : C 185 : D=MIN ( D , UN-1.D-05 ) 186 : ELSE 187 : D=XZERO 188 : END IF 189 : C 190 : C on teste la croissance de D 191 : C 192 : D=MAX ( D , DINI ) 193 : C 194 : C on le stocke dans les variables internes finales 195 : C 196 : VARF(2)= D 197 : C 198 : C on calcule les contraintes finales 199 : C 200 : DO 500 ISTRS=1,NSTRS 201 : SIGF(ISTRS)=SIG0(ISTRS)*(UN-D) 202 : 500 CONTINUE 203 : 204 : * Et on stocke la différence entre contraintes réelles et contraintes plastiques 205 : 206 : DO ISTRS=1,NSTRS 207 : VARF(ISTRS+2) = SIG0(ISTRS)-SIGF(ISTRS) 208 : END DO 209 : 210 : END IF 211 : 212 : RETURN 213 : END 214 : 215 :
© Cast3M 2003 - All rights reserved.
Disclaimer