Télécharger coul1.eso

Retour à la liste

coul1
C COUL1     SOURCE    JK148537  23/03/27    21:15:03     11641          
      SUBROUTINE COUL1 ( IB,IGAU,NSTRS,SIG00,EPST00,EPIN00,VAR0,NVARI,
     .     DEPST0,IFOUR,XMAT,NMATT,IVAL,DD,SIGF1,DEFP1,VARF,KERRE)
C-----------------------------------------------------------------------
C
C         PLASTICITE MODELE 3D MOHR COULOMB ELEMENTS JOINTS JOI1
C
C-----------------------------------------------------------------------
C
C         PLASTICITE MODELE 3D MOHR COULOMB ELEMENTS JOINTS JOI1
C
C ENTREES
C     NSTRS           = NOMBRE DE COMPOSANTES DE CONTRAINTES
C     SIG00           = CONTRAINTES INITIALES (AU DEBUT DU PAS)
C     VAR0(NVARI)     = VARIABLES INTERNES    (AU DEBUT DU PAS)
C       VAR0(1)       = 'EPSE' DEFORMATION PLASTIQUE EQUIVALENTE
C       VAR0(2)       = 'EPOU' DEF. PLASTIQUE DUE A L'OUVERTURE SEULE
C       VAR0(3)       = 'STAT' ETAT DU JOINT
C       VAR0(4)       = 'LAM1' MULTIPLICATEUR PLASTIQUE EN TRACTION
C       VAR0(5)       = 'PLA1' DEF. PLAST. DS. LA DIR. DE GLISSEMENT 1
C       VAR0(6)       = 'PLA2' DEF. PLAST. DS. LA DIR. DE GLISSEMENT 2
C       VAR0(7)       = 'PLA3' DEF. PLAST. DS. LA DIR. NORMALE AUX
C                              DIRECTIONS DE GLISSEMENT
C       VAR0(8)       = 'EPCO' DEF. PLASTIQUE DUE A LA COMPRESSION SEULE
C       VAR0(9)       = 'LAM3' MULTIPLICATEUR PLASTIQUE EN COMPRESSION
C
C     NMATT           = NOMBRE DE COMPOSANTES DE MATERIAU
C     XMAT(NMATT)     = COMPOSANTES DU MATERIAU
C     XMAT ( 1)       = 'V1X '                          V1=(V1X,V1Y,V1Z)
C     XMAT ( 2)       = 'V1Y '
C     XMAT ( 3)       = 'V1Z '
C     XMAT ( 4)       = 'V2X '                          V2=(V2X,V2Y,V2Z)
C     XMAT ( 5)       = 'V2Y '
C     XMAT ( 6)       = 'V2Z '                                  V3=V1^V2
C     XMAT ( 7)       = 'KN  ' RAIDEUR EN TRANSLATION DS LA DIRECTION V1
C     XMAT ( 8)       = 'KS1 ' RAIDEUR EN TRANSLATION DS LA DIRECTION V2
C     XMAT ( 9)       = 'KS2 ' RAIDEUR EN TRANSLATION DS LA DIRECTION V3
C     XMAT (10)       = 'QN  ' RAIDEUR EN ROTATION AUTOUR DE V1
C     XMAT (11)       = 'QS1 ' RAIDEUR EN ROTATION AUTOUR DE V2
C     XMAT (12)       = 'QS2 ' RAIDEUR EN ROTATION AUTOUR DE V3
C     XMAT (13)       = 'QT  ' RAIDEUR TANGENTE NORMALE AU PLAN DE GLISS.
C     XMAT (14)       = 'FNE ' LIMITE ELASTIQUE DE L'EFFORT NORMAL
C     XMAT (15)       = 'COHE' COHESION
C     XMAT (16)       = 'FRIC' ANGLE DU CONE EN °
C     XMAT (17)       = 'TYPE' NORMALE AU PLAN DE GLISSEMENT
C
C     EPST00(NSTRS)   = DEFORMATION TOTALE INITIALE (AU DEBUT DU PAS)
C     DEPST(NSTRS)    = INCREMENT DE DEFORMATION TOTALE
C     DD(NSTRS,NSTRS) = MATRICE DE HOOK
C
C SORTIES
C     STRESS(NSTRS)   = CONTRAINTES FINALES
C     STATEV(NSTATV)  = VARIABLES INTERNES FINALES
C     KINC           = INDICE D'ERREUR
C                     = 60 PAS DE CONVERGENCE APRES 10 ITERATIONS
C                          INTERNES (VARIABLE ISSU DANS COULM)
C                     = 61
C
C-----------------------------------------------------------------------
C
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
-INC CCREEL
C
      DIMENSION SIG00(*),EPST00(*),EPIN00(*),VAR0(*),DEPST0(*),XMAT(*)
      DIMENSION DD(NSTRS,NSTRS),IDI1(6),IDI2(6),TRAID(6,6)
      DIMENSION SIGF1(*),DEFP1(*),VARF(*)
C    Mes arguments
      DIMENSION SIG0(3),EPST0(3),EPIN0(3),DEPST(3)
      DIMENSION SIGF(3),DEFP(3)
C     ------------------------------------------------------------------
      DATA IDI1 / 2 , 3 , 1 , 5 , 6 , 4 /
      DATA IDI2 / 3 , 1 , 2 , 6 , 4 , 5 /
C     ------------------------------------------------------------------
C
C -------IDEBMA POINTE SUR LA FIN DES 2 VECTEURS DU REPERE DE L'EF JOI1
      IDEBMA=6
C
C------CHARGEMENT DE TOUTES LES CARACTERISTIQUES MECANIQUES
      CALL ZERO(TRAID,6,6)
      DO 100 I=1,NSTRS
        TRAID(I,I)=XMAT(IDEBMA+I)
 100  CONTINUE
C ----------DIRECTION NORMALE AU PLAN DE GLISSEMENT
      IF (INT(XMAT(IDEBMA+NSTRS+5)).EQ.1) THEN
      IDIRN=1
      ELSE IF (INT(XMAT(IDEBMA+NSTRS+5)).EQ.2) THEN
      IDIRN=4
      ENDIF
C
C ----------DIRECTIONS DE GLISSEMENT
      IDIR1=IDI1(IDIRN)
      IDIR2=IDI2(IDIRN)
C
C ----------CONTRAINTES AU DEBUT DU PAS POUR LES 3 DIRECTIONS
      SIG0(1)=SIG00(IDIR1)
      SIG0(2)=SIG00(IDIR2)
      SIG0(3)=SIG00(IDIRN)
C
C ----------DEF. PLASTIQUES AU DEBUT DU PAS POUR LES 3 DIRECTIONS
      IDEBLA=4
      DO I=1,3
        EPIN0(I)=VAR0(IDEBLA+I)
      END DO
C
C  ----------INCREMENT DES DEFORMATIONS TOTALES POUR LES 3 DIRECTIONS
      DEPST(1)=DEPST0(IDIR1)
      DEPST(2)=DEPST0(IDIR2)
      DEPST(3)=DEPST0(IDIRN)
C
C ----------PARAMETRES POUR LE TEST DE NULLITE DES CONTRAINTES ET
C ----------DEFORMATIONS
C
      PRESP = +1.D-5
      PRESM = -1.D-5
      PREEP = +1.D-10
*     PREEP = 0.D0
      PREEM = -1.D-10
C
C
C ----------PARAMETRE D'IMPPRESSION
C
      IMPRI = 0
C
C ----------NOMBRE MAXIMUM D'ITERATIONS INTERNES AUTORISEES
C
      MAXITE = 1000
C
C
C------CHARGEMENT DU EPSILON DE PRECISION DES ITERATIONS INTERNES
C
      PRECIS = 1.D-6
      PRECIN = -1.D-6
C
C------CHARGEMENT DES CARACTERISTIQUES MECANIQUES
C------DANS LA DIECTION NORMALE ET LE PLAN DE GLISSEMENT
C
      EC     = XMAT(IDEBMA+IDIRN)
      G1     = XMAT(IDEBMA+IDIR1)
      G2     = XMAT(IDEBMA+IDIR2)
*- cs deb
*      EF     = XMAT(IDEBMA+NSTRS+1)
      QT     = XMAT(IDEBMA+NSTRS+1)
      FNE    = XMAT(IDEBMA+NSTRS+2)
      IF (ABS(EC-QT).LT.(XZPREC*MAX(EC,QT))) THEN
        write(6,*) 'coul1.eso : EC,QT =',EC,QT
        CALL ERREUR(26)
        RETURN
      ELSE
        HH     = EC*QT/(EC-QT)
      ENDIF 
*      ECN    = XMAT(IDEBMA+NSTRS+2)
*- cs fin
      C      = XMAT(IDEBMA+NSTRS+3)
      PHI    = XMAT(IDEBMA+NSTRS+4)
      PHI    = (PHI*XPI)/180.D0
      ZMU    = 0.0D0
      ZMU    = (ZMU*XPI)/180.D0
*      FTRC   = XMAT(IDEBMA+NSTRS+15)
      FTRC   = 0.
C  AM 28/7/95
C  FTRC EST EN 4 ET NON EN 9 A CAUSE DE RHO ET ALPH
C  VOIR COMMENTAIRE DANS ECOUL2
C
C   LA LIMITE EN TRACTION EST MISE BRUTALEMENT A 0
C   DES QU'ELLE EST DEPASSEE.
C
C
C------CHARGEMENT DE LA MATRICE DE HOOK INVERSE
C------POUR LE CALCUL DE L'INCREMENT DE DEFORMATION
C------A PARTIR DE L'INCREMENT DE CONTRAINTE
C------         DEPSN = D**-1 DSIGT
C                 ( D**-1 = DD )
C
      DD(1,1) = 1.0D0/G1
      DD(1,2) = 0.0D0
      DD(1,3) = 0.0D0
      DD(2,1) = 0.0D0
      DD(2,2) = 1.0D0/G2
      DD(2,3) = 0.0D0
      DD(3,1) = 0.0D0
      DD(3,2) = 0.0D0
      DD(3,3) = 1.0D0/EC
C
C-----------INCREMENT DE DEFORMATION
C
      DEPS1 = DEPST(1)
      DEPS2 = DEPST(2)
      DEPSN = DEPST(3)
C
C             TRAITEMENT DU PASSAGE DIRECT ENTRE UN ETAT OUVERT
C             (I.E. EXISTENCE D'UN JEU AU TEMPS N) A UN ETAT
C             COMPRIME AU TEMPS N+1
C
      DEFPN=0.0D0
      INDIC=0
      IFIN =0
      STOTA1 = 0.0D0
      STOTA2 = 0.0D0
      DEFP11 = 0.0D0
      DEFP12 = 0.0D0
      DEFP21 = 0.0D0
      DEFP22 = 0.0D0
      ALP    = 1.0D0
      S1INT  = 0.0D0
      S2INT  = 0.0D0
*      ZJEU = EPIN0(3)
 
C          DEFORMATION PLASTIQUE DUE A LA COMPRESSION SEULE
C
      EPCO=VAR0(8)
      DEFPNC= EPCO
      ZJEU = EPIN0(3)-EPCO
C
C  AM  CORRECTION POUR JOINT AVEC RESISTANCE EN TRACTION
C                                CAS DU RAPPROCHEMENT DES LEVRES
*- cs deb
        IF (DEPSN.LT.-1.D-4) THEN
*      IF (DEPSN.LT.-1.D-4.AND.SIG0(3).GT.PRESP) THEN
*- cs fin
C                                IL EXISTE UN JEU
         IF (ZJEU.GT.1.D-3) THEN
C                                PASSAGE DIRECT D'UN ETAT OUVERT
C                                A UN ETAT COMPRIME
cc               IF (ABS(DEPSN).GT.ZJEU) THEN
            valdif=(ABS(DEPSN) - ZJEU)/ABS(DEPSN)
ccc               IF (valdif.gt.3.05D-6) THEN
            IF (valdif.gt.1.D-3) THEN
C
               ALP = ZJEU/(ABS(DEPSN))
               STOTA1= DEPS1
               DEPS1 = ALP * DEPS1
               STOTA2= DEPS2
               DEPS2 = ALP * DEPS2
               STOCK = DEPSN
               DEPSN = ALP * DEPSN
               INDIC = 1
               GO TO 99
            ENDIF
C
         ENDIF
C
      ENDIF
C
 98   CONTINUE
C
C                                IL RESTE A PARCOURIR LE TRAJET
C                                (1-ALP) * DEPSN
C
      IF (INDIC.EQ.1) THEN
         DEPS1=(1.0D0-ALP)*STOTA1
         DEPS2=(1.0D0-ALP)*STOTA2
         DEPSN=(1.0D0-ALP)*STOCK
         INDIC=0
         IFIN =1
      ENDIF
C
 99   CONTINUE
C
C
C ----CAS GENERAL SANS TRANSITION DIRECTE D'UN ETAT OUVERT A UN ETAT
C ----A UN ETAT COMPRIME
C
C
C
C
      IF ((INDIC.EQ.1).OR.(IFIN.EQ.0)) THEN
C
C--------------------CALCUL DE L ETAT DE CONTRAINTES, DE DEFORMATIONS
C--------------------DES VARIABLES INTERNES AU PAS PRECEDENT
C
C------CHARGEMENT DES VALEURS AU PAS PRECEDENT
C
C
C----------DEFORMATIONS
C
C             DEFORMATIONS INELASTIQUES
C
         EPSS1 = EPIN0(1)
         EPSS2 = EPIN0(2)
         EPSN1 = EPIN0(3)
C
C             DEFORMATIONS ELASTIQUES
C
         DEFS10 = DD(1,1)*SIG0(1)
         DEFS20 = DD(2,2)*SIG0(2)
c
         if (sig0(3).le.presp) then
            DEFN0 = DD(3,3)*SIG0(3)
         else if (sig0(3).gt.presp) then
C
C  AM  CORRECTION POUR JOINT AVEC RESISTANCE EN TRACTION
C
C           DEFN0 = 0.0D0
            DEFN0 = DD(3,3)*SIG0(3)
         endif
C
C             DEFORMATIONS TOTALES AU PAS PRECEDENT N
C
         EPS1N  = EPSS1 + DEFS10
         EPS2N  = EPSS2 + DEFS20
         EPSNN  = EPSN1 + DEFN0
C
      ENDIF
C
C
C
C
C ----CAS DU 2ND PASSAGE : D'UN ETAT FERME ON PASSE A UN ETAT COMPRIME
C
C
C
C
      IF (IFIN.EQ.1) THEN
C
C--------------------CALCUL DE L ETAT DE CONTRAINTES, DE DEFORMATIONS
C--------------------DES VARIABLES INTERNES AU PAS PRECEDENT
C
C------CHARGEMENT DES VALEURS AU PAS PRECEDENT
C
C----------DEFORMATIONS
C
C             DEFORMATIONS INELASTIQUES
C
         EPSS1 = DEFP11
         EPSS2 = DEFP21
*- cs deb
*         EPSN1 = 0.0D0
         EPSN1 = EPCOF
*- cs fin
C
C             DEFORMATIONS ELASTIQUES
C
         DEFS10 = DD(1,1)*S1INT
         DEFS20 = DD(2,2)*S2INT
         DEFN0  = 0.0D0
C
C             DEFORMATIONS TOTALES AU PAS PRECEDENT N
C
         EPS1N  = EPSS1 + DEFS10
         EPS2N  = EPSS2 + DEFS20
         EPSNN  = EPSN1 + DEFN0
C
      ENDIF
C
C             ESTIMATION DE LA DEFO. TOTALE AU PAS N+1
C
      EP1NP1 = EPS1N + DEPS1
      EP2NP1 = EPS2N + DEPS2
      EPNNP1 = EPSNN + DEPSN
      eptmco = epnnp1-EPCO
C
C----------EPSILON ELASTIQUE (1 2 N) AU PAS N+1 A L'ITERATION 0
C
      EE1P10 = DEFS10 + DEPS1
      EE2P10 = DEFS20 + DEPS2
      EENP10 = DEFN0  + DEPSN
C
C----------SIGMA (1 2 N) AU PAS N+1 A L'ITERATION 0
C
      S1NP10 = G1 * EE1P10
      S2NP10 = G2 * EE2P10
C
C
C  AM  CORRECTION POUR JOINT AVEC RESISTANCE EN TRACTION
C
C        if  (depsn.lt.preem) then
      IF  (DEPSN.LT.PREEM.AND.SIG0(3).LE.PRESP) THEN
c
c -----------fermeture du joint
c
*         if (epnnp1.lt.preep) then
         if (eptmco.lt.preep) then
c
c --------------fermeture du joint avec compression des levres
c
*            if (eenp10.le.preep) then
c
c ----------------compression : pente EC
c
               SNNP10 = EC * EENP10
*            endif
C
*            if (ABS(epnnp1).le.preep) then
            if (ABS(eptmco).le.preep) then
c
c ----------------levres jointives
c
               SNNP10 = 0.0D0
            endif
         endif
c
*         if ((epnnp1.ge.preep).or.(indic.eq.1)) then
         if ((eptmco.ge.preep).or.(indic.eq.1)) then
c
c --------------rapprochement des levres sans compression
c
            SNNP10 = 0.0D0
c
         endif
c
      ELSE
c
c -----------ouverture du joint OU FERMETURE SI JOINT PAS CASSE
c
         SNNP10 = EC * EENP10
c
      endif
C
C-----------VARIABLES INTERNES
C
C          DEFORMATION PLASTIQUE EQUIVALENTE
C
      EPSEQ=VAR0(1)
C
C          DEFORMATION PLASTIQUE DUE A L'OUVERTURE SEULE
C
      EPOUN=VAR0(2)
C
C          ETAT DU JOINT
C
C            STAT   = 0 LEVRES EN CONTACT, AVEC OU SANS CISAILLEMENT
C            STAT   = 1 JOINT EN OUVERTURE PENTE EC
C            STAT   = 2 JOINT EN FERMETURE PENTE EC
C            STAT   = 3 JOINT EN FERMETURE PLASTICITE
C
C
C
      STAT =VAR0(3)
C
C           MULTIPLICATEUR PLASTIQUE EN TRACTION
C
      XLAM1=VAR0(4)
C
C      MISE A 0 DE LA LIMITE EN TRACTION SI DEJA CASSE OU SI JEU
C
      RTRAC=FTRC
      IF(XLAM1.GT.0.D0) RTRAC=0.D0
C
C
C          DEFORMATION PLASTIQUE DUE A LA COMPRESSION SEULE
C
*      EPCO=VAR0(8)
      EPCOF=EPCO
      VARF(8)=EPCOF
C
C
C           MULTIPLICATEUR PLASTIQUE EN COMPRESSION
C
      XLAM3=VAR0(9)
      XLAM3F=XLAM3
      VARF(9)=XLAM3F
C
C
C
C-----------ETUDE DE L'ECOULEMENT SELON LES DIFFERENTS CRITERES
C
      CRICIS = CALG(S1NP10,S2NP10,SNNP10,PHI,C)
      CRINOR = SNNP10 - RTRAC
      CRICIP = CALN(S1NP10,S2NP10,SNNP10,PHI,C)
      ABSTAU = SQRT( (S1NP10 ** 2) + (S2NP10 ** 2) )
C
C----------AIGUILLAGE SUIVANT LA ZONE 0, 1, 2, 3
C
      IF ((CRICIS.LE.PRECIS).AND.(CRINOR.LE.PRECIS)) THEN
         GO TO 10
      ENDIF
C     
*      IF ((CRICIS.GT.PRECIS).AND.(CRICIP.LE.PRECIS)) THEN
      IF ((CRICIS.GT.PRECIS).AND.(CRINOR.LE.PRECIS)) THEN      
         GO TO 40
      ENDIF
C     
*      IF ((CRICIP.GT.PRECIS).AND.(ABSTAU.GT.C)) THEN
*      IF ((CRICIP.GT.PRECIS).AND.(ABSTAU.GT.PRECIS)) THEN
      IF ((ABSTAU.GT.PRECIS).AND.(CRINOR.GT.PRECIN)) THEN
         GO TO 20
      ENDIF
C     
*      IF ((ABSTAU.LE.C).AND.(CRINOR.GT.PRECIN)) THEN
      IF ((ABSTAU.LE.PRECIS).AND.(CRINOR.GT.PRECIN)) THEN
         GO TO 30
      ENDIF
C     
C-----------CAS DE LA CONVERGENCE A L ITERATION 0 :
C
 10   CONTINUE
C
C----------CHARGEMENT DES VALEURS DE CONTRAINTES, DEFORMATIONS,
C----------VARIABLES INTERNES APRES CONVERGENCE INTERNE
C
C
      IF (IMPRI.EQ.1) THEN
         WRITE(6,*) 'ETIQUETTE 10 - LES DEUX CRITERES SONT VERIFIES'
      ENDIF
C
*- cs deb
      SNNP11=SNNP10
*- cs fin
C
      IF (INDIC.EQ.1) THEN
C                              FERMETURE DU JOINT SANS COMPRESSION
*         print *,'10 ok1'
         S1INT = S1NP10
         DEFP11= 0.0D0
         S2INT = S2NP10
         DEFP21= 0.0D0
         DEFPN = -(epin0(3)-DEFPNC)
         GO TO 98
C
      ELSE
C
         IF (IFIN.EQ.0) THEN
C
C
c
            defp(1)=0.0d0
            defp(2)=0.0d0
C
            if ( abs(snnp10).le.presp ) then
*            print *,'10 ok2'
c
c            deformation plastique du pas (n+1) moins
c            deformation plastique du pas (n)
c
               defp(3) = epnnp1 - epsn1
c
            else
*            print *,'10 ok3'
c
c            comportement elastoplastique du joint en fermeture
c
*- cs deb
*             DEFP(3) = 0.0D0
             CALL PLAQP(EC,QT,FNE,SNNP10,EPCO,XLAM3,SNNP11,EPCOF,XLAM3F)
             DEFP(3) = EPCOF-EPCO
*- cs fin
c
            endif
c
         ELSE
*         print *,'10 ok4'
c
c             cas du passage direct d'un etat ouvert a un etat
c             comprime
c
            DEFP12  = 0.0D0
            DEFP(1) = 0.0d0
            DEFP22  = 0.0D0
            DEFP(2) = 0.0d0
cc                DEFP(3) = 0.0d0
*- cs deb
*            DEFP(3) = defpn
             CALL PLAQP(EC,QT,FNE,SNNP10,EPCO,XLAM3,SNNP11,EPCOF,XLAM3F)
             DEFP(3) = DEFPN+(EPCOF-EPCO)
 
*- cs fin
            IFIN = 0
         ENDIF
C
C     CALCUL DES VARIABLES INTERNES : DEPEQ  STAT  EPOUN
C
c
c      calcul du delta epsilon plastique entre
c      le pas actuel et le pas precedent
c
         delep1 = defp(1)
         delep2 = defp(2)
         delepn = defp(3)
c
c      produit contracte pour le calcul de la defo. plasti. cumulee
c
         pro1=(delep1 * delep1) + (delep2 * delep2) + (delepn * delepn)
*- cs         depeq = sqrt( (2.0D0/3.0D0) * pro1 )
         depeq = SQRT(pro1)
         VARF(1) = EPSEQ + DEPEQ
         VARF(4) = XLAM1
         VARF(9) = XLAM3F
c
c      calcul de la deformation normale due a l'ouverture seule
c
         if ( abs(snnp11).le.presp ) then
*- cs deb
            epoun = defp(3)-(EPCOF-EPCO)
*- cs fin
            VARF(2) = EPOUN
         else
*            epoun = 0.0d0
            epoun = DEFPN
            VARF(2) = EPOUN
         endif
C
*- cs deb
C      deformation normale plastique due a la compression
C
         VARF(8)=EPCOF
C
*- cs fin
c
c     etat du joint  stat = 0   cisaillement pur, levres en contact
c                    stat = 1   ouverture du joint contrainte nulle
c                    stat = 2   compression du joint elasticite pente EC
c                    stat = 3   compression du joint plasticite
c
         stat = 2.D0
         eptmco = epnnp1-EPCOF
         if ( abs(snnp11).le.presp) then
            stat = 1.D0
            if ( abs(eptmco).le.1.D-8 ) stat = 0.D0
         endif
         F3=CRIPQ(EC,QT,FNE,SNNP11,XLAM3F)
         if (ABS(F3).LE.PRECIS) stat = 3.D0
         VARF(3) = STAT
         iori = 1
C
         SIGF(1) = S1NP10
         SIGF(2) = S2NP10
         SIGF(3) = SNNP11
C
         GOTO 9900
c
      ENDIF
C
C----------CRITERES DE L'EFFORT NORMAL ET DE CISAILLEMENT
C----------TOUS LES DEUX VIOLES EN MEME TEMPS
C
 20   CONTINUE
C
C            ECOULEMENT SELON LES DEUX CRITERES A L'ITERATION ITER
C
      ITER = 1
      DELAM1=0.D0
C
 2    CONTINUE
C
      IF (ITER.GT.MAXITE) THEN
         WRITE(6,*) ' LES DEUX CRITERES SONT VIOLES '
         WRITE(6,*) ' DIVERGENCE APRES 1000 ITERATIONS INTERNES'
         STOP
      ENDIF
C
      IF (IMPRI.EQ.1) THEN
         write(6,*) 'ETIQUETTE 20 - LES DEUX CRITERES SONT VIOLES'
      ENDIF
C
C
C            CALCUL DES CORRECTEURS PLASTIQUES
C
C
C                DERIVEE DU CRITERE F SUIVANT SIGMA
C
      DFIDS1 = 0.0D0
      DFIDS2 = 0.0D0
      DFIDSN = 1.0D0
C
C                DERIVEE DU CRITERE G SUIVANT SIGMA
C
      DENO = SQRT ( (S1NP10**2) + (S2NP10**2) )
      DGIDS1 = S1NP10 / DENO
      DGIDS2 = S2NP10 / DENO
      DGIDSN = TAN(PHI)
      DZGIDS = TAN(ZMU)
C
C                MATRICE DE HOOK * DFI/DSIGMA
C
      HODF1 = G1 * DFIDS1
      HODF2 = G2 * DFIDS2
      HODFN = EC * DFIDSN
C
C                MATRICE DE HOOK * DGI/DSIGMA
C
      HODG1 = G1 * DGIDS1
      HODG2 = G2 * DGIDS2
      HODGN = EC * DZGIDS
C
C                FIFI = DFI/DSIGMA * MATRICE DE HOOK * DFI/DSIGMA
C
      FIFI = (DFIDS1 * HODF1) +  (DFIDS2 * HODF2) + (DFIDSN * HODFN)
C
C                FIGI = DFI/DSIGMA * MATRICE DE HOOK * DGI/DSIGMA
C
      FIGI = (DFIDS1 * HODG1) +  (DFIDS2 * HODG2) + (DFIDSN * HODGN)
C
C                GIFI = DGI/DSIGMA * MATRICE DE HOOK * DFI/DSIGMA
C
      GIFI = (DGIDS1 * HODF1) +  (DGIDS2 * HODF2) + (DGIDSN * HODFN)
C
C                GIFI = DGI/DSIGMA * MATRICE DE HOOK * DFI/DSIGMA
C
      GIGI = (DGIDS1 * HODG1) +  (DGIDS2 * HODG2) + (DGIDSN * HODGN)
C
C                RESOLUTION DU SYSTEME POUR TROUVER LES CORRECTEURS
C                PLASTIQUES DLAM1 ET DLAM2
C
      DET = (FIFI * GIGI) - (FIGI * GIFI)
C
      IF (ABS(DET).LE.PRECIS) THEN
         WRITE(6,*) ' SYSTEME INDETERMINE (F G): ITERATION ',ITER
         STOP
      ENDIF
C
      DLAM1 = (CRINOR * GIGI) - (CRICIS * FIGI)
      DLAM1 = DLAM1 / DET
      DELAM1 = DELAM1 + DLAM1
      RTRAC=0.D0
C
      DLAM2 = (CRICIS * FIFI) - (CRINOR * GIFI)
      DLAM2 = DLAM2 / DET
C
C                CALCUL DU NOUVEL ETAT DE CONTRAINTES
C
      S1NP11=S1NP10 - (G1 * ( (DLAM1*DFIDS1) + (DLAM2*DGIDS1)))
      S2NP11=S2NP10 - (G2 * ( (DLAM1*DFIDS2) + (DLAM2*DGIDS2)))
      SNNP11=SNNP10 - (EC * ( (DLAM1*DFIDSN) + (DLAM2*DZGIDS)))
C
C                ALLER EN 52 POUR LE TEST DE LA CONVERGENCE
      GOTO 52
C
C
 52   CONTINUE
C
C                CONVERGENCE CHECK
C
      CRICIS = CALG(S1NP11,S2NP11,SNNP11,PHI,C)
      CRINOR = SNNP11 - RTRAC
C
      IF ((ABS(CRICIS).LE.PRECIS).AND.(ABS(CRINOR).LE.PRECIS)) THEN
C
C----------LES DEUX CRITERES SONT VERIFIES
C
C
         IF (IMPRI.EQ.1) THEN
            write(6,*) 'PASSAGE N. 5 - 2'
            write(6,*) 'LES DEUX CRITERES SONT VERIFIES'
            write(6,*) 'CONVERGENCE A L ITERATION',ITER
         ENDIF
C
C----------CHARGEMENT DES VALEURS DE CONTRAINTES, DEFORMATIONS,
C----------VARIABLES INTERNES APRES CONVERGENCE INTERNE
C
C
         SIGF(1) = S1NP11
         SIGF(2) = S2NP11
         SIGF(3) = SNNP11
C
         DEFP(1) = ( EP1NP1 - ( DD(1,1) * SIGF(1) ) ) - EPSS1
         DEFP(2) = ( EP2NP1 - ( DD(2,2) * SIGF(2) ) ) - EPSS2
c
         if ( abs(sigf(3)) .le. presp ) then
c
c            deformation plastique du pas (n+1) moins
c            deformation plastique du pas (n)
c
            DEFP(3) = EPNNP1 - EPSN1
         else
c
c        le joint se deforme toujours elastiquement en fermeture
c
            DEFP(3) = 0.0D0
c
         endif
C
C     CALCUL DES VARIABLES INTERNES : DEPEQ  STAT  EPOUN
C
c
c      calcul du delta epsilon plastique entre
c      le pas actuel et le pas precedent
c
         delep1 = defp(1)
         delep2 = defp(2)
         delepn = defp(3)
c
c      produit contracte pour le calcul de la defo. plasti. cumulee
c
         pro1=(delep1 * delep1) + (delep2 * delep2) + (delepn * delepn)
*- cs         depeq = sqrt( (2.0D0/3.0D0) * pro1 )
         depeq = SQRT(pro1)
         VARF(1) = EPSEQ + DEPEQ
         VARF(4) = XLAM1 + DELAM1
c
c      calcul de la deformation normale due a l'ouverture seule
c
         if ( abs(snnp11).le.presp ) then
            epoun = delepn
            VARF(2) = EPOUN
         else
            epoun = 0.0d0
            VARF(2) = EPOUN
         endif
c
c     etat du joint  stat = 0   cisaillement pur, levres en contact
c                    stat = 1   ouverture du joint contrainte nulle
c                    stat = 2   compression du joint elasticite pente EC
c                    stat = 3   compression du joint plasticite
c
         stat = 2.D0
         eptmco = epnnp1-EPCOF
         if ( abs(snnp11).le.presp) then
            stat = 1.D0
            if ( abs(eptmco).le.1.D-8 ) stat = 0.D0
         endif
         F3=CRIPQ(EC,QT,FNE,SNNP11,XLAM3F)
         if (ABS(F3).LE.PRECIS) stat = 3.D0
         VARF(3) = STAT
         iori = 2
c
         GOTO 9900
      ELSE
C
         S1NP10=S1NP11
         S2NP10=S2NP11
         SNNP10=SNNP11
         CRICIS = CALG(S1NP10,S2NP10,SNNP10,PHI,C)
         CRINOR = SNNP10 - RTRAC
         ITER = ITER + 1
C
         GO TO 2
      ENDIF
C
C
C____________________________________________________________________
C
C
C----------CRITERE DE L'EFFORT NORMAL VIOLE MAIS
C----------CRITERE DE CISAILLEMENT VERIFIE
C
 30   CONTINUE
C
C            ECOULEMENT SELON LE CRITERE F A L'ITERATION ITER
C
      ITER = 1
      DELAM1=0.D0
C
 3    CONTINUE
C
      IF (ITER.GT.MAXITE) THEN
         WRITE(6,*) ' LE CRITERE NORMAL EST VIOLE '
         WRITE(6,*) ' DIVERGENCE APRES 1000 ITERATIONS INTERNES'
         STOP
      ENDIF
C
      IF (IMPRI.EQ.1) THEN
         write(6,*) 'ETIQUETTE 30 - CRITERE DE L EFFORT NORMAL VIOLE'
         write(6,*) 'ITERATION',ITER
      ENDIF
C
C
C
C            CALCUL DU CORRECTEUR PLASTIQUE
C
C
C                DERIVEE DU CRITERE F SUIVANT SIGMA
C
      DFIDS1 = 0.0D0
      DFIDS2 = 0.0D0
      DFIDSN = 1.0D0
C
C                MATRICE DE HOOK * DGI/DSIGMA
C
      HODF1 = G1 * DFIDS1
      HODF2 = G2 * DFIDS2
      HODFN = EC * DFIDSN
C
C                FIFI = DFI/DSIGMA * MATRICE DE HOOK * DFI/DSIGMA
C
      FIFI = (DFIDS1 * HODF1) +  (DFIDS2 * HODF2) + (DFIDSN * HODFN)
C
C                RESOLUTION DE L'EQUATION POUR TROUVER LE CORRECTEUR
C                PLASTIQUE DLAM1
C
      IF (ABS(FIFI).LE.PRECIS) THEN
         WRITE(6,*) ' SYSTEME INDETERMINE (F) : ITERATION ',ITER
         STOP
      ENDIF
C
      DLAM1 = CRINOR / FIFI
      DELAM1 = DELAM1 + DLAM1
      RTRAC = 0.D0
C
C                CALCUL DU NOUVEL ETAT DE CONTRAINTES
C
      S1NP11=S1NP10 - (G1 *  (DLAM1*DFIDS1) )
      S2NP11=S2NP10 - (G2 *  (DLAM1*DFIDS2) )
      SNNP11=SNNP10 - (EC *  (DLAM1*DFIDSN) )
C
C                ALLER EN 53 POUR LE TEST DE LA CONVERGENCE
      GOTO 53
C
C
C
 53   CONTINUE
C
C                CONVERGENCE CHECK
C
      CRICIS = CALG(S1NP11,S2NP11,SNNP11,PHI,C)
      CRINOR = SNNP11 - RTRAC
C
      IF (ABS(CRINOR).LE.PRECIS) THEN
C
C----------LES DEUX CRITERES SONT VERIFIES
C
C
*         IF (IMPRI.EQ.1) THEN
*            write(6,*) 'PASSAGE N. 5 - 3'
*            write(6,*) 'LES DEUX CRITERES SONT VERIFIES'
*            write(6,*) 'CONVERGENCE A L ITERATION',ITER
*         ENDIF
C
C----------CHARGEMENT DES VALEURS DE CONTRAINTES, DEFORMATIONS,
C----------VARIABLES INTERNES APRES CONVERGENCE INTERNE
C
C
         SIGF(1) = S1NP11
         SIGF(2) = S2NP11
         SIGF(3) = SNNP11
C
         DEFP(1) = ( EP1NP1 - ( DD(1,1) * SIGF(1) ) ) - EPSS1
         DEFP(2) = ( EP2NP1 - ( DD(2,2) * SIGF(2) ) ) - EPSS2
c
         if ( abs(sigf(3)) .le. presp ) then
c
c            deformation plastique du pas (n+1) moins
c            deformation plastique du pas (n)
c
            DEFP(3) = EPNNP1 - EPSN1
         else
c
c        le joint se deforme toujours elastiquement en fermeture
c
            DEFP(3) = 0.0D0
c
         endif
C
C     CALCUL DES VARIABLES INTERNES : DEPEQ  STAT  EPOUN
C
c
c      calcul du delta epsilon plastique entre
c      le pas actuel et le pas precedent
c
         delep1 = defp(1)
         delep2 = defp(2)
         delepn = defp(3)
c
c      produit contracte pour le calcul de la defo. plasti. cumulee
c
         pro1=(delep1 * delep1) + (delep2 * delep2) + (delepn * delepn)
*- cs         depeq = sqrt( (2.0D0/3.0D0) * pro1 )
         depeq = SQRT(pro1)
         VARF(1) = EPSEQ + DEPEQ
         VARF(4) = XLAM1 + DELAM1
c
c      calcul de la deformation normale due a l'ouverture seule
c
         if ( abs(snnp11).le.presp ) then
            epoun = delepn
            VARF(2) = EPOUN
         else
            epoun = 0.0d0
            VARF(2) = EPOUN
         endif
c
c     etat du joint  stat = 0   cisaillement pur, levres en contact
c                    stat = 1   ouverture du joint contrainte nulle
c                    stat = 2   compression du joint elasticite pente EC
c                    stat = 3   compression du joint plasticite
c
         stat = 2.D0
         eptmco = epnnp1-EPCOF
         if ( abs(snnp11).le.presp) then
            stat = 1.D0
            if ( abs(eptmco).le.1.D-8 ) stat = 0.D0
         endif
         F3=CRIPQ(EC,QT,FNE,SNNP11,XLAM3F)
         if (ABS(F3).LE.PRECIS) stat = 3.D0
         VARF(3) = STAT
         iori = 3
c
         GOTO 9900
      ELSE
C
         S1NP10=S1NP11
         S2NP10=S2NP11
         SNNP10=SNNP11
         CRICIS = CALG(S1NP10,S2NP10,SNNP10,PHI,C)
         CRINOR = SNNP10 - RTRAC
         ITER = ITER + 1
C
         GO TO 3
      ENDIF
C
C
C____________________________________________________________________
C
C
C
C----------CRITERES DE CISAILLEMENT VIOLE MAIS
C----------CRITERE DE L'EFFORT NORMAL VERIFIE
C
 40   CONTINUE
C
C            ECOULEMENT SELON LE CRITERE G A L'ITERATION ITER
C
      ITER = 1
      DELAM3=0.D0
C
 4    CONTINUE
C
      IF (ITER.GT.MAXITE) THEN
         WRITE(6,*) ' LE CRITERE DE CISAILLEMENT EST VIOLE '
         WRITE(6,*) ' DIVERGENCE APRES 1000 ITERATIONS INTERNES'
         STOP
      ENDIF
C
      IF (IMPRI.EQ.1) THEN
         write(6,*) 'ETIQUETTE 40 - CRITERE DE CISAILLEMENT VIOLE'
      ENDIF
C
C
C       LE CRITERE DE PLASTICITE EN COMPRESSION EST-IL VIOLE ?
C
      F3=CRIPQ(EC,QT,FNE,SNNP10,XLAM3)
C
      IF (F3.LE.PRECIS) THEN
C-----------------------------------------------------------------------
C     CALCUL DU CORRECTEUR PLASTIQUE POUR LE CRITERE DE CISAILLEMENT
C-----------------------------------------------------------------------
C
C           DERIVEE DU CRITERE DE COULOMB G SUIVANT SIGMA
C
      DENO = SQRT ( (S1NP10**2) + (S2NP10**2) )
      DGIDS1 = S1NP10 / DENO
      DGIDS2 = S2NP10 / DENO
      DZGIDS = TAN(ZMU)
      DGIDSN = TAN(PHI)
C
C                MATRICE DE HOOK * DGI/DSIGMA
C
      HODG1 = G1 * DGIDS1
      HODG2 = G2 * DGIDS2
      HODGN = EC * DZGIDS
C
C                GIGI = DGI/DSIGMA * MATRICE DE HOOK * DGI/DSIGMA
C
      GIGI = (DGIDS1 * HODG1) +  (DGIDS2 * HODG2) + (DGIDSN * HODGN)
C
C
C                RESOLUTION DE L'EQUATION POUR TROUVER LE CORRECTEUR
C                PLASTIQUE DLAM2
C
      IF (ABS(GIGI).LE.PRECIS) THEN
         WRITE(6,*) ' SYSTEME INDETERMINE (G) : ITERATION ',ITER
         STOP
      ENDIF
C
C     Cisaillement
      DLAM2 = CRICIS / GIGI
C
C
C     Plasticité en compression
      DLAM3 = 0.
      DELAM3 = DELAM3 + DLAM3
      XLAM3F=XLAM3+DELAM3
C
C                CALCUL DU NOUVEL ETAT DE CONTRAINTES
C
CC       SNNP11=SNNP10 - (EC *  (DLAM2*DGIDSN) )
C
      S1NP11=S1NP10 - (G1 *  (DLAM2*DGIDS1) )
      S2NP11=S2NP10 - (G2 *  (DLAM2*DGIDS2) )
      SNNP11=SNNP10 - (EC *  (DLAM2*DZGIDS) )
C
      ELSE
C-----------------------------------------------------------------------
C     CALCUL DES CORRECTEURS PLAST. POUR LE CRITERE DE CISAILLEMENT ET
C                  DE PLASTICITE EN COMPRESSION
C-----------------------------------------------------------------------
 
C           DERIVEE DU CRITERE DE PLASTICITE F3 SUIVANT SIGMA ET LE
C                  MULTIPLICATEUR PLASTIQUE XLAM3
C
C
      DFIDS1 = 0.0D0
      DFIDS2 = 0.0D0
      DFIDSN = SNNP10/ABS(SNNP10)
      DFIDLA = -HH
C
C           DERIVEE DU CRITERE DE COULOMB G SUIVANT SIGMA
C
      DENO = SQRT ( (S1NP10**2) + (S2NP10**2) )
      DGIDS1 = S1NP10 / DENO
      DGIDS2 = S2NP10 / DENO
      DZGIDS = TAN(ZMU)
      DGIDSN = TAN(PHI)
C
C                MATRICE DE HOOK * DFI/DSIGMA
C
      HODF1 = G1 * DFIDS1
      HODF2 = G2 * DFIDS2
      HODFN = EC * DFIDSN
C
C                MATRICE DE HOOK * DGI/DSIGMA
C
      HODG1 = G1 * DGIDS1
      HODG2 = G2 * DGIDS2
      HODGN = EC * DZGIDS
C
C
C     FIFI = DFI/DSIGMA * MATRICE DE HOOK * DFI/DSIGMA - DFI/DLAM2
C
      FIFI = (DFIDS1*HODF1)+(DFIDS2*HODF2)+(DFIDSN*HODFN)-DFIDLA
C
C                FIGI = DFI/DSIGMA * MATRICE DE HOOK * DGI/DSIGMA
C
      FIGI = (DFIDS1 * HODG1) +  (DFIDS2 * HODG2) + (DFIDSN * HODGN)
C
C                GIFI = DGI/DSIGMA * MATRICE DE HOOK * DFI/DSIGMA
C
      GIFI = (DGIDS1 * HODF1) +  (DGIDS2 * HODF2) + (DGIDSN * HODFN)
C
C                GIGI = DGI/DSIGMA * MATRICE DE HOOK * DGI/DSIGMA
C
      GIGI = (DGIDS1 * HODG1) +  (DGIDS2 * HODG2) + (DGIDSN * HODGN)
C
C                RESOLUTION DU SYSTEME POUR TROUVER LES CORRECTEURS
C                PLASTIQUES DLAM2 ET DLAM3
C
      DET = (FIFI * GIGI) - (FIGI * GIFI)
C
      IF (ABS(DET).LE.PRECIS) THEN
         WRITE(6,*) ' SYSTEME INDETERMINE (F G): ITERATION ',ITER
         STOP
      ENDIF
C
C     Cisaillement
      DLAM2 = (CRICIS * FIFI) - (F3 * GIFI)
      DLAM2 = DLAM2 / DET
C
C     Plasticité en compression
      DLAM3 = (F3 * GIGI) - (CRICIS * FIGI)
      DLAM3 = DLAM3 / DET
      DELAM3 = DELAM3 + DLAM3
      XLAM3F=XLAM3+DELAM3
      EPCOF=EPCOF+DLAM3*DFIDSN
C
C
C                CALCUL DU NOUVEL ETAT DE CONTRAINTES
C
      S1NP11=S1NP10 - (G1 * ( (DLAM2*DGIDS1) + (DLAM3*DFIDS1)))
      S2NP11=S2NP10 - (G2 * ( (DLAM2*DGIDS2) + (DLAM3*DFIDS2)))
      SNNP11=SNNP10 - (EC * ( (DLAM2*DZGIDS) + (DLAM3*DFIDSN)))
C
      ENDIF
C-----------------------------------------------------------------------
C      FIN DE CALCUL DU OU DES CORRECTEURS PLASTIQUES
C-----------------------------------------------------------------------
C
      if  (SNNP11.ge.0.0d0)  SNNP11 = 0.0d0
C
C
 54   CONTINUE
C
C                CONVERGENCE CHECK
C
      CRICIS = CALG(S1NP11,S2NP11,SNNP11,PHI,C)
      F3=CRIPQ(EC,QT,FNE,SNNP11,XLAM3F)
C
      IF ((ABS(CRICIS).LE.PRECIS.AND.(F3.LE.PRECIS))) THEN
C
         IF (INDIC.EQ.1) THEN
C
C----------1ER TRAJET : ON PARCOURT ALP * DEPSN
C
            DEFPN = -(epin0(3)-DEFPNC)
            S1INT = S1NP11
            S2INT = S2NP11
            DEFP11 = EP1NP1 - ( DD(1,1) * S1NP11 )
            DEFP21 = EP2NP1 - ( DD(2,2) * S2NP11 )
            GO TO 98
         ELSE
            IF (IFIN.EQ.0) THEN
C
C
C----------LES DEUX CRITERES SONT VERIFIES
C
C
               IF (IMPRI.EQ.1) THEN
                  write(6,*) 'CONVERGENCE A L ITERATION',ITER
               ENDIF
C
C----------CHARGEMENT DES VALEURS DE CONTRAINTES, DEFORMATIONS,
C----------VARIABLES INTERNES APRES CONVERGENCE INTERNE
C
C
               SIGF(1) = S1NP11
               SIGF(2) = S2NP11
               SIGF(3) = SNNP11
C
               DEFP(1) = ( EP1NP1 - ( DD(1,1) * SIGF(1) ) ) - EPSS1
               DEFP(2) = ( EP2NP1 - ( DD(2,2) * SIGF(2) ) ) - EPSS2
c
               if ( abs(sigf(3)) .le. presp ) then
c
c            deformation plastique du pas (n+1) moins
c            deformation plastique du pas (n)
c
                  DEFP(3) = EPNNP1 - EPSN1
               else
c
c        le joint se deforme plastiquement en fermeture
c
C                  DEFP(3) = 0.0D0
                  DEFP(3) = EPCOF-EPCO
c
               endif
            ELSE
C
C----------CAS DU PASSAGE DIRECT D'UN ETAT OUVERT A UN ETAT COMPRIME
C
               DEFP12  = EP1NP1 - (DD(1,1)*S1NP11)
               DEFP22  = EP2NP1 - (DD(2,2)*S2NP11)
               DEFP(1) =  DEFP12 - epin0(1)
               DEFP(2) =  DEFP22 - epin0(2)
*               DEFP(3) = -epin0(3)
               DEFP(3) = DEFPN+(EPCOF-EPCO)
               SIGF(1) = S1NP11
               SIGF(2) = S2NP11
               SIGF(3) = SNNP11
               IFIN    = 0
            ENDIF
C
         ENDIF
C
C     CALCUL DES VARIABLES INTERNES : DEPEQ  STAT  EPOUN
C
c
c      calcul du delta epsilon plastique entre
c      le pas actuel et le pas precedent
c
         delep1 = defp(1)
         delep2 = defp(2)
         delepn = defp(3)
c
c      produit contracte pour le calcul de la defo. plasti. cumulee
c
         pro1=(delep1 * delep1) + (delep2 * delep2) + (delepn * delepn)
*- cs         depeq = sqrt( (2.0D0/3.0D0) * pro1 )
         depeq = SQRT(pro1)
 
         VARF(1) = EPSEQ + DEPEQ
         VARF(4) = XLAM1
         VARF(8) = EPCOF
         VARF(9) = XLAM3F
c
c      calcul de la deformation normale due a l'ouverture seule
c
         if ( abs(snnp11).le.presp ) then
*            epoun = delepn
            epoun = delepn-(EPCOF-EPCO)
            VARF(2) = EPOUN
         else
*            epoun = 0.0d0
            epoun = DEFPN
            VARF(2) = EPOUN
         endif
c
c     etat du joint
c                stat = 0  levres en contact, avec ou sans cisaillement
c                stat = 1  ouverture du joint pente EC
c                stat = 2  compression du joint elastique pente EC
c                stat = 3  compression du joint plastique
c
         stat = 2.D0
*         eptmco = epnnp1-EPCOF
         if ( abs(SIGF(3)).le.presp) then
            stat = 0.D0
            if (SIG0(3).GT.XZERO) stat = 1.D0
            if ( abs(eptmco).gt.1.D-8 ) stat = 1.D0
         endif         
         F3=CRIPQ(EC,QT,FNE,SIGF(3),XLAM3F)
         if (DABS(F3).LE.PRECIS) stat = 3.D0
         VARF(3) = STAT
         iori = 4
c
         GOTO 9900
      ELSE
C
         S1NP10=S1NP11
         S2NP10=S2NP11
         SNNP10=SNNP11
*         CRICIS = CALG(S1NP10,S2NP10,SNNP10,PHI,C)
         ITER = ITER + 1
C
         GO TO 4
      ENDIF
C
C
C____________________________________________________________________
C
 9900 CONTINUE
C ----------CONTRAINTES A LA FIN DU PAS POUR LES 3 DIRECTIONS
      SIGF1(IDIR1)=SIGF(1)
      SIGF1(IDIR2)=SIGF(2)
      SIGF1(IDIRN)=SIGF(3)
C
      I1=4
      I2=6
      IF (IDIRN.GE.4) THEN
        I1=1
        I2=3
      ENDIF
      DO 350 I=I1,I2
        SIGF1(I)=XZERO
        IF (VARF(3).EQ.1.D0.AND.IORI.NE.4) GOTO 350
        DO 300 J=I1,I2
          SIGF1(I)=SIGF1(I)+TRAID(I,J)*(EPST00(J)+DEPST0(J))
 300    CONTINUE
 350  CONTINUE
C
C ----------DEF. PLASTIQUES A LA FIN DU PAS POUR LES 3 DIRECTIONS
      DO I=1,3
        VARF(IDEBLA+I)=EPIN0(I)+DEFP(I)
      END DO
C
      DEFP1(IDIR1)=VARF(IDEBLA+1)
      DEFP1(IDIR2)=VARF(IDEBLA+2)
      DEFP1(IDIRN)=VARF(IDEBLA+3)
C
      RETURN
      END
 
 
 
 
 
 
 
 
 
 
 
 

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