Sources Esope | Cast3M

avan1
C AVAN1     SOURCE    CHAT      05/01/12    21:29:47     5004
C----------------------------------------------------------------
       SUBROUTINE AVAN1(TAU,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,
     +              EPSVPT,VARPT,XPM,DSPT,
     +              XMAT,NSTRS0,NVARI,NCOMAT,IFOUR,IVTEST,
     +              MFR)
C
      IMPLICIT INTEGER(I-N)
       IMPLICIT REAL*8(A-H,O-Z)
       DIMENSION SIG(*),EPSV(*),VAR(*),SIG1(*),EPSV1(*),VAR1(*)
       DIMENSION EPSVPT(*),VARPT(*),XPM(*),DSPT(*),XMAT(*)
       DIMENSION XX(6),EPSV0(6),EPSV00(6)
C----------------------------------------------------------------
C  ====> connaissant l'{tat, c'est @ dire SIG, EPSV et VAR
C        ainsi que les deriv{es EPSVPT, VARPT et DSPT
C        on calcule l'{tat SIG1, EPSV1, VAR1 au temps T + TAU
C----------------------------------------------------------------
C
       NSTRS=NSTRS0
C
C------ AFFECTATION DE QUELQUES VARIABLES
       YOU = XMAT(1)
       XNU = XMAT(2)
       XM  = XMAT(10)
       AN  = XMAT(17)
C
C------- MAJ DE EPSV
       DO 30 I =1,NSTRS
          EPSV1(I) = EPSV(I) + TAU*EPSVPT(I)
  30   CONTINUE
C
C------- MAJ DE SIG
C        SIG1 = SIG + T.(E:EPSPT - E:EPSVPT ) , DSPT = E:EPSPT
C      cas des contraintes planes
       IF (IFOUR.EQ.-2) THEN
          SIG(3) = 0.0
          X2MU = YOU/(1.+XNU)
          XCO  = X2MU/(1.-XNU)
       SIG1(1) =SIG(1) + TAU*(DSPT(1)- XCO*(EPSVPT(1)+XNU*EPSVPT(2)))
       SIG1(2) =SIG(2) + TAU*(DSPT(2)- XCO*(EPSVPT(2)+XNU*EPSVPT(1)))
       SIG1(4) =SIG(4) + TAU*(DSPT(4)- X2MU*EPSVPT(4))
          SIG1(3) = 0.0
C      cas tri dim
       ELSEIF (IFOUR.EQ.2.OR.IFOUR.EQ.0.OR.IFOUR.EQ.-1.OR.
     &         IFOUR.EQ.1.OR.IFOUR.EQ.-3) THEN
       IF (MFR.EQ.5) THEN
         X2MU=YOU/(1.0+XNU)
         XCO = X2MU/(1.0-XNU)
         SIG1(1)=SIG(1)+TAU*(DSPT(1)-XCO*(EPSVPT(1)+XNU*EPSVPT(2)))
         SIG1(2)=SIG(2)+TAU*(DSPT(2)-XCO*(EPSVPT(2)+XNU*EPSVPT(1)))
         SIG1(3)=SIG(3)+TAU*(DSPT(3)-X2MU*EPSVPT(3))
         SIG1(4)=SIG(4)+TAU*(DSPT(4)-X2MU*EPSVPT(4))
         SIG1(5)=SIG(5)+TAU*(DSPT(5)-X2MU*EPSVPT(5))
       ELSE
         TRACE = EPSVPT(1)+EPSVPT(2)+EPSVPT(3)
         X2MU = YOU/(1+XNU)
         XLAMB = X2MU*XNU/(1-2*XNU)
         DO 31 I=1,NSTRS
           A = 0.0
           IF (I.LE.3) A=1.0
       SIG1(I) = SIG(I) + TAU*(DSPT(I)-XLAMB*A*TRACE-X2MU*EPSVPT(I))
  31     CONTINUE
       ENDIF
C= Modes de calcul 1D contraintes planes suivant z (DYCZ et GYCZ)
       ELSE IF (IFOUR.EQ.4.OR.IFOUR.EQ.8) THEN
         SIG(3) = 0.
         XCO = YOU/(1.-XNU*XNU)
         SIG1(1) = SIG(1) + TAU*(DSPT(1)-XCO*(EPSVPT(1)+XNU*EPSVPT(2)))
         SIG1(2) = SIG(2) + TAU*(DSPT(2)-XCO*(EPSVPT(2)+XNU*EPSVPT(1)))
         SIG1(3) = 0.
C= Modes de calcul 1D contraintes planes suivant y (CYDZ et CYGZ) et
C=                 1D axisymetrique contraintes axiales nulles (AXCZ)
       ELSE IF (IFOUR.EQ.5.OR.IFOUR.EQ.10.OR.IFOUR.EQ.13) THEN
         SIG(2) = 0.
         XCO = YOU/(1.-XNU*XNU)
         SIG1(1) = SIG(1) + TAU*(DSPT(1)-XCO*(EPSVPT(1)+XNU*EPSVPT(3)))
         SIG1(2) = 0.
         SIG1(3) = SIG(3) + TAU*(DSPT(3)-XCO*(EPSVPT(3)+XNU*EPSVPT(1)))
C= Modes de calcul 1D contraintes planes suivant y et z (CYCZ)
       ELSE IF (IFOUR.EQ.6) THEN
         SIG(2) = 0.
         SIG(3) = 0.
         SIG1(1) = SIG(1) + TAU*(DSPT(1)-YOU*EPSVPT(1))
         SIG1(2) = 0.
         SIG1(3) = 0.
C= Autres modes de calcul 1D deformations planes (DYDZ GYDZ DYGZ GYGZ)
C=                        1D axisymetrique (AXDZ AXGZ) et 1D spherique
       ELSE IF (IFOUR.EQ. 3.OR.IFOUR.EQ. 7.OR.IFOUR.EQ. 9.OR.IFOUR.EQ.11
     .      .OR.IFOUR.EQ.12.OR.IFOUR.EQ.14.OR.IFOUR.EQ.15) THEN
         X2MU = YOU/(1.+XNU)
         TRACE = EPSVPT(1)+EPSVPT(2)+EPSVPT(3)
         XCO = X2MU*XNU*TRACE/(1.-XNU)
         SIG1(1) = SIG(1) + TAU*(DSPT(1)-X2MU*EPSVPT(1)-XCO)
         SIG1(2) = SIG(2) + TAU*(DSPT(2)-X2MU*EPSVPT(2)-XCO)
         SIG1(3) = SIG(3) + TAU*(DSPT(3)-X2MU*EPSVPT(3)-XCO)
       ELSE
         WRITE(*,*) 'CAS NON PREVU DANS AVAN1 (GUIONNET)'
         STOP
       ENDIF
C
C------ SAUVEGARDE DE VAR, MEMOIRE COMPRISE
       DO 32 I=1,NVARI,1
         VAR1(I) = VAR(I)
  32   CONTINUE
C
C il n'y a que 4+NSTRS valeurs dans VARPT,
C
C------ MAJ DE VAR(1..4)
       IF (MFR.EQ.5) THEN
          NSTRS=NSTRS0+1
          EPSV0(1)=EPSV1(1)
          EPSV0(2)=EPSV1(2)
          EPSV0(3)=-1.*(EPSV0(1)+EPSV0(2))
          EPSV0(4)=EPSV1(3)
          EPSV0(5)=EPSV1(4)
          EPSV0(6)=EPSV1(5)
          EPSV00(1)=EPSV(1)
          EPSV00(2)=EPSV(2)
          EPSV00(3)=-1.*(EPSV00(1)+EPSV00(2))
          EPSV00(4)=EPSV(3)
          EPSV00(5)=EPSV(4)
          EPSV00(6)=EPSV(5)
       ELSE
          NSTRS=NSTRS0
          DO 01 I=1,NSTRS
           EPSV0(I)=EPSV1(I)
           EPSV00(I)=EPSV(I)
 01       CONTINUE
       ENDIF
       VAR1(1) = VAR(1) + TAU*VARPT(1)
       VAR1(2) = VAR(2) + TAU*VARPT(2)
C      Chargement cyclique ?
       IF (VAR(2).LT.VAR(5+NSTRS)) THEN
          VAR1(3) = VAR(3) + TAU*VARPT(3)
          VAR1(4) = VAR(4) + TAU*VARPT(4)
       ELSE
         DPN = VAR1(2)**AN - VAR(2)**AN
         VAR1(3) = VAR(3) + DPN
         VAR1(4) = VAR(4) + DPN
       ENDIF
C
C------ MAJ DE VAR(4+1..NSTRS)
C         D X = XPM  * D P**XM + XPT * D T
C         P1 = VAR(2)
       DPM = VAR1(2)**XM - VAR(2)**XM
       DO 34 I=1,NSTRS,1
         VAR1(4+I) = VAR(4+I) + XPM(I)*DPM + VARPT(4+I)*TAU
  34   CONTINUE
C
C     CRITERE D'INVERSION : EXTREMUM D'UNE DES COMPOSANTES DU TENSEUR
C                           DES DEFORMATIONS VISCOPLASTIQUES
C
      CALL GIOCRT(EPSV0,EPSV00,NSTRS,VAR1,NVARI,IVTEST)
C
       RETURN
       END