1 : C LCGLDM SOURCE BR232186 12/10/04 21:15:18 7520 2 : SUBROUTINE LCGLDM (XMAT,VAR0,VARF,SIG0,SIGF,DEPST,XCAR) 3 : C 4 : C====&===1=========2=========3=========4=========5=========6=========7== 5 : C Commentaires : Subroutine permettant de mettre en oeuvre le 6 : C modele LGCR_DM (EDF R&D/AMA) pour representer 7 : C le comportement cyclique du béton armé sous 8 : C sollicitations cycliques 9 : C 10 : C Traits : - Endommagement scalaire 11 : C - 2D (DKT/COQUES/PLAQUES MINCES) pré-intégrée 12 : C - Monotone et cyclique 13 : C - Effet unilateral (retour à l'origine) 14 : C 15 : C Auteur : B. Richard (Dr - Ing.) - CEA/DEN/DANS/DM2S/SEMT/EMSI 16 : C====&===1=========2=========3=========4=========5=========6=========7== 17 : C 18 : C----DECLARATION GENERALES---------------------------------------------- 19 : C 20 : IMPLICIT INTEGER(I-N) 21 : IMPLICIT REAL*8(A-H,O-Z) 22 : C 23 : C----DECLARATION SEGMENTS----------------------------------------------- 24 : C 25 : REAL*8 EPSIM(3,3),EPSIF(3,3),UNIT(2,2),SIGMAM(3,3) 26 : REAL*8 SIGMAF(3,3),VAPM(3),VEPM(3,3),VAPF(3),VEPF(3,3) 27 : REAL*8 XMAT(*),SIG0(*),SIGF(*),VAR0(*),VARF(*),DEPST(*),XCAR(*) 28 : REAL*8 XLAMM,XLAMF,XMUM,XMUF,XD1,XD2,XTOL 29 : REAL*8 XTRAV1(3,3),XTRAV2(3,3),VEPMT(3,3),VEPFT(3,3) 30 : REAL*8 XTRAV3(3,3),XTRAV4(3,3),XIMI(2),XIFI(2),BETAF(3,3) 31 : REAL*8 XTRAV5(2,2),XTRAV6(2,2),BETAM(3,3),SGEFF(3) 32 : C 33 : C----DECLARATION PARAMETRES--------------------------------------------- 34 : C 35 : PARAMETER (UN=1.D0) 36 : PARAMETER (XTOL=1.D-8) 37 : 37 : CALL ZERO(EPSIM,3,3) 37 : CALL ZERO(EPSIF,3,3) 38 : 39 : DO I=1,2 40 : DO J=1,2 41 : IF (I.EQ.J) THEN 42 : UNIT(I,J) = 1.0D0 43 : ELSE 44 : UNIT(I,J) = 0.0D0 45 : ENDIF 46 : ENDDO 47 : ENDDO 48 : C 49 : C----MISE EN DONNEES---------------------------------------------------- 50 : C 51 : YOUN = XMAT(1) 52 : XNU = XMAT(2) 53 : 54 : XYOUNM = XMAT(1) 55 : XNUM = XMAT(2) 56 : XYOUNF = XMAT(5) 57 : XNUF = XMAT(6) 58 : 59 : EPAI = XCAR(1) 60 : 61 : XLAMM = XNUM*XYOUNM*EPAI/((1.0D0+XNUM)*(1.0D0-2.0D0*XNUM)) 62 : XMUM = XYOUNM*EPAI/(2.0D0*(1.0D0+XNUM)) 63 : 64 : XLAMF = (XNUF*XYOUNF*(EPAI**(3.0D0)))/ 65 : & (12.0D0*(1.0D0-(XNUF**(2.0D0)))) 66 : XMUF = (XYOUNF*(EPAI**(3.0D0)))/ 67 : & (24.0D0*(1.0D0+XNUF)) 68 : 69 : XGAMMT = XMAT(7) 70 : XGAMMC = XMAT(8) 71 : XGAMF = XMAT(9) 72 : XSEUI0 = XMAT(10) 73 : XALFA = XMAT(11) 74 : C 75 : C----VARIABLES INTERNES------------------------------------------------- 76 : C 77 : EPSIM(1,1) = VAR0(1) 78 : EPSIM(2,2) = VAR0(2) 79 : EPSIM(1,2) = VAR0(3)/2.0D0 80 : EPSIM(2,1) = VAR0(3)/2.0D0 81 : 82 : EPSIF(1,1) = VAR0(4) 83 : EPSIF(2,2) = VAR0(5) 84 : EPSIF(1,2) = VAR0(6)/2.0D0 85 : EPSIF(2,1) = VAR0(6)/2.0D0 86 : 87 : XD1 = VAR0(9) 88 : XD2 = VAR0(10) 89 : 90 : XD10 = VAR0(9) 91 : XD20 = VAR0(10) 92 : 93 : C print*,'DEPST',(DEPST(I),I=1,3) 94 : 95 : C 96 : C----ACTUALISATION DEFORMATION------------------------------------------ 97 : C 98 : EPSIM(1,1) = EPSIM(1,1) + DEPST(1) 99 : EPSIM(2,2) = EPSIM(2,2) + DEPST(2) 100 : EPSIM(1,2) = EPSIM(1,2) + DEPST(3)/2.0D0 101 : EPSIM(2,1) = EPSIM(2,1) + DEPST(3)/2.0D0 102 : 103 : EPSIF(1,1) = EPSIF(1,1) + DEPST(4) 104 : EPSIF(2,2) = EPSIF(2,2) + DEPST(5) 105 : EPSIF(1,2) = EPSIF(1,2) + DEPST(6)/2.0D0 106 : EPSIF(2,1) = EPSIF(2,1) + DEPST(6)/2.0D0 107 : C 108 : C----PREDICTION ELASTIQUE----------------------------------------------- 109 : C 110 : 111 : C Deformation hors plan 111 : CALL CALCEP33(EPSIM,XGAMMT,XGAMMC,XD1,XD2,XMUM,XLAMM,XEPS33) 112 : C 113 : C----CALCUL DES SEUILS YM1 ET YM2--------------------------------------- 114 : C 114 : CALL CALCYMI(EPSIM,XGAMMT,XGAMMC,XLAMM,XMUM,XD1,XEPS33,YM1) 114 : CALL CALCYMI(EPSIM,XGAMMT,XGAMMC,XLAMM,XMUM,XD2,XEPS33,YM2) 115 : C 116 : C----CALCUL DES SEUILS YF1 ET YF2--------------------------------------- 117 : C 118 : 119 : C Calcul de YF01 et YF02 119 : CALL CALCYF01(EPSIF,XGAMF,XLAMF,XMUF,YF01) 119 : CALL CALCYF02(EPSIF,XGAMF,XLAMF,XMUF,YF02) 120 : 121 : C Calcul de YF1 et YF2 122 : YF1 = (XALFA/(XD1+XALFA)**2.0D0)*YF01 123 : YF2 = (XALFA/(XD2+XALFA)**2.0D0)*YF02 124 : C 125 : C----CALCUL DES SEUILS Y1 ET Y2----------------------------------------- 126 : C 127 : XSEUIL1 = YM1 + YF1 - XSEUI0 128 : XSEUIL2 = YM2 + YF2 - XSEUI0 129 : C 130 : C----TEST SEUIL 1------------------------------------------------------- 131 : C 132 : IF (XSEUIL1.GT.0.0D0) THEN 133 : 134 : C Debut des iterations internes 135 : DO NINTER=1,100 136 : 137 : C Calcul du residu 1 137 : CALL CALCYMI(EPSIM,XGAMMT,XGAMMC,XLAMM,XMUM,XD1,XEPS33,YM1) 138 : 138 : CALL CALCYF01(EPSIF,XGAMF,XLAMF,XMUF,YF01) 139 : 140 : YF1 = (XALFA/(XD1+XALFA)**2.0D0)*YF01 141 : 142 : XRES1 = YM1 + YF1 - XSEUI0 143 : 144 : C Calcul de l'operateur tangent 1 144 : CALL CALCRTG1(EPSIM,EPSIF,XGAMMT,XGAMMC,XGAMF,XALFA,XLAMM, 145 : & XMUM,XLAMF,XMUF,XD1,XD2,XEPS33,DR1DD1) 146 : 147 : C Corrections de l'endommagement 148 : XDD1 = -1.0D0*(1.0D0/DR1DD1)*XRES1 149 : 150 : C Critere d arret des iterations internes 151 : XCRI1 = DABS(XDD1*DR1DD1/XSEUI0) 152 : 153 : IF (XCRI1.LE.XTOL) THEN 154 : 155 : GOTO 111 156 : 157 : ELSE 158 : 159 : XD1 = XD1 + XDD1 160 : 160 : CALL CALCEP33(EPSIM,XGAMMT,XGAMMC,XD1,XD2,XMUM,XLAMM, 161 : & XEPS33) 162 : 163 : ENDIF 164 : 165 : ENDDO 166 : ENDIF 167 : 168 : 111 CONTINUE 169 : 170 : C 171 : C----TEST SEUIL 2------------------------------------------------------- 172 : C 173 : IF (XSEUIL2.GT.0.0D0) THEN 174 : 175 : C Debut des iterations internes 176 : DO NINTER=1,100 177 : 178 : C Calcul du residu 1 178 : CALL CALCYMI(EPSIM,XGAMMT,XGAMMC,XLAMM,XMUM,XD2,XEPS33,YM2) 179 : 179 : CALL CALCYF02(EPSIF,XGAMF,XLAMF,XMUF,YF02) 180 : 181 : YF2 = (XALFA/(XD2+XALFA)**2.0D0)*YF02 182 : 183 : XRES2 = YM2 + YF2 - XSEUI0 184 : 185 : C Calcul de l'operateur tangent 1 185 : CALL CALCRTG2(EPSIM,EPSIF,XGAMMT,XGAMMC,XGAMF,XALFA,XLAMM, 186 : & XMUM,XLAMF,XMUF,XD1,XD2,XEPS33,DR2DD2) 187 : 188 : C Corrections de l'endommagement 189 : XDD2 = -1.0D0*(1.0D0/DR2DD2)*XRES2 190 : 191 : C Critere d arret des iterations internes 192 : XCRI2 = DABS(XDD2*DR2DD2/XSEUI0) 193 : 194 : IF (XCRI2.LE.XTOL) THEN 195 : 196 : GOTO 112 197 : 198 : ELSE 199 : 200 : XD2 = XD2 + XDD2 201 : 201 : CALL CALCEP33(EPSIM,XGAMMT,XGAMMC,XD1,XD2,XMUM,XLAMM, 202 : & XEPS33) 203 : 204 : ENDIF 205 : 206 : ENDDO 207 : ENDIF 208 : 209 : 112 CONTINUE 210 : 211 : XTRA = EPSIM(1,1)+EPSIM(2,2) 211 : CALL CALCXIM (XTRA,XGAMMT,XGAMMC,XD1,XD2,XIMF) 212 : C 213 : C-----CALCUL DE CONTRAINTES--------------------------------------------- 214 : C 215 : 216 : C Test sur les endommagements pour obtenir une dissipation positive 217 : IF (XD1.LE.XD10) THEN 218 : XD1 = XD10 219 : ENDIF 220 : 221 : IF (XD2.LE.XD20) THEN 222 : XD2 = XD20 223 : ENDIF 224 : 225 : C Cacul du repere propre 225 : CALL JACOB4(EPSIM,2,VAPM,VEPM) 226 : C print*,'EPSIM=',((EPSIM(I,J),I=1,2),J=1,2) 227 : C print*,'VAPM=',(VAPM(I),I=1,3) 228 : C print*,'VEPM=',((VEPM(I,J),I=1,2),J=1,2) 229 : 229 : CALL JACOB4(EPSIF,2,VAPF,VEPF) 230 : 231 : C print*,'EPSIF=',((EPSIF(I,J),I=1,3),J=1,3) 232 : C print*,'VAPF=',(VAPF(I),I=1,3) 233 : C print*,'VEPF=',((VEPF(I,J),I=1,2),J=1,2) 234 : 235 : C Calcul de l'endommagement xim, xim1 et xim2 236 : XTRA = EPSIM(1,1)+EPSIM(2,2) 236 : CALL CALCXIM(XTRA,XGAMMT,XGAMMC,XD1,XD2,XIM) 236 : CALL CALCXIM(VAPM(1),XGAMMT,XGAMMC,XD1,XD2,XIM1) 236 : CALL CALCXIM(VAPM(2),XGAMMT,XGAMMC,XD1,XD2,XIM2) 237 : 238 : XIMI(1) = XIM1 239 : XIMI(2) = XIM2 240 : 241 : C Calcul de l'endommagement xif, xif1 et xif2 242 : XTRAK = EPSIF(1,1)+EPSIF(2,2) 242 : CALL CALCXIF(XTRAK,XALFA,XGAMF,XD1,XD2,XIF) 242 : CALL CALCXIF(VAPF(1),XALFA,XGAMF,XD1,XD2,XIF1) 242 : CALL CALCXIF(VAPF(2),XALFA,XGAMF,XD1,XD2,XIF2) 243 : 244 : XIFI(1) = XIF1 245 : XIFI(2) = XIF2 246 : 247 : C Calcul du repere propre des deformations generalisees 247 : CALL TRSPOD(VEPM,3,3,VEPMT) 247 : CALL TRSPOD(VEPF,2,2,VEPFT) 248 : 249 : C Calcul des contraintes generalisees dans le repere propre 250 : XTRA = EPSIM(1,1)+EPSIM(2,2) 251 : XTRAK = EPSIF(1,1)+EPSIF(2,2) 252 : 252 : CALL ZERO(SIGMAM,3,3) 252 : CALL ZERO(SIGMAF,3,3) 253 : 254 : XTRAK = EPSIF(1,1) + EPSIF(2,2) 255 : 256 : DO I = 1,2 257 : SIGMAM(I,I) = XIM*XLAMM*(XTRA+XEPS33)+ 258 : & 2.0D0*XMUM*VAPM(I)*XIMI(I) 259 : 260 : C Abandon de la non linearite en partie flexion 261 : C SIGMAF(I,I) = XIF*XLAMF*XTRAK+ 262 : C & 2.0D0*XMUF*VAPF(I)*XIFI(I) 263 : ENDDO 264 : 265 : C Calcul des contraintes generalisees dans le repere global 265 : CALL ZERO(XTRAV5,3,3) 265 : CALL ZERO(XTRAV6,3,3) 266 : 266 : CALL MATMAT(VEPM,SIGMAM,3,3,3,XTRAV5) 266 : CALL ZERO(SIGMAM,3,3) 266 : CALL MATMAT(XTRAV5,VEPMT,3,3,3,SIGMAM) 267 : 268 : C Abandon de la non linearite en partie flexion 268 : C CALL MATMAT(VEPF,SIGMAF,3,3,3,XTRAV6) 268 : C CALL ZERO(SIGMAF,3,3) 268 : C CALL MATMAT(XTRAV6,VEPFT,3,3,3,SIGMAF) 269 : C 270 : C----CONTRAINTES------------------------------------------------------- 271 : C 272 : SIGF(1) = SIGMAM(1,1) 273 : SIGF(2) = SIGMAM(2,2) 274 : SIGF(3) = SIGMAM(1,2) 275 : 276 : 277 : FAC = (EPAI**3)/12.D0 278 : AUX = FAC*XYOUNF/(1.D0-XNUF*XNUF) 279 : AUX1 = FAC*XYOUNF*.5D0/(1.D0+XNUF) 280 : AUX2 = EPAI*XYOUNF*.5D0/(1.d0+XNUF)/1.2d0 281 : 282 : SIGF(4) = SIG0(4)+AUX*(DEPST(4)+XNUF*DEPST(5)) 283 : SIGF(5) = SIG0(5)+AUX*(DEPST(5)+XNUF*DEPST(4)) 284 : SIGF(6) = SIG0(6)+AUX1*DEPST(6) 285 : 286 : C print*,'SIGF(I)',(SIGF(I),I=4,6) 287 : 288 : C Cisaillement élastique 289 : SIGF(7) = SIG0(7)+AUX2*DEPST(7) 290 : SIGF(8) = SIG0(8)+AUX2*DEPST(8) 291 : 292 : C 293 : C----VARIABLES INTERNES------------------------------------------------- 294 : C 295 : VARF(1) = EPSIM(1,1) 296 : VARF(2) = EPSIM(2,2) 297 : VARF(3) = EPSIM(1,2)*2.0D0 298 : 299 : VARF(4) = EPSIF(1,1) 300 : VARF(5) = EPSIF(2,2) 301 : VARF(6) = EPSIF(1,2)*2.0D0 302 : 303 : VARF(9) = XD1 304 : VARF(10) = XD2 305 : 306 : VARF(11) = 1.0D0-XIM 307 : VARF(12) = 1.0D0-XIF 308 : 309 : RETURN 310 : END 311 : 312 : 313 : 314 : 315 : 316 : 317 : 318 : 319 : 320 : 321 : 322 : 323 : 324 : 325 : 326 : 327 : 328 : 329 :
© Cast3M 2003 - All rights reserved.
Disclaimer