Télécharger umagat.eso

Retour à la liste

umagat
C UMAGAT    SOURCE    STRU      07/05/31    21:15:29     5744
      SUBROUTINE UMAGAT ( STRESS, STATEV, DDSDDE,
     &                  STRAN, DSTRAN, DTIME,
     &                  TEMP, DTEMP, PREDEF, DPRED,
     &                  NTENS, NSTATV,
     &                  PROPS, NPROPS, COORDS,
     &                  NOEL, NPT, KINC )
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
C-----------------------------------------------------------------------
C
C     DESCRIPTION FONCTIONNELLE :
C     -------------------------
C
C     Dependance du module Utilisateur UMAT
C     Appelee dans le cas suivant :
C
C     ILOI =  7 : Modele 'VISCOPLASTIQUE' 'GATT_MONERIE' de CAST3M,
C                 INPLAS=107
C                 RESTRICTION A LA FORMULATION MASSIVE
C
C     Composantes de materiau
C     -----------------------
C     PROPS ( 1) : 'YOUN'
C     PROPS ( 2) : 'NU  '
C     PROPS ( 3) : 'RHO '
C     PROPS ( 4) : 'ALPH'
C     PROPS ( 5) : 'R   '
C     PROPS ( 6) : 'DG  '
C     PROPS ( 7) : 'K1  '
C     PROPS ( 8) : 'M1  '
C     PROPS ( 9) : 'Q1  '
C     PROPS (10) : 'N1  '
C     PROPS (11) : 'K2  '
C     PROPS (12) : 'M2  '
C     PROPS (13) : 'Q2  '
C     PROPS (14) : 'N2  '
C     PROPS (15) : 'OMEG'
C     PROPS (16) : 'Q   '
C     PROPS (17) : 'H   '
C     PROPS (18) : 'BETA'
C     PROPS (19) : 'K   '
C     PROPS (20) : 'A   '
C     PROPS (21) : 'Q3  '
C     PROPS (22) : 'N3  '
C     PROPS (23) : 'ADEN'
C     PROPS (24) : 'KGON'
C     PROPS (25) : 'POR0'
C     PROPS (26) : 'BUMI'
C     PROPS (27) : 'EFIS'
C     PROPS (28) : 'DG0 '
C     PROPS (29) : 'CR1 '
C     PROPS (30) : 'CR2 '
C     PROPS (31) : 'CR3 '
C     PROPS (32) : 'CR  '
C     PROPS (33) : 'KP  '
C     PROPS (34) : 'AP  '
C     PROPS (35) : 'BP  '
C     PROPS (36) : 'QP  '
C     PROPS (37) : 'VISQ' (non utilise mais c'est la premiere comp. fac. pour gatt.eso)
C     PROPS (38) : 'TYPE'
C     PROPS (39) : 'COMP'
C     PROPS (40) : 'DYN '
C
C     Variables internes
C     ------------------
C     STATEV( 1) : 'EPSE'
C     STATEV( 2) : 'PORO'
C     STATEV( 3) : 'BU  '
C     STATEV( 4) : 'TETA'
C     STATEV( 5) : 'EPSD'
C     STATEV( 6) : 'EPSG'
C
C     STATEV( 7) : 'ENXX'
C     STATEV( 8) : 'ENYY'
C     STATEV( 9) : 'ENZZ'
C     STATEV(10) : 'GNXY'
C     STATEV(11) : 'GNXZ'
C     STATEV(12) : 'GNYZ'
C
C     Parametres externes
C     -------------------
C     TEMP       : 'T   '
C     PREDEF(1)  : 'DFIS'
C
C-----------------------------------------------------------------------
C     Arguments de l'interface
C
      INTEGER       NTENS, NSTATV, NPROPS,
     &              NOEL, NPT, KINC
C
      REAL*8        STRESS(NTENS), STATEV(*),
     &              DDSDDE(NTENS,NTENS),
     &              STRAN(NTENS), DSTRAN(NTENS),
     &              DTIME,
     &              TEMP, DTEMP, PREDEF(*), DPRED(*),
     &              PROPS(NPROPS),
     &              COORDS(3)
C
C     Variables locales
C
      INTEGER    NVARI
      PARAMETER (NVARI=6)
C
C     Tableaux de travail
C
      REAL*8     VAR(NVARI),VAR1(NVARI),VARP1(NVARI),VARP2(NVARI),
     &           VAR12(NVARI),VARP3(NVARI),VARP4(NVARI),VAR13(NVARI),
     &           SIG(8),SIG1(8),SIG12(8),SIG13(8),
     &           DSPT(8),XX(8),
     &           EVP1(8),EVP2(8),EVP3(8),EVP4(8),
     &           EPSV(8),EPSV1(8),EPSV12(8),EPSV13(8)
C
      REAL*8     EPIN0(8),DSIGT(8),VARF(NVARI),SIGF(8),EPINF(8)
C
C     Variables scalaires
C
      LOGICAL    DTLIBR
      REAL*8     DT, DTLEFT, BORNE, RMAX, RMIN, DIV, FAC,
     &           T0, TF, FI0, FIF,
     &           YOU, XNU, AUX0, AUX, AUX1, GEGE,
     &           XMAX,
     &           TI0, TI1, TPOINT, DELTAT, TI12,
     &           FII0, FII1, FPOINT, DELTAF, FII12,
     &           VARPBU,
     &           TAU, ASIG, ERRABS, TAU2, RA, SQRA, RA1
      INTEGER    MSOUPA, I, J, ITERO, NSSINC, NITERA,
     &           INPLAS, MFR1, IFOURB
C
      REAL*8     UN, DEUX, UNDEMI
      PARAMETER (UN=1.D0,DEUX=2.D0,UNDEMI=.5D0)
C
      REAL*8     PRECIS
C
      REAL*8     PROCON
      EXTERNAL   PROCON
C
C------------------- Debut du code executable --------------------------
C
C=======================================================================
C 1 - INITIALISATIONS
C     Pas de temps, parametres de pilotage des iterations internes,
C     parametres du modele
C=======================================================================
C
      PRECIS=1.d-8
      DTLIBR=.TRUE.
C
      DT = DTIME
      IF (DT.LT.0.0) THEN
         KINC = -71
         RETURN
      ENDIF
      IF (DT.EQ.0.0) DT = 1.e-20
C
      DTLEFT = DT
      BORNE  = 2.0
      RMAX   = 1.3
      RMIN   = 0.7
      DIV    = 7.0
      FAC    = 3.0
C
      MSOUPA = 200
C
      INPLAS = 107
      MFR1   = 1
      IFOURB = 2
C
      T0  = TEMP
      TF  = T0 + DTEMP
      FI0 = PREDEF(1)
      FIF = FI0 + DPRED(1)
C
C=======================================================================
C 2 - Prediction elastique de l'increment de contraintes
C     La sequence ci-dessous se substitue a l'appel a CALSIG
C     Le calcul des coefficients de la matrice de Hooke reprend la
C     sequence d'operations de HOOKIS
C=======================================================================
C
      YOU = PROPS(1)
      XNU = PROPS(2)
C
      AUX0 = UN/((UN+XNU)*(UN-DEUX*XNU))
      AUX  = YOU*AUX0*(UN-XNU)
      AUX1 = YOU*AUX0*XNU
      GEGE = YOU*UNDEMI/(UN+XNU)
C
      DO 10 J=1,NTENS
         DO 11 I=1,NTENS
            DDSDDE(I,J)=0.0D0
  11     CONTINUE
  10  CONTINUE
C
      DDSDDE(1,1) = AUX
      DDSDDE(1,2) = AUX1
      DDSDDE(1,3) = AUX1
C
      DDSDDE(2,1) = AUX1
      DDSDDE(2,2) = AUX
      DDSDDE(2,3) = AUX1
C
      DDSDDE(3,1) = AUX1
      DDSDDE(3,2) = AUX1
      DDSDDE(3,3) = AUX
C
      DDSDDE(4,4) = GEGE
      DDSDDE(5,5) = GEGE
      DDSDDE(6,6) = GEGE
C
      DO 20 I=1,NTENS
         DSIGT(I)=0.0D0
  20  CONTINUE
      DO 21 J=1,NTENS
         DO 22 I=1,NTENS
            DSIGT(I)=DSIGT(I)+DDSDDE(I,J)*DSTRAN(J)
  22     CONTINUE
  21  CONTINUE
C
C=======================================================================
C 3 - Pre-traitement des deformations inelastiques au debut du pas
C     Passage gamma -> epsilon pour les termes extradiagonaux
C=======================================================================
C
      DO 34 J=1,3
         EPIN0(J)=STATEV(NVARI+J)
  34  CONTINUE
      DO 35 J=4,6
         EPIN0(J)=0.5D0*STATEV(NVARI+J)
  35  CONTINUE
C
C=======================================================================
C 4 - AFFECTATION DE XMAX
C=======================================================================
C
      XMAX=PROPS(1)*1.D-3
C
C=======================================================================
C 5 - INITIALISATIONS AVANT ITERATIONS EN SOUS-INCREMENTS
C=======================================================================
C
      ITERO = 0
 6543 CONTINUE
C
      ITERO = 1 + ITERO
      IF ( ITERO.NE.1) THEN
         DTLIBR = .TRUE.
         PRECIS = PRECIS * 7.D0
         IF (ITERO.GT.3) THEN
            KINC = -72
            RETURN
         ENDIF
      ENDIF
C
      DTLEFT = DT
      TAU    = DTLEFT
      TI0=T0
      TI1=TF
      TPOINT=(TF-T0)/DT
      FII0=FI0
      FII1=FIF
      FPOINT=(FIF-FI0)/DT
      ASIG = SQRT(PROCON(STRESS,STRESS,NTENS))
      ERRABS = PRECIS*ASIG
      IF (XMAX.GT.ASIG) ERRABS = PRECIS*XMAX
C
      DO 40 I=1,NTENS
         SIG(I) = STRESS(I)
         EPSV(I) = EPIN0(I)
         DSPT(I) = DSIGT(I)/DT
  40  CONTINUE
C
      DO 50 I=1,NVARI
         VAR(I)=STATEV(I)
  50  CONTINUE
C
C=======================================================================
C 6 - ITERATIONS EN SSINCREMENTS  /FIN SI DTLEFT = 0
C=======================================================================
C
      NSSINC = 0
      NITERA = 0
  60  CONTINUE
C
         NSSINC = NSSINC + 1
         IF (NSSINC.GT.MSOUPA) THEN
            DTLIBR=.FALSE.
            GOTO 6543
         ENDIF
C
         DELTAT=TPOINT*TAU
         TI1=TI0+DELTAT
         DELTAF=FPOINT*TAU
         FII1=FII0+DELTAF
         CALL INCRE1(TAU,SIG,EPSV,VAR,EVP1,VARP1,PROPS,NTENS,NVARI,
     &               INPLAS,NPROPS,MFR1,FII0,FII1,TI1)
C
C        ---------------------------------------------------------------
C        DEBUT DES ITERATIONS SUR TAU OPTIMAL /FIN SI RA PETIT
C
         NITERA  = 0
  70     CONTINUE
C
            NITERA  = NITERA + 1
C
C           ____________________________________________________________
            TAU2=TAU*0.5D0
            CALL ADVAN1(TAU2,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,DSPT,EVP1,
     &                  VARP1,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
C
C           ____________________________________________________________
            DELTAT=TPOINT*TAU
            TI1=TI0+DELTAT
            DELTAF=FPOINT*TAU
            FII1=FII0+DELTAF
            VAR1(3)=VAR(3)
            CALL INCRE1(TAU,SIG1,EPSV1,VAR1,EVP2,VARP2,PROPS,NTENS,
     &                  NVARI,INPLAS,NPROPS,MFR1,FII0,FII1,TI1)
            DO 701 I=1,NTENS
               EVP2(I) = 0.5D0*(EVP1(I)+EVP2(I))
 701        CONTINUE
            VARPBU=VARP2(3)
            DO 702 I=1,NVARI
               VARP2(I)= 0.5D0*(VARP1(I) + VARP2(I))
 702        CONTINUE
C
C           ____________________________________________________________
            DELTAT=TPOINT*TAU2
            TI12=TI0+DELTAT
            DELTAF=FPOINT*TAU2
            FII12=FII0+DELTAF
            CALL ADVAN1(TAU2,SIG,EPSV,VAR,SIG12,EPSV12,VAR12,DSPT,EVP2,
     &                  VARP2,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
            VAR12(3)=VAR(3)
            CALL INCRE1(TAU2,SIG12,EPSV12,VAR12,EVP3,VARP3,PROPS,NTENS,
     &                  NVARI,INPLAS,NPROPS,MFR1,FII0,FII12,TI12)
C
C           ____________________________________________________________
            TI1=TI12+DELTAT
            FII1=FII12+DELTAF
            CALL ADVAN1(TAU2,SIG12,EPSV12,VAR12,SIG13,EPSV13,VAR13,DSPT,
     &                  EVP3,VARP3,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
            CALL INCRE1(TAU2,SIG13,EPSV13,VAR13,EVP4,VARP4,PROPS,NTENS,
     &                  NVARI,INPLAS,NPROPS,MFR1,FII12,FII1,TI1)
            DO 711 I=1,NTENS
               EVP4(I) = 0.5D0*(EVP3(I)+EVP4(I))
 711        CONTINUE
            DO 712 I=1,NVARI
               VARP4(I)= 0.5D0*(VARP3(I) + VARP4(I))
 712        CONTINUE
            CALL ADVAN1(TAU2,SIG12,EPSV12,VAR12,SIGF,EPINF,VARF,DSPT,
     &                  EVP4,VARP4,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
C
C           ____________________________________________________________
            DELTAT=TPOINT*TAU2
            TI12=TI0+DELTAT
            DELTAF=FPOINT*TAU2
            FII12=FII0+DELTAF
            VARF(3)=VAR(3)
            CALL INCRE1(TAU2,SIGF,EPINF,VARF,EVP4,VARP4,PROPS,NTENS,
     &                  NVARI,INPLAS,NPROPS,MFR1,FII0,FII12,TI12)
            DO 721 I=1,NTENS
               EVP2(I) = (EVP1(I)+EVP4(I))/6.D0+EVP3(I)*2.D0/3.D0
 721        CONTINUE
            DO 722 I=1,NVARI
               VARP2(I)=(VARP1(I) + VARP4(I))/6.D0+VARP3(I)*2.D0/3.D0
 722        CONTINUE
            VARP2(3)=VARPBU
            CALL ADVAN1(TAU,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,DSPT,EVP2,
     &                  VARP2,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
            VARF(3)=VAR1(3)
C
C           ____________________________________________________________
C           CALCUL DU RAPPORT : ERREUR CALCULEE / ERREUR ADMISE
C
            DO 730 I=1,NTENS
               XX(I) = SIGF(I)-SIG1(I)
 730        CONTINUE
            RA = SQRT(PROCON(XX,XX,NTENS))/(ERRABS)
            SQRA = SQRT(RA)
C
C           ____________________________________________________________
C           TEST DE FIN D'ITERATIONS /  MISE A JOUR DE TAU
C           DIV =7   BORNE = 2
C           SI   SQRA>7   TAU = TAU/7    ET NOUVEL ESSAI
C           SI 2<RA<7*7   ON VISE RA = 1  ET NOUVEL ESSAI
C
            IF (DTLIBR) THEN
C
C               Petite ruse pour dejouer l'optimisation
                RA1=RA*1.D0
C
                IF ((RA.GT.DIV*DIV).OR.(RA.NE.RA1)) THEN
                   TAU = TAU/DIV
                   DELTAT=TPOINT*TAU
                   TI1=TI0+DELTAT
                  GOTO 70
                ELSE IF ( RA.GT.(BORNE)) THEN
                   TAU = TAU/SQRA
                   DELTAT=TPOINT*TAU
                   TI1=TI0+DELTAT
                  GOTO 70
               ENDIF
            ENDIF
C           ____________________________________________________________
C
C        FIN D'ITERATIONS SUR TAU OPTIMAL / MISE A JOUR DES VARIABLES
C        Ici RA < BORNE
C        On avance en temps
C
         DO 80 I=1,NTENS,1
            SIG(I)  = SIGF(I)
            EPSV(I) = EPINF(I)
  80     CONTINUE
         DO 90 I=1,NVARI,1
            VAR(I)  = VARF(I)
  90     CONTINUE
C        ---------------------------------------------------------------
C
C     TEST DE FIN DES ITERATIONS EN SSINCREMENTS / MISE A JOUR DE TAU
C     * si   SQRA<1/3        TAU = TAU*3
C     * si   1/3<SQRA<RMIN   on vise RA = 1
C     * si   RMIN<SQRA<RMAX  TAU inchang{
C     * si   SQRA>RMAX       on vise RA = 1
C     Fin des boucles en ss increments  si tau = DTLEFT
C
      IF ( TAU.LT.DTLEFT ) THEN
         DTLEFT = DTLEFT - TAU
         IF ( FAC*SQRA.LT.1.D0) THEN
            TAU=TAU*FAC
         ELSE IF ( (SQRA.LT.RMIN).OR.(SQRA.GT.RMAX) ) THEN
            TAU=TAU/SQRA
         ENDIF
         IF (TAU.GT.DTLEFT) TAU = DTLEFT
         TI0=TI1
         DELTAT=TPOINT*TAU
         TI1=TI0+DELTAT
         FII0=FII1
         DELTAF=FPOINT*TAU
         FII1=FII0+DELTAF
         GOTO 60
      ENDIF
C
      IF (ABS(TAU-DTLEFT).GT.(TAU/1000.)) THEN
         KINC = -73
         RETURN
      ENDIF
C
C=======================================================================
C 7 - STOCKAGE DES RESULTATS AVANT SORTIE
C=======================================================================
C
C     Contraintes
C
      DO 100 J=1,NTENS
         STRESS(J)=SIGF(J)
 100  CONTINUE
C
C     Variables internes
C
      DO 110 J=1,NVARI
         STATEV(J)=VARF(J)
 110  CONTINUE
C
C     Deformations inelastiques
C     Passage epsilon -> gamma pour les termes extradiagonaux
C
      DO 120 J=1,3
         STATEV(NVARI+J)=EPINF(J)
 120  CONTINUE
      DO 121 J=4,6
         STATEV(NVARI+J)=2.0D0*EPINF(J)
 121  CONTINUE
C
      RETURN
      END
 
 

© Cast3M 2003 - Tous droits réservés.
Mentions légales