Download ricbet3d.eso

Back to the list

   1 : C RICBET3D  SOURCE    CHAT      12/08/06    21:15:30     7479           
   2 :           SUBROUTINE RICBET3D(XMAT,SIG0,SIGF,VAR0,VARF,DEPST,ISTEP)
   3 : C
   4 : C====&===1=========2=========3=========4=========5=========6=========7==
   5 : C Commentaires : Subroutine permettant de mettre en oeuvre le 
   6 : C                modele RICBET (B. Richard) pour representer
   7 : C                le comportement 3D d'un béton ordinaire sous                 
   8 : C                sollicitations cycliques
   9 : C
  10 : C Traits       : - Endommagement scalaire 
  11 : C                - Glissement frottant
  12 : C                - 3D
  13 : C                - Monotone et cyclique
  14 : C                - Effet unilatteral
  15 : C                - Comportement lineaire en compression
  16 : C                - Non localite portant sur le seuil d'endommagement (Fd)
  17 : C
  18 : C Auteurs      : B. Richard (Dr - Ing.) & F. Ragueneau (Pr) - LMT/ENS
  19 : C Date         : Mars 2011
  20 : C====&===1=========2=========3=========4=========5=========6=========7==
  21 : C
  22 : C-----DECLARATION GENERALE----------------------------------------------
  23 : C
  24 :       IMPLICIT REAL*8(A-H,O-Z)
  25 : C
  26 : C-----DECLARATION GENERALE----------------------------------------------
  27 : C
  28 :       REAL*8 X(3,3),UNIT(3,3),EPSPI(3,3),EPSI(3,3),SIGMA(3,3)
  29 :       REAL*8 EPSDEV(3,3),SIGMA1(3,3),EPSPIDEV(3,3)
  30 :       REAL*8 EPSIPOS(3,3),EPSINEG(3,3),SIGPI(3,3),DEPSN(3,3)
  31 :       REAL*8 TERME1(3,3),TERME2(3,3),TERME3(3,3),TERME4(3,3)
  32 :       REAL*8 SIGMOX(3,3),TRAV(3,3),SIG1DEV(3,3),TERME5(3,3)
  33 :       REAL*8 XDT,XTRAV2(3,3),XTRAV3(3,3),DEFN(3,3),EPSIT(3,3)     
  34 :       REAL*8 XMAT(*),SIG0(*),SIGF(*),VAR0(*),VARF(*),DEPST(*)
  35 :       REAL*8 XDFP(3,3),XDGP(3,3),XDJ2P(3,3),DEVGP(3,3)
  36 : C
  37 : C-----DEFINITION DE PARAMETRES------------------------------------------
  38 : C
  39 :       PARAMETER (IMAX = 1000)
  40 : C
  41 : C-----MISE EN DONNEE----------------------------------------------------
  42 : C       
  43 :       XE = XMAT(1)
  44 :       XNU = XMAT(2)
  45 :       XFT = XMAT(5)
  46 :       ALDIR = XMAT(6)
  47 :       GAMMA1 = XMAT(7)
  48 :       A1 = XMAT(8)
  49 :       SIGRF = XMAT(9)
  50 : 
  51 :       XAF = XMAT(10)
  52 :       XBF = XMAT(11)
  53 :       XAG = XMAT(12)
  54 :       XBG = XMAT(13)
  55 :       XAC = XMAT(14)
  56 :       XBC = XMAT(15)
  57 :       XSIGU = XMAT(16)
  58 :       XFC = XMAT(17)
  59 : 
  60 :       XK = XE/((1.0D0-2.D0*XNU)*3.0D0)
  61 :       XG = XE/(2.0D0*(1.D0+XNU))
  62 : C
  63 : C-----DEFINITION DE LA MATRICE UNITE D ORDRE 2------------------------------
  64 : C
  65 :           DO I=1,3
  66 :          DO J=1,3
  67 :             IF (I.EQ.J) THEN
  68 :                UNIT(I,J) = 1.0D0
  69 :             ELSE
  70 :                UNIT(I,J) = 0.0D0
  71 :             ENDIF
  72 :          ENDDO
  73 :       ENDDO
  74 : C
  75 : C-----VARIABLES INTERNES----------------------------------------------------
  76 : C
  77 : 
  78 : C     Tenseur de glissement interne
  79 :       EPSPI(1,1) = VAR0(2)
  80 :       EPSPI(2,2) = VAR0(3)
  81 :       EPSPI(3,3) = VAR0(4)
  82 :       EPSPI(1,2) = VAR0(5)
  83 :       EPSPI(1,3) = VAR0(6)
  84 :       EPSPI(2,3) = VAR0(7)
  85 :       EPSPI(2,1) = VAR0(5)
  86 :       EPSPI(3,1) = VAR0(6)
  87 :       EPSPI(3,2) = VAR0(7)
  88 : 
  89 : C     Tenseur de la contrainte de rappel
  90 :       X(1,1) = VAR0(8)
  91 :       X(2,2) = VAR0(9)
  92 :       X(3,3) = VAR0(10)
  93 :       X(1,2) = VAR0(11)
  94 :       X(1,3) = VAR0(12)
  95 :       X(2,3) = VAR0(13)
  96 :       X(2,1) = VAR0(11)
  97 :       X(3,1) = VAR0(12)
  98 :       X(3,2) = VAR0(13)
  99 : 
 100 : C     Endommagement
 101 :       XDT = VAR0(14)
 102 : 
 103 : C     Force associee a l ecrouissage isotrope
 104 :       XFFT = VAR0(15)
 105 : 
 106 : C     Deformations totales      
 107 :       EPSI(1,1) = VAR0(16)
 108 :       EPSI(2,2) = VAR0(17)
 109 :       EPSI(3,3) = VAR0(18)
 110 :       EPSI(1,2) = VAR0(19)
 111 :       EPSI(1,3) = VAR0(20)
 112 :       EPSI(2,3) = VAR0(21)
 113 :       EPSI(2,1) = VAR0(19)
 114 :       EPSI(3,1) = VAR0(20)
 115 :       EPSI(3,2) = VAR0(21)
 116 : 
 117 : C     Tenseur des deformations plastiques
 118 :       DEFN(1,1) = VAR0(35)
 119 :       DEFN(2,2) = VAR0(36)
 120 :       DEFN(3,3) = VAR0(37)
 121 :       DEFN(1,2) = VAR0(38)
 122 :       DEFN(1,3) = VAR0(39)
 123 :       DEFN(2,3) = VAR0(40)
 124 :       DEFN(2,1) = VAR0(38)
 125 :       DEFN(3,1) = VAR0(39)
 126 :       DEFN(3,2) = VAR0(40)
 127 : 
 128 : C     Ecrouissage plastique
 129 :       XECO = VAR0(41)
 130 : 
 131 : C     Deformation plastique cumulee
 132 :       XDPP = VAR0(42)
 133 : 
 134 : C     Contraintes tests      
 135 :       SIGMA1(1,1) = VAR0(22)
 136 :       SIGMA1(2,2) = VAR0(23)
 137 :       SIGMA1(3,3) = VAR0(24)
 138 :       SIGMA1(1,2) = VAR0(25)
 139 :       SIGMA1(1,3) = VAR0(26)
 140 :       SIGMA1(2,3) = VAR0(27)
 141 :       SIGMA1(2,1) = VAR0(25)
 142 :       SIGMA1(3,1) = VAR0(26)
 143 :       SIGMA1(3,2) = VAR0(27)
 144 :       
 145 : C     Indicateur
 146 :       XBEN = VAR0(22)+VAR0(23)+VAR0(24)
 147 : 
 148 : C     Contrainte initiales
 149 :       SIGMA(1,1) = SIG0(1)
 150 :       SIGMA(2,2) = SIG0(2)
 151 :       SIGMA(3,3) = SIG0(3)
 152 :       SIGMA(1,2) = SIG0(4)
 153 :       SIGMA(1,3) = SIG0(5)
 154 :       SIGMA(2,3) = SIG0(6)
 155 :       SIGMA(2,1) = SIG0(4)
 156 :       SIGMA(3,1) = SIG0(5)
 157 :       SIGMA(3,2) = SIG0(6)
 158 : 
 159 : C
 160 : C-----ACTUALISATION DES DEFORMATIONS------------------------------------
 161 : C       
 162 :       EPSI(1,1) = EPSI(1,1) + DEPST(1)
 163 :       EPSI(2,2) = EPSI(2,2) + DEPST(2)
 164 :       EPSI(3,3) = EPSI(3,3) + DEPST(3)
 165 :       EPSI(1,2) = EPSI(1,2) + (DEPST(4)/2.0D0)
 166 :       EPSI(1,3) = EPSI(1,3) + (DEPST(5)/2.0D0)
 167 :       EPSI(2,3) = EPSI(2,3) + (DEPST(6)/2.0D0)
 168 :       EPSI(2,1) = EPSI(2,1) + (DEPST(4)/2.0D0)
 169 :       EPSI(3,1) = EPSI(3,1) + (DEPST(5)/2.0D0)
 170 :       EPSI(3,2) = EPSI(3,2) + (DEPST(6)/2.0D0)
 171 :           
 172 : C
 173 : C-----CALCUL DE QUELQUES QUANTITES--------------------------------------
 174 : C
 175 : 
 176 : C     Calcul de la trace de EPSI
 177 :       TRA = EPSI(1,1)+EPSI(2,2)+EPSI(3,3)
 178 : 
 179 : C     Calcul du deviateur de EPSI
 179 :       CALL DEVIAT(EPSI,EPSDEV)
 180 :           
 181 : C     Calcul de la trace de EPSPI
 182 :       TRAPI = EPSPI(1,1)+EPSPI(2,2)+EPSPI(3,3)
 183 :       
 184 : C     Calcul du deviateur de EPSPI
 184 :       CALL DEVIAT(EPSPI,EPSPIDEV)
 185 :       
 186 : C     Calcul du deviateur de DEFN
 186 :       CALL DEVIAT(DEFN,DEPSN)
 187 : 
 188 : C     Calcul de la trace de DEFN
 189 :       TRADN = DEFN(1,1)+DEFN(2,2)+DEFN(3,3)
 190 :       
 191 : C
 192 : C-----PARTIE POSITIVES ET NEGATIVES DE (EPS)-----------------------------
 193 : C
 194 :       DO I=1,3
 195 :          DO J=1,3
 196 :             EPSIT(I,J) = EPSI(I,J)-DEFN(I,J)
 197 :          ENDDO
 198 :       ENDDO
 199 : 
 199 :       CALL POSNEG(EPSIT,EPSIPOS,EPSINEG)
 200 : C
 201 : C-----SEUIL ENDOMMAGEMENT------------------------------------------------
 202 : C
 203 : 
 204 : C     Calcul de l enerigie Y+
 205 :       YP = 0.0D0
 206 :       
 207 :       DO I=1,3
 208 :          DO J=1,3
 209 :             YP = YP + EPSIPOS(I,J)*EPSIPOS(I,J)
 210 :          ENDDO
 211 :       ENDDO
 212 : 
 213 :       YP = 0.5D0*YP*XE
 214 : 
 215 : C
 216 : C-----NON LOCAL----------------------------------------------------------
 217 : C
 218 :       IF (ISTEP.EQ.0) THEN
 219 :          SEUIL1=YP
 220 :          VARF(1)=YP
 221 :       ELSE IF (ISTEP.EQ.1) THEN
 222 :          VARF(1)=YP
 223 :          DO I=2,42
 224 :             VARF(I)=VAR0(I)
 225 :          ENDDO
 226 :          GOTO 2000
 227 :       ELSE IF (ISTEP.EQ.2) THEN
 228 :          SEUIL1=VAR0(1)
 229 :          VARF(1)=SEUIL1
 230 :       ENDIF 
 231 : C
 232 : C-----EVOLUTION DE L ENDOMMAGEMENT---------------------------------------
 233 : C
 234 : 
 235 : C     Calcul du seuil
 236 :       XFD = SEUIL1 - ((XFT**2)/(2.0D0*XE) + XFFT)
 237 : 
 238 : C     Evolution explicite de l'endommagement
 239 :       IF ((XFD.GT.0.0D0).AND.(XBEN.GE.0.0D0)) THEN
 240 :          XDT  = 1.0D0-1.0D0/(1.0D0+ALDIR*(SEUIL1-((XFT**2)/(2.0D0*XE))))
 241 :          XFFT = SEUIL1 - ((XFT**2)/(2.0D0*XE))
 242 :       ENDIF   
 243 : C
 244 : C-----SEUIL GLISSEMENT INTERNE-------------------------------------------
 245 : C
 246 : 
 247 : C     Calcul de sigma_pi
 248 :       DO I=1,3
 249 :          DO J=1,3
 250 :             SIGPI(I,J) = XDT*XK*(TRA-TRAPI-TRADN)*UNIT(I,J)+
 251 :      &      2.0D0*XG*XDT*(EPSDEV(I,J)-EPSPIDEV(I,J)-DEPSN(I,J))
 252 :          ENDDO
 253 :       ENDDO
 254 : 
 255 : C     Calcul de la quantite (sigma_pi - X)
 256 :       DO I=1,3
 257 :          DO J=1,3
 258 :             SIGMOX(I,J) = SIGPI(I,J)-X(I,J)
 259 :          ENDDO
 260 :       ENDDO
 261 : 
 262 : C     Calcul de I2(sigma_pi-X)
 262 :       CALL I2(SIGMOX,XFPI)
 263 : 
 264 : C     Sigma_test = Sigma_finale
 265 :       DO I=1,3
 266 :          DO J=1,3
 267 :             SIGMA(I,J)=SIGMA1(I,J)
 268 :          ENDDO
 269 :       ENDDO
 270 : 
 271 : C     Mise a 1 du multiplicateur plastique
 272 :       XLAM = 1.0D0
 273 : 
 274 : C     CAlcul de la trace de Sigma_finale
 275 :       TRASIG1 = SIGMA1(1,1)+SIGMA1(2,2)+SIGMA1(3,3)      
 276 : 
 277 : C     Test sur la positivite du seuil de glissement interne
 278 :       IF ((XFPI.GT.0.0D0).AND.(XBEN.GE.0.D0)) THEN      
 279 :          
 280 : C         On retient la valeur initiale du seuil
 281 :           XFPI0 = XFPI
 282 : 
 283 : C         On itere en bloquant le nombre d iteration max a IMAX
 284 :           DO NI=1,IMAX
 285 : 
 286 : C             On realise un test de sortie
 287 :               IF ((DABS(XFPI/XFPI0 ).LE.1.0D-8).OR.
 288 :      &            (DABS(XLAM).LE.(1.0D-10))) THEN
 289 : 
 290 :                  GOTO 400
 291 : 
 292 :               ELSE
 293 : 
 294 : C                Calcul des derivees dfpi/dsigpi et dphipi/dsigpi
 294 :                  CALL DI2(SIGMOX,TERME1)
 294 :                  CALL DI2(SIGMOX,TERME2)
 295 : 
 296 : C                Calcul des derivees dfpi/dX et dphipi/dX
 297 :                  DO I=1,3
 298 :                     DO J=1,3
 299 :                        TERME3(I,J) = 0.0D0
 300 :                        TERME4(I,J) = 0.0D0
 301 : 
 302 :                        TERME3(I,J) = -1.0D0*TERME2(I,J)
 303 :                        TERME4(I,J) = -1.0D0*TERME2(I,J)+A1*X(I,J)
 304 :                     ENDDO
 305 :                  ENDDO
 306 : 
 307 : C                On calcul le terme C:dfpi/dsigpi
 308 :                  TRA_SE = TERME2(1,1)+TERME2(2,2)+TERME2(3,3)
 309 :                  
 309 :                  CALL DEVIAT(TERME2,TRAV)
 310 : 
 311 :                  DO I=1,3
 312 :                     DO J=1,3
 313 :                        TERME5(I,J) = 0.0D0
 314 :                        TERME5(I,J) = XK*TRA_SE*UNIT(I,J)+
 315 :      &                 2.0D0*XG*TRAV(I,J)
 316 :                     ENDDO
 317 :                  ENDDO
 318 : 
 319 : C                On calcul le denominateur du multiplicateur plastique
 320 :                  XT1 = 0.0D0
 321 :                  DO I=1,3
 322 :                     DO J=1,3
 323 :                        XT1 = XT1 + TERME1(I,J)*TERME5(I,J)
 324 :                     ENDDO
 325 :                  ENDDO
 326 :                  
 327 :                  XT2 = 0.0D0                 
 328 :                  DO I=1,3
 329 :                     DO J=1,3
 330 :                        XT2 = XT2 + TERME3(I,J)*TERME4(I,J)
 331 :                     ENDDO
 332 :                  ENDDO
 333 : 
 334 : C                Calcul du multiplicateur plastique
 335 :                  XLAM = XFPI/(XDT*XT1 + GAMMA1*XT2)
 336 : 
 337 : C                On actualisatise X et sigpi
 338 :                  DO I=1,3
 339 :                     DO J=1,3
 340 :                        SIGPI(I,J)=SIGPI(I,J)-XDT*XLAM*TERME5(I,J)
 341 :                        X(I,J)=X(I,J)-GAMMA1*XLAM*TERME4(I,J)
 342 :                     ENDDO
 343 :                  ENDDO
 344 : 
 345 : C                On recalcule le nouveau seuil
 346 :                  DO I=1,3
 347 :                     DO J=1,3
 348 :                        SIGMOX(I,J) = 0.0D0
 349 :                        SIGMOX(I,J) = SIGPI(I,J)-X(I,J)
 350 :                     ENDDO
 351 :                  ENDDO
 352 : 
 352 :                  CALL I2(SIGMOX,XFPI)
 353 : 
 354 :               ENDIF
 355 : 
 356 :           ENDDO
 357 :           
 358 :  400      CONTINUE 
 359 : 
 360 : 
 361 : C         On calcul epspi
 362 :           TRA_SE=SIGPI(1,1)+SIGPI(2,2)+SIGPI(3,3)
 362 :           CALL DEVIAT(SIGPI,TRAV)
 363 : 
 364 :           TRAPI = (TRA-TRADN) - (TRA_SE/(XDT*3.0D0*XK))
 365 : 
 366 :           DO I=1,3
 367 :              DO J=1,3
 368 :                 EPSPIDEV(I,J) = (EPSDEV(I,J)-DEPSN(I,J))-
 369 :      &          (1.0D0/(2.0D0*XG*XDT))*TRAV(I,J)
 370 :              ENDDO
 371 :           ENDDO
 372 : 
 373 :           DO I=1,3
 374 :              DO J=1,3
 375 :                 EPSPI(I,J) = (TRAPI/3.0D0)*UNIT(I,J)+EPSPIDEV(I,J)
 376 :              ENDDO
 377 :           ENDDO
 378 : 
 379 :       ENDIF
 380 : C
 381 : C-----CONTRAINTES ACTUALISEE
 382 : C
 383 :       TRA1 = SIGMA1(1,1)+SIGMA1(2,2)+SIGMA1(3,3)
 384 :       
 385 :       IF (XBEN.GE.0.0D0) THEN
 386 :          DO I=1,3
 387 :             DO J=1,3
 388 :                SIGMA1(I,J)=XK*(TRA-XDT*TRAPI-TRADN)*UNIT(I,J)+
 389 :      &          2.0D0*XG*
 390 :      &          (EPSDEV(I,J)-XDT*EPSPIDEV(I,J)-DEPSN(I,J))
 391 :             ENDDO
 392 :          ENDDO
 393 :             
 394 :       ELSEIF (XBEN.LE.SIGRF) THEN
 395 :       
 396 :             
 397 : C        La contrainte effective
 398 :          DO I=1,3
 399 :             DO J=1,3
 400 :                SIGMA1(I,J)=XK*(TRA-TRADN)*UNIT(I,J)+
 401 :      &         2.0D0*XG*(EPSDEV(I,J)-DEPSN(I,J))
 402 :             ENDDO
 403 :          ENDDO
 404 : 
 405 : C        On verifie le seuile de plasticite
 405 :          CALL J2(SIGMA1,TEMP1)
 406 : 
 407 :          TEMP2=SIGMA1(1,1)+SIGMA1(2,2)+SIGMA1(3,3)
 408 : 
 409 :          TEMP2 = TEMP2/3.00
 410 :  
 411 :          XFP = XAF*TEMP1+XBF*TEMP2 - (XFC + XECO)
 412 : 
 413 : C        On verifie la positivite du seuil
 414 :          IF (XFP.GT.0.0D0) THEN
 415 :          
 416 : C           On retient la valeur initiale pour le critere relatif
 417 :             XFP0 = XFP
 418 : 
 419 : C           On commence les choses serieuses
 420 :             DO NP=1,IMAX
 421 :             
 422 : C           La derivee du critere
 422 :             CALL DI2(SIGMA1,XDJ2P)
 423 :             
 424 :             DO I=1,3
 425 :             DO J=1,3
 426 :             XDFP(I,J) = XAF*XDJ2P(I,J)+XBF*UNIT(I,J)/3.0D0
 427 :             XDGP(I,J) = XAG*XDJ2P(I,J)+XBG*UNIT(I,J)/3.0D0
 428 :             ENDDO
 429 :             ENDDO
 430 : 
 431 : C           Le hessien elasique             
 431 :             CALL DEVIAT(XDGP,DEVGP)
 432 : 
 433 :             XTGP = XDGP(1,1)+XDGP(2,2)+XDGP(3,3)
 434 :             
 435 :             DO I=1,3
 436 :             DO J=1,3
 437 :             TRAV(I,J)=0.0D0
 438 :             TRAV(I,J)=XK*XTGP*UNIT(I,J)+2.0D0*XG*DEVGP(I,J)
 439 :             ENDDO
 440 :             ENDDO
 441 : 
 442 : C           La derivee de l ecrouissage            
 443 :             XTEMP= DEXP(-1.0D0*XBC*XDPP)*(XAC-(XAC*XDPP+XFC+XSIGU)*XBC)
 444 : 
 445 : C           Le critere actuel 
 446 :            
 446 :             CALL J2(SIGMA1,TEMP1)
 447 : 
 448 :             TEMP2=SIGMA1(1,1)+SIGMA1(2,2)+SIGMA1(3,3)
 449 : 
 450 :             TEMP2 = TEMP2/3.00
 451 :  
 452 :             XFP = XAF*TEMP1+XBF*TEMP2 - (XFC + XECO)
 453 : 
 454 : C           Le multiplicateur plastique
 454 :             CALL DBLECONT(XDFP,TRAV,TEMP9)
 455 : 
 456 :             XLAMP = XFP/(TEMP9+XTEMP)
 457 : 
 458 : C           Actualisation des variables internes
 459 :             XECO = XECO + XLAMP*XTEMP
 460 :             
 461 :             XDPP = XDPP + XLAMP
 462 :             
 463 :             DO I=1,3
 464 :             DO J=1,3
 465 :             DEFN(I,J)=DEFN(I,J) + XLAMP*XDGP(I,J)
 466 :             ENDDO
 467 :             ENDDO
 468 : 
 469 :             DO I=1,3
 470 :             DO J=1,3
 471 :             SIGMA1(I,J)=SIGMA1(I,J) - XLAMP*TRAV(I,J)
 472 :             ENDDO
 473 :             ENDDO
 474 :             
 475 : C           On calcul le seuil pour un eventuel arret
 475 :             CALL J2(SIGMA1,TEMP1)
 476 : 
 477 :             TEMP2=SIGMA1(1,1)+SIGMA1(2,2)+SIGMA1(3,3)
 478 : 
 479 :             TEMP2 = TEMP2/3.00
 480 :  
 481 :             XFP = XAF*TEMP1+XBF*TEMP2 - (XFC + XECO)
 482 :             
 483 : C           Le critere d arret des iterations internes            
 484 :             IF ((DABS(XFP/XFP0).LE.(1.0D-8)).OR.(NP.GE.IMAX)) THEN
 485 :                
 486 :                GOTO 294
 487 :                
 488 :             ENDIF
 489 :             
 490 :             
 491 :             ENDDO
 492 :          
 493 :  294     CONTINUE        
 494 :          
 495 :          ENDIF
 496 :          
 497 :       ELSE
 498 :       
 499 : C        Calcul du deviateur de DEFN
 499 :          CALL DEVIAT(DEFN,DEPSN)
 500 : 
 501 : C        Calcul de la trace de DEFN
 502 :          TRADN = DEFN(1,1)+DEFN(2,2)+DEFN(3,3)
 503 : 
 504 :          TRASIG1 = 3.0D0*XK*(TRA-XDT*TRAPI-TRADN)/
 505 :      &             (1.0D0-(3.0D0*XK/SIGRF)*XDT*TRAPI)
 506 :          
 507 :          DO I=1,3
 508 :             DO J=1,3
 509 :                SIG1DEV(I,J)=2.00D0*XG*(EPSDEV(I,J)-XDT*EPSPIDEV(I,J)*
 510 :      &         (1.0D0-(TRASIG1/SIGRF))-DEPSN(I,J))
 511 :             ENDDO
 512 :          ENDDO
 513 : 
 514 :          DO I=1,3
 515 :             DO J=1,3
 516 :                SIGMA1(I,J)=TRASIG1*UNIT(I,J)/3.0D0+SIG1DEV(I,J)
 517 :             ENDDO
 518 :          ENDDO
 519 : 
 520 :       ENDIF 
 521 :           
 522 : C
 523 : C-----CALCUL DE L ENDOMMAGEMENT EFFECTIF (GRANDEUR INDICATIVE)
 524 : C-----PAS DE SENS THERMO
 525 : C
 526 :       XFSIG = 1.0D0-(SIGMA1(1,1)+SIGMA1(2,2)+SIGMA1(3,3))/SIGRF
 527 :       
 528 :       IF ((XFSIG.LE.1.0D0).AND.(XFSIG.GE.0.0D0)) THEN
 529 :           XDEFF = XDT*XFSIG
 530 :       ELSE
 531 :           XDEFF = XDT
 532 :       ENDIF  
 533 :       
 534 : C
 535 : C-----STOCKAGE DES CONTRAINTES AVEC SYMMETRISATION POUR ROBUSTESSE
 536 : C
 537 :       SIGF(1)= SIGMA1(1,1)
 538 :       SIGF(2)= SIGMA1(2,2)
 539 :       SIGF(3)= SIGMA1(3,3)
 540 :       SIGF(4)= (SIGMA1(1,2)+SIGMA1(2,1))/2.0D0
 541 :       SIGF(5)= (SIGMA1(1,3)+SIGMA1(3,1))/2.0D0
 542 :       SIGF(6)= (SIGMA1(2,3)+SIGMA1(3,2))/2.0D0
 543 : C
 544 : C-----STOCKAGE DES DEFORMATIONS DE GLISSEMENT
 545 : C
 546 :       VARF(2)=EPSPI(1,1)
 547 :       VARF(3)=EPSPI(2,2)
 548 :       VARF(4)=EPSPI(3,3)
 549 :       VARF(5)=EPSPI(1,2)
 550 :       VARF(6)=EPSPI(1,3)
 551 :       VARF(7)=EPSPI(2,3)
 552 : C
 553 : C-----STOCKAGE DE L ECROUISSAGE CINEMATIQUE
 554 : C
 555 :       VARF(8)=X(1,1)
 556 :       VARF(9)=X(2,2)
 557 :       VARF(10)=X(3,3)
 558 :       VARF(11)=X(1,2)
 559 :       VARF(12)=X(1,3)
 560 :       VARF(13)=X(2,3)
 561 : C
 562 : C-----STOCKAGE DES CONTRAINTES TESTS
 563 : C
 564 :       VARF(22)=SIGMA1(1,1)
 565 :       VARF(23)=SIGMA1(2,2)
 566 :       VARF(24)=SIGMA1(3,3)
 567 :       VARF(25)=SIGMA1(1,2)
 568 :       VARF(26)=SIGMA1(1,3)
 569 :       VARF(27)=SIGMA1(2,3)
 570 : C
 571 : C-----TENSEUR DES DEFORMATIONS PLASTIQUES
 572 : C     
 573 :       VARF(35) = DEFN(1,1)
 574 :       VARF(36) = DEFN(2,2)
 575 :       VARF(37) = DEFN(3,3)
 576 :       VARF(38) = DEFN(1,2)
 577 :       VARF(39) = DEFN(1,3)
 578 :       VARF(40) = DEFN(2,3)
 579 : C      
 580 : C-----Ecrouissage plastique
 581 : C
 582 :       VARF(41) = XECO
 583 : C
 584 : C-----Deformation plastique cumulee
 585 : C
 586 :       VARF(42) = XDPP      
 587 : C
 588 : C-----STOCKAGE DES CONTRAINTES TESTS
 589 : C
 590 :       VARF(29)=SIGPI(1,1)
 591 :       VARF(30)=SIGPI(2,2)
 592 :       VARF(31)=SIGPI(3,3)
 593 :       VARF(32)=SIGPI(1,2)
 594 :       VARF(33)=SIGPI(1,3)
 595 :       VARF(34)=SIGPI(2,3)
 596 : C
 597 : C-----ENDOMMAGEMENT EFFECTIF
 598 : C
 599 :       VARF(28) = XDEFF
 600 : C
 601 : C-----STOCKAGE DE L ENDOMMAGEMENT
 602 : C
 603 :       VARF(14) = XDT
 604 : C
 605 : C-----STOCKAGE DES FORCES THERMODYNAMIQUES ASSOCIEES A L ECROUISSAGE ISOTROPE
 606 : C
 607 :       VARF(15) = XFFT
 608 : C
 609 : C-----STOCKAGE DES DEFORMATIONS TOTALES
 610 : C
 611 :       VARF(16)=EPSI(1,1)
 612 :       VARF(17)=EPSI(2,2)
 613 :       VARF(18)=EPSI(3,3)
 614 :       VARF(19)=EPSI(1,2)
 615 :       VARF(20)=EPSI(1,3)
 616 :       VARF(21)=EPSI(2,3)
 617 : C
 618 : C-----BALISE DE SORTIE POUR LE NON LOCAL
 619 : C
 620 :  2000 CONTINUE
 621 : 
 622 :       RETURN
 623 :           END
 624 : 
 625 : 
 626 : 
 627 : 
 628 : 
 629 : 
 630 : 
 631 : 
 632 : 
 633 : 
 634 : 
 635 : 
 636 : 
 637 : 
 638 : 
 639 : 
 640 : 
 641 : 
 642 : 
 643 :  
 644 :  

© Cast3M 2003 - All rights reserved.
Disclaimer