C ADVAN1    SOURCE    OF166741  25/11/04    21:15:02     12349          

C=======================================================================
C=                                                                     =
C=  THIS ROUTINE UPDATES THE STRAINS ,STRESSES & INTERNAL VARIABLES.   =
C=                                                                     =
C=======================================================================

      SUBROUTINE ADVAN1 (TAU,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,DSPT,EPSVPT,
     .                   VARPT,XMAT,NSTRS,NVARI,IFOUR,INPLAS,MFR)

      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)

      DIMENSION SIG(*),EPSV(*),VAR(*),SIG1(*),EPSV1(*),VAR1(*),
     .          DSPT(*),EPSVPT(*),VARPT(*),XMAT(*)

      DATA UN / 1.D0 /

C  STEP-1 : UPDATING STRAINS
C ----------------------------
      YOU=XMAT(1)
      XNU=XMAT(2)
      IF (INPLAS.EQ.107) THEN
C= EPSVPT contient les vit de def viscopl. + densification + gonflement
         A=VARPT(5)+VARPT(6)
         DO 10 I=1,NSTRS
           IF (I.GT.3) A=0.D0
           EPSV1(I)=EPSV(I)+TAU*(EPSVPT(I)-A)
  10     CONTINUE
      ELSE
         DO 15 I=1,NSTRS
           EPSV1(I)=EPSV(I)+TAU*EPSVPT(I)
  15     CONTINUE
      ENDIF
C  STEP-2 : UPDATING STESSES
C ----------------------------
C= Mode de calcul 2D contraintes planes
      IF (IFOUR.EQ.-2) THEN
        SIG(3)=0.D0
        X2MU = YOU/(UN+XNU)
        XCO = X2MU/(UN-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.D0
        SIG1(4) = SIG(4) + TAU*(DSPT(4)-X2MU*EPSVPT(4))
C= Modes de calcul 3D, 2D deformations planes (GENE ou non),
C=                 2D axisymetrique, 2D series de Fourier
      ELSE IF (IFOUR.EQ.2.OR.IFOUR.EQ.0.OR.IFOUR.EQ.1.OR.
     .         IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
        X2MU = YOU/(UN+XNU)
        ALMB = X2MU*XNU/(UN-2*XNU)
        IF (MFR.EQ.5) THEN
          X2MU = YOU/(UN+XNU)
          XCO = X2MU/(UN-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)
          CC=ALMB*TRACE
          DO 20 I=1,3
            SIG1(I) = SIG(I) + TAU*(DSPT(I)-CC-X2MU*EPSVPT(I))
  20      CONTINUE
          DO 21 I=4,NSTRS
            SIG1(I) = SIG(I) + TAU*(DSPT(I)-X2MU*EPSVPT(I))
  21      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.D0
         XCO = YOU/(UN-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.D0
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.D0
         XCO = YOU/(UN-XNU*XNU)
         SIG1(1) = SIG(1) + TAU*(DSPT(1)-XCO*(EPSVPT(1)+XNU*EPSVPT(3)))
         SIG1(2) = 0.D0
         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.D0
         SIG(3) = 0.D0
         SIG1(1) = SIG(1) + TAU*(DSPT(1)-YOU*EPSVPT(1))
         SIG1(2) = 0.D0
         SIG1(3)=0.D0
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/(UN+XNU)
         TRACE = EPSVPT(1)+EPSVPT(2)+EPSVPT(3)
         XCO = X2MU*XNU*TRACE/(UN-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)
       ENDIF

C  STEP-3 : UPDATING THE INTERNAL VARIABLES
C------------------------------------------
       IF (MFR.EQ.5) THEN
         NSTRS0=NSTRS+1
       ELSE
         NSTRS0=NSTRS
       ENDIF
       IF ((INPLAS.GE.19.AND.INPLAS.LE.23).OR.INPLAS.EQ.61) THEN
         DO 30 I=1,2*NSTRS0+3
           VAR1(I)=VARPT(I)
  30     CONTINUE
         NSUB=2*NSTRS0+4
         VAR1(NSUB)=VAR(NSUB)+TAU*VARPT(NSUB)
         IF (INPLAS.EQ.20.OR.INPLAS.EQ.22.OR.INPLAS.EQ.23.OR.
     .       INPLAS.EQ.61) THEN
           DO 40 I=1,NSTRS0
             VAR1(NSUB+I)=VAR(NSUB+I)+TAU*VARPT(NSUB+I)
  40       CONTINUE
         ENDIF
       ELSE IF (INPLAS.EQ.24) THEN
         VAR1(1) = VAR(1) + TAU*VARPT(1)
         VAR1(2) = VAR(2) + TAU*VARPT(2)
       ELSE IF (INPLAS.EQ.25.OR.INPLAS.EQ.53.OR.INPLAS.EQ.70.OR.
     .          INPLAS.EQ.76.OR.INPLAS.EQ.77.OR.INPLAS.EQ.165.OR.
     .          INPLAS.EQ.192) THEN
         DO 50 I=1,NVARI
           VAR1(I)=VAR(I)+TAU*VARPT(I)
  50     CONTINUE
       ELSE IF (INPLAS.EQ.107) THEN
         DO 60 I=1,NVARI
           VAR1(I)=VAR(I)+TAU*VARPT(I)
  60     CONTINUE
         IF (VAR1(2).LT.0.D0) THEN
           VAR1(2) = 0.D0
         ELSE IF (VAR1(2).GT.0.3D0) THEN
           VAR1(2) = 0.3D0
         ENDIF
       ELSE IF(INPLAS.EQ.130) THEN
         DO 70 I=1,NVARI
           VAR1(I)=VAR(I)+TAU*VARPT(I)
  70     CONTINUE
         VAR1(3) = MIN( VAR1(3), 0.98D0 )
         VAR1(4) = MIN( VAR1(4), 0.98D0 )
       ENDIF

       RETURN
       END