Télécharger la plaquette Cast3M


fr en

Télécharger cbetoc.eso

Retour à la liste

Numérotation des lignes : 

cbetoc

  1. C CBETOC    SOURCE    PV        17/12/08    21:15:11     9660           
  2. C BETOCY    SOURCE    AM        95/12/27    21:15:06     1962
  3.       SUBROUTINE CBETOC(wrk52,wrk53,wrk54,WRK2,NCURVT,NCURVC)
  4.       IMPLICIT INTEGER(I-N)
  5.       IMPLICIT REAL*8(A-H,O-Z)
  6. C
  7. C MODELE DE MACONNERIE ENDOMMAGEABLE EN CYCLIQUE
  8. C
  9. C=======================================================================
  10. C CETTE ROUTINE EST APPELE DANS ECOUL2
  11. C
  12. C ENTREES:
  13. C -------
  14. C SIG0      = CONTRAINTES AU DEBUT DU PAS D'INTEGRATION
  15. C DSIGT     = INCREMENT DE CONT. CALCULE ELASTIQUEMENT A PARTIR DE
  16. C             L'INCREMENT DES DEFORMATIONS TOTALES
  17. C VAR0      = VARIABLES INTERNES AU DEBUT DU PAS D'INTEGRATION
  18. C
  19. C XMAT      = CARACTERISTIQUES MECANIQUES DU MATERIAU
  20. C
  21. C SORTIES:
  22. C -------
  23. C SIGF      = CONTR. A LA FIN DU PAS D'INTEGRATION
  24. C VARF      = VARIABLES INTERNES A LA FIN DU PAS D'INTEGRATION
  25. C DEFP      = INCREMENT DES DEFORM. PLASTIQUES A LA FIN DU PAS
  26. C             D'INTEGRATION
  27. C KERRE      = INDICE QUI REGIT LES ERREURS
  28. C
  29. C=========================================================================
  30. C      VARIABLES NECESSAIRES
  31. C=========================================================================
  32.  
  33. -INC PPARAM
  34. -INC CCOPTIO
  35. -INC DECHE
  36. C
  37. C  Segment pour le Broyden
  38. C
  39.       SEGMENT QUASIN
  40.         REAL*8 XH0(MN,MN)
  41.         REAL*8 WUP(MN,ITM),YUP(MN,ITM)
  42.         INTEGER IT
  43.         REAL*8 R(MN)
  44.         REAL*8 D(MN),DES(MN)
  45.       ENDSEGMENT
  46. C  Segment des courbes d'ecrouissage
  47.        SEGMENT WRK2
  48.          REAL*8 TRAC(LTRAC)
  49.        ENDSEGMENT
  50.  
  51.        DIMENSION  SIGT(4),SIGI(4), DSIGI(4), DSIGIDI(4)
  52. C  Matrice de compliance (C-1)
  53.        DIMENSION CINV(3)
  54. C  Variables d'ecrouissage isotrope
  55.        DIMENSION XKISO(2)
  56. C  Travaux plastiques cycliques et Deformations plastiques equivalentes
  57. C  associes au mecanisme de compression  de la grande surf.
  58.        REAL*8 WORK
  59.        DIMENSION EPSI(2),EPSI0(2),EPSII(2)
  60. C  Contraintes principales et variables d'ecrouissage cinematique
  61.        DIMENSION SIGTP(2), SIGHP(2),SIGHP0(2)
  62. C  Cos et Sin  de l'angle definissant la premiere direction principale
  63.        REAL*8 CPHI, CPHI0, CPHI00, SPHI, SPHI0, SPHI00
  64. C  Indicateur du contact
  65.        DIMENSION KONTAC(4),KONTA0(4)
  66. C  Vecteurs contrainte transitoire
  67.        DIMENSION SITPRO(4)
  68. C  Jacobien initial
  69.        DIMENSION XJ00(3,3)
  70. C  Variable phi apparaissant au denominateur (calcul de la normale et de A)
  71.        REAL*8  PHI
  72. C  Matrices contenant A-1
  73. C  AINV --> Inverse de A pour la grande surface
  74. C                      [ A B 0] A=AINV(1)
  75. C  AINVGRANDESURFACE = [-B A 0] B=AINV(2)
  76. C                      [ 0 0 C] C=AINV(3)
  77.        DIMENSION AINV(3)
  78. C  Vecteurs normaux aux surface
  79. C         XNOR1 --> 1eme mecanisme
  80. C         XNOR2 --> 2eme mecanisme
  81.        DIMENSION XNOR1(4), XNOR2(4)
  82. C
  83. C  Multiplicateurs plastiques 2 mecanismes + 2 surfaces dans GAMBDA
  84. C  Valeurs des criteres dans FFF
  85. C  Valeurs des limites dans RESIST
  86. C Ordre pour les 3 tableaux:
  87. C            1 --> Traction petite surface direction 1
  88. C            2 --> Traction petite surface direction 2
  89. C            3 --> Compression petite surface direction 1
  90. C            4 --> Compression petite surface direction 2
  91. C            5 --> Traction grande surface
  92. C            6 --> Compression grande surface
  93.        DIMENSION GAMBDA(6), FFF(6), RESIST(6)
  94. C Valeur necessaire pour le calcul des residus de la grande surface
  95.        DIMENSION XKISO0(2),XKISOI(2)
  96. C
  97. C  Separation de l'increment en deux ou trois parties dans le cas du contact
  98. C
  99.        REAL*8 SPLITT
  100.        REAL*8 XFFFS
  101. C
  102. C  Pour differencier les differents cas
  103. C     ICAS --> NUMERO DU MULTIPLICATEUR PLASTIQUE  A METTRE A JOUR
  104. C     IHARD --> NUMERO DE LA VARIABLE D'ECROUISSAGE CINEMATIQUE
  105. C     ISGNTC --> Signe indiquant si on est en traction ou compression
  106. C     ISGN2  --> Signe indiquant si on prend la plus grande ou la plus
  107. C                petite contrainte principale
  108. C     NCORNER --> Revient-on sur la pointe ?
  109.       INTEGER ICAS, IHARD, ISGNTC, ISGN2
  110.       INTEGER ICASB, IHARDB, ISGNTCB, ISGN2B
  111.       INTEGER NCORNER
  112. C
  113. C  Variable intermediaire valant la resistance de la grande surface
  114. C
  115.       DIMENSION FG(2)
  116. C
  117. C  Critere de convergence
  118. C
  119.       DIMENSION EPSIL(6)
  120. C========================================================================
  121. C Parametres utiles
  122. C========================================================================
  123.       PARAMETER (XZER=0.D0,UN=1.D0,DEUX=2.D0,UNDEMI=.5D0)
  124.       PARAMETER (UNRADE=0.707106781D0)
  125.       PARAMETER (EPSILO=1.D-8,EPSILO2=-1.D-10)
  126. C========================================================================
  127. C
  128. C     CALCUL DE DSIGT
  129. C
  130. C     calcul de la matrice elastique
  131. C
  132.       CMATE = 'ISOTROPE'
  133.       KCAS=1
  134.       CALL DOHMAS(XMAT,CMATE,IFOUR,4,KCAS,DDHOOK,IRTD)
  135.       DO IE1=1,4
  136.         DSIGT(IE1)=XZER
  137.         DO IE2=1,4
  138.           DSIGT(IE1)=DSIGT(IE1)+DDHOOK(IE1,IE2)*DEPST(IE2)
  139.         ENDDO
  140.       ENDDO
  141.       DSIGT(3)=XZER
  142. C=========================================================================
  143. C      LECTURE DES PARAMETRES DU MATERIAU
  144. C
  145. C=========================================================================
  146. C Module d'Young YOUN, Poisson XNU, Shear GG
  147. C Ecrouissage cinematique XHH1 ,
  148.       YOUN = XMAT(1)
  149.       XNU  = XMAT(2)
  150.       GGG  = YOUN/(2*(UN+XNU))
  151.       XHH1 = XMAT(5)
  152.       CPI  = YOUN/(UN-XNU)
  153. C
  154.       RESIST(1) = XMAT(6)
  155.       RESIST(2) = XMAT(6)
  156.       RESIST(3) = XMAT(7)
  157.       RESIST(4) = XMAT(7)
  158.       RESIST(5) = XMAT(8)
  159.       RESIST(6) = XMAT(9)
  160. C Deformation plastique equivalente et travail plastique de ref.
  161. C pour la dissipation cyclique
  162.       EPS0      = XMAT(10)
  163.       WOR0      = XMAT(11)
  164.       DEGRAD    = EPS0/WOR0
  165.       XLCAT      = XMAT(14)
  166.       XLCAC      = XMAT(15)
  167. C
  168. C
  169. C  Evolution des variables d'ecrouissage isotrope
  170. C
  171. C  Voila les positions des points d'entree des deux courbes dans TRAC
  172.       IPT    = 1
  173.       IPC    = IPT+2*NCURVT
  174. C     NCURVT : nombre de point de la courbe de traction
  175. C     NCURVC : nombre de point de la courbe de compression
  176. C===========================================================
  177. C   Fin lecture des caracteristiques
  178. C   On passe aux contraintes et variables d'etat
  179. C===========================================================
  180. C
  181. C  CONTRAINTES ELASTIQUES
  182. C
  183.       SIGI(1) = SIG0(1)
  184.       SIGI(2) = SIG0(2)
  185.       SIGI(4) = SIG0(4)
  186. C DELTA SIGMA TRIAL
  187.       DSIGI(1) = DSIGT(1)
  188.       DSIGI(2) = DSIGT(2)
  189.       DSIGI(4) = DSIGT(4)
  190. C
  191. C  ECROUISSAGE CINEMATIQUE
  192. C
  193.       SIGHP0(1) = VAR0(1)
  194.       SIGHP0(2) = VAR0(2)
  195. C
  196. C   ECROUISSAGE ISOTROPE
  197. C
  198.       EPSI0(1) = VAR0(3)
  199.       EPSI0(2) = VAR0(4)
  200. C
  201.       CPHI00 = VAR0(5) + 1.D0
  202.       CPHI0 = CPHI00
  203.       SPHI00 = VAR0(6)
  204.       SPHI0 = SPHI00
  205. C
  206.       KONTA0(1) = nint(VAR0( 7))
  207.       KONTA0(2) = nint(VAR0( 8))
  208.       KONTA0(3) = nint(VAR0( 9))
  209.       KONTA0(4) = nint(VAR0(10))
  210. C
  211. C  Variables de contact
  212. C
  213.       SPLITT = 0.D0
  214.       NCON = 0
  215.       IDICHO = 0
  216.       NDICHO = 0
  217.       NPASS = 0
  218.       DSIGIDI(1) = XZER
  219.       DSIGIDI(2) = XZER
  220.       DSIGIDI(4) = XZER
  221. C===============================================
  222. C  Fin de la lecture des contraintes et des variables d'etat
  223. C================================================
  224. C
  225. C  Premier calcul des residus puis DRIVER
  226. C================================================
  227. C
  228. C   Entree quand on revient des routines, dernier test avant de finir
  229. C
  230.  10   CONTINUE
  231.       SIGT(1) = SIGI(1) + DSIGI(1)
  232.       SIGT(2) = SIGI(2) + DSIGI(2)
  233.       SIGT(4) = SIGI(4) + DSIGI(4)
  234.       CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  235.       EPSI(1) = EPSI0(1)
  236.       EPSI(2) = EPSI0(2)
  237.       CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  238.      &XKISO(1),DUMM,KERRE)
  239.       CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  240.      &XKISO(2),DUMM,KERRE)
  241.       SIGHP(1) = SIGHP0(1)
  242.       SIGHP(2) = SIGHP0(2)
  243.       NPASS2 = 0
  244.       KONTAC(1) = KONTA0(1)
  245.       KONTAC(2) = KONTA0(2)
  246.       KONTAC(3) = KONTA0(3)
  247.       KONTAC(4) = KONTA0(4)
  248. C
  249. C Calcul des residus initiaux
  250. C
  251.       NMEC   = 0
  252.       FFF(1) = XZER
  253.       FFF(2) = XZER
  254.       FFF(3) = XZER
  255.       FFF(4) = XZER
  256.       FFF(5) = XZER
  257.       FFF(6) = XZER
  258.       FINT1  = XZER
  259.       FINT2  = XZER
  260.       GAMBDA(1) = XZER
  261.       GAMBDA(2) = XZER
  262.       GAMBDA(3) = XZER
  263.       GAMBDA(4) = XZER
  264.       GAMBDA(5) = XZER
  265.       GAMBDA(6) = XZER
  266.       ICAS = 0
  267.       ICASB = 0
  268. C
  269. C Critere de convergence
  270. C
  271.       EPSIL(1) = ABS(RESIST(1)*EPSILO)
  272.       EPSIL(2) = ABS(RESIST(2)*EPSILO)
  273.       EPSIL(3) = ABS(RESIST(3)*EPSILO)
  274.       EPSIL(4) = ABS(RESIST(4)*EPSILO)
  275.       EPSIL(5) = ABS(RESIST(5)*EPSILO)
  276.       EPSIL(6) = ABS(RESIST(6)*EPSILO)
  277. C
  278. C
  279.       IF (KONTAC(1).EQ.0) THEN
  280.         CALL CRITPE (SIGTP(1), SIGHP0(1), RESIST(1), 1, FFF(1))
  281.         FINT1 = 0.
  282.         IF (FFF(1).GE.EPSIL(1)) THEN
  283.          NMEC = NMEC + 1
  284.          IF (NMEC.EQ.1) THEN
  285.            ICAS = 1
  286.            ISGNTC = 1
  287.            IHARD = 1
  288.          ELSE
  289.            ICASB = 1
  290.            ISGNTCB = 1
  291.            IHARDB = 1
  292.          ENDIF
  293.         ENDIF
  294.       ELSE
  295.         FG(1) = XKISO(1)*RESIST(5)
  296.         CALL CRITPE (SIGTP(1), XZER, FG(1), 1, FINT1)
  297.       ENDIF
  298. C
  299.       IF (KONTAC(2).EQ.0) THEN
  300.         CALL CRITPE (SIGTP(2), SIGHP0(2), RESIST(2), 1, FFF(2))
  301.         FINT2 = 0.
  302.         IF (FFF(2).GE.EPSIL(2)) THEN
  303.          NMEC = NMEC + 1
  304.          IF (NMEC.EQ.1) THEN
  305.            ICAS = 2
  306.            ISGNTC = 1
  307.            IHARD = 2
  308.          ELSE
  309.            ICASB = 2
  310.            ISGNTCB = 1
  311.            IHARDB = 2
  312.          ENDIF
  313.         ENDIF
  314.       ELSE
  315.         FG(1) = XKISO(1)*RESIST(5)
  316.         CALL CRITPE (SIGTP(2), XZER, FG(1), 1, FINT2)
  317.       ENDIF
  318. C
  319.       IF ((FINT1.GE.EPSIL(5)).OR.(FINT2.GE.EPSIL(5))) THEN
  320.        NMEC = NMEC + 1
  321.        IF (FINT1.GE.FINT2) THEN
  322.          FFF(5) = FINT1
  323.          IF (NMEC.EQ.1) THEN
  324.            IHARD = 1
  325.            ICAS = 5
  326.            ISGNTC = 1
  327.          ELSE
  328.            IHARDB = 1
  329.            ICASB = 5
  330.            ISGNTCB = 1
  331.          ENDIF
  332.        ELSE
  333.          FFF(5) = FINT2
  334.          IF (NMEC.EQ.1) THEN
  335.            IHARD = 2
  336.            ICAS = 5
  337.            ISGNTC = 1
  338.          ELSE
  339.            IHARDB = 2
  340.            ICASB = 5
  341.            ISGNTCB = 1
  342.          ENDIF
  343.        ENDIF
  344.       ENDIF
  345. C
  346.       IF (KONTAC(3).EQ.0) THEN
  347.         CALL CRITPE (SIGTP(1), SIGHP0(1), RESIST(3), -1, FFF(3))
  348.         FINT1 = 0.
  349.         IF (FFF(3).GE.EPSIL(3)) THEN
  350.          NMEC = NMEC + 1
  351.          IF (NMEC.EQ.1) THEN
  352.            ICAS = 3
  353.            ISGNTC = -1
  354.            IHARD = 1
  355.          ELSE
  356.            ICASB = 3
  357.            ISGNTCB = -1
  358.            IHARDB = 1
  359.          ENDIF
  360.         ENDIF
  361.       ELSE
  362.         FG(1) = XKISO(2)*RESIST(6)
  363.         CALL CRITPE (SIGTP(1), XZER, FG(1), -1, FINT1)
  364.       ENDIF
  365. C
  366.       IF (KONTAC(4).EQ.0) THEN
  367.         CALL CRITPE (SIGTP(2), SIGHP0(2), RESIST(4), -1, FFF(4))
  368.         FINT2 = 0.
  369.         IF (FFF(4).GE.EPSIL(4)) THEN
  370.          NMEC = NMEC + 1
  371.          IF (NMEC.EQ.1) THEN
  372.            ICAS = 4
  373.            ISGNTC = -1
  374.            IHARD = 2
  375.          ELSE
  376.            ICASB = 4
  377.            ISGNTCB = -1
  378.            IHARDB = 2
  379.          ENDIF
  380.         ENDIF
  381.       ELSE
  382.         FG(1) = XKISO(2)*RESIST(6)
  383.         CALL CRITPE (SIGTP(2), XZER, FG(1), -1, FINT2)
  384.       ENDIF
  385. C
  386.       IF ((FINT1.GE.EPSIL(6)).OR.(FINT2.GE.EPSIL(6))) THEN
  387.         NMEC = NMEC + 1
  388.         IF (FINT1.GE.FINT2) THEN
  389.           FFF(6) = FINT1
  390.           IF (NMEC.EQ.1) THEN
  391.             IHARD = 1
  392.             ICAS = 6
  393.             ISGNTC = -1
  394.           ELSE
  395.             IHARDB = 1
  396.             ICASB = 6
  397.             ISGNTCB = -1
  398.           ENDIF
  399.         ELSE
  400.           FFF(6) = FINT2
  401.           IF (NMEC.EQ.1) THEN
  402.             IHARD = 2
  403.             ICAS = 6
  404.             ISGNTC = -1
  405.           ELSE
  406.             IHARDB = 2
  407.             ICASB = 6
  408.             ISGNTCB = -1
  409.           ENDIF
  410.         ENDIF
  411.       ENDIF
  412. C
  413. C
  414. C Le tableau FFF contient les residus
  415. C
  416. C================================================
  417. C               DRIVER
  418. C================================================
  419.       IF (NMEC.EQ.0) GOTO 9999
  420. C
  421.       IF (NMEC.EQ.1) GOTO 1000
  422. C
  423.       IF (NMEC.EQ.2) GOTO 4000
  424. C
  425.       GOTO 9999
  426. C==================================================
  427. C               CAS 1:
  428. C      Un seul mecanisme  est active
  429. C==================================================
  430.  1000 CONTINUE
  431. C
  432.         CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  433. C
  434. C Mis a jour de isgn2
  435. C
  436.         IF (IHARD.EQ.1) THEN
  437.           IF (SIGTP(1).GE.SIGTP(2)) THEN
  438.             ISGN2 = 1
  439.           ELSE
  440.             ISGN2 = -1
  441.           ENDIF
  442.           HARDC = 2
  443.         ELSE
  444.           IF (IHARD.EQ.2) THEN
  445.             IF (SIGTP(2).GE.SIGTP(1)) THEN
  446.               ISGN2 = 1
  447.             ELSE
  448.               ISGN2 = -1
  449.             ENDIF
  450.           ENDIF
  451.           HARDC = 1
  452.         ENDIF
  453. C
  454.         IF (IHARD.EQ.1) THEN
  455.           IHARDC = 2
  456.         ELSE
  457.           IHARDC = 1
  458.         ENDIF
  459. C
  460.         NCORNER = 0
  461. C
  462.  1010 CONTINUE
  463. C
  464. C  Definition du segment
  465. C
  466.        ITM = 30
  467.        MN = 1
  468.        SEGINI,QUASIN
  469.        IT = -1
  470. C
  471. C  Calcul de xhh2 et xkiso
  472. C
  473.        IF (ICAS.GE.5) THEN
  474.          IF (ICAS.EQ.5) THEN
  475.           CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  476.      &XKISO(1),XHH2,KERRE)
  477.           CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  478.      &XKISO(2),DUMM,KERRE)
  479.          ELSE
  480.           CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  481.      &XKISO(1),DUMM,KERRE)
  482.           CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  483.      &XKISO(2),XHH2,KERRE)
  484.          ENDIF
  485.          EPSI(1) = EPSI0(1)
  486.          XKISO0(1) = XKISO(1)
  487.          EPSI(2) = EPSI0(2)
  488.          XKISO0(2) = XKISO(2)
  489.          XHH2I = XHH2
  490.          EPSII(1) = EPSI0(1)
  491.          XKISOI(1) = XKISO(1)
  492.          EPSII(2) = EPSI0(2)
  493.          XKISOI(2) = XKISO(2)
  494.        ENDIF
  495. C
  496.        DO 1011,IE1=1,6
  497.          FFF(IE1) = XZER
  498.  1011  CONTINUE
  499. C
  500.         IF (ICAS.LE.4) THEN
  501.           CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD),RESIST(ICAS),
  502.      &               ISGNTC,FFF(ICAS))
  503.         ELSE
  504.           FG(1) = XKISOI(ICAS-4)*RESIST(ICAS)
  505.           CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FINT1)
  506.           CALL CRITPE (SIGTP(IHARDC),XZER,FG(1),ISGNTC,FINT2)
  507.           FFF(ICAS) = MAX(FINT1,FINT2)
  508.         ENDIF
  509. C
  510. C Calcul de la normale et du jacobien
  511. C
  512.        IF (ICAS.LE.4) THEN
  513.          CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  514.          XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)-XHH1
  515.        ELSE
  516.          CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,2,XNOR1)
  517.          XJ00(1,1) = -XXCCYY(XNOR1,XNOR1,YOUN,XNU)
  518.      &            -ISGNTC*RESIST(ICAS)*XHH2
  519.        ENDIF
  520. C
  521.        R(1) = + FFF(ICAS)
  522.        GAMBDA(1)=XZER
  523.        GAMBDA(2)=XZER
  524.        GAMBDA(3)=XZER
  525.        GAMBDA(4)=XZER
  526.        GAMBDA(5)=XZER
  527.        GAMBDA(6)=XZER
  528. C
  529.        XJ0INV = UN/XJ00(1,1)
  530. C
  531.        CALL ZERO(D,1,MN)
  532.        XH0(1,1) = XJ0INV
  533. C
  534. C  On rentre dans les iterations internes
  535. C
  536.       DO I=0,ITM
  537. C
  538. C  Appel de Broyden
  539. C
  540.         CALL BROYDE (QUASIN)
  541. C
  542.         GAMBDA(ICAS) = D(1)
  543. C
  544. C  Calcul du phi et des termes de la matrice A et A-1
  545. C
  546.         SITPRO(1) = SIGI(1) + DSIGI(1)
  547.      &       - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  548.      &       + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  549.         SITPRO(2) = SIGI(2) + DSIGI(2)
  550.      &       - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  551.      &       + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  552.         SITPRO(4) = SIGI(4) + DSIGI(4)
  553. C
  554. C Mis a jour des variables d'ecrouissage
  555. C
  556.         SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  557.         SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  558.         EPSI(1) = EPSI0(1) + GAMBDA(5)
  559.         IF (ICAS.LE.4) THEN
  560.           WORK = GAMBDA(ICAS)*RESIST(ICAS)*ISGNTC
  561.         ELSE
  562.           WORK = XZER
  563.         ENDIF
  564.         EPSI(2) = EPSI0(2) + GAMBDA(6) + DEGRAD*WORK
  565. C
  566.         IF (ICAS.LE.4) THEN
  567.            CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  568.      &XKISO(1),DUMM,KERRE)
  569.            CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  570.      &XKISO(2),DUMM,KERRE)
  571.         ELSE
  572.           IF (ICAS.EQ.5) THEN
  573.             CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  574.      &XKISO(1),XHH2,KERRE)
  575.             CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  576.      &XKISO(2),DUMM,KERRE)
  577.           ELSE
  578.             CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  579.      &XKISO(1),DUMM,KERRE)
  580.             CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  581.      &XKISO(2),XHH2,KERRE)
  582.           ENDIF
  583.         ENDIF
  584. C
  585. C Calcul des phi
  586. C
  587.         IF (ICAS.LE.4) THEN
  588.           PHI = ISGNTC*(-UNDEMI*(SITPRO(1) + SITPRO(2)) + SIGHP(IHARD)
  589.      &      +  RESIST(ICAS))
  590.         ELSE
  591.           PHI = ISGNTC*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  592.      &         +(XKISOI(ICAS - 4)+(EPSI(ICAS-4) - EPSII(ICAS-4))*XHH2I)
  593.      &         *RESIST(ICAS))
  594.         ENDIF
  595. C
  596. C  Calcul de AINV
  597. C
  598.         SOMLAM = GAMBDA(1) + GAMBDA(2) + GAMBDA(3) + GAMBDA(4)
  599.      & + GAMBDA(5) + GAMBDA(6)
  600.         AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  601.         AINV(2) = GGG*SOMLAM/(DEUX*PHI + DEUX*GGG*SOMLAM)
  602.         AINV(3) = PHI/(PHI + (GGG*SOMLAM))
  603. C
  604. C   Calcul des contraintes a l'aide de AINV
  605. C
  606.         SIGT(1) = AINV(1)*SITPRO(1) + AINV(2)*SITPRO(2)
  607.         SIGT(2) = AINV(2)*SITPRO(1) + AINV(1)*SITPRO(2)
  608.         SIGT(4) = AINV(3)*SITPRO(4)
  609. C
  610. C   Nouvelles contraintes principales
  611. C
  612.         CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  613. C
  614. C
  615.         IF (ICAS.LE.4) THEN
  616.           CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD),RESIST(ICAS),
  617.      &               ISGNTC,FFF(ICAS))
  618.           R(1)= + FFF(ICAS)
  619.         ELSE
  620.            FG(1) = (XKISOI(ICAS-4)+(EPSI(ICAS-4)-EPSII(ICAS-4))*XHH2I)
  621.      &             *RESIST(ICAS)
  622. C
  623.           CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FINT1)
  624.           CALL CRITPE (SIGTP(IHARDC),XZER,FG(1),ISGNTC,FINT2)
  625.           FFF(ICAS) = MAX (FINT1,FINT2)
  626. C
  627.           R(1) = + FFF(ICAS)
  628.      &         - ISGNTC*RESIST(ICAS)
  629.      &         * (XKISO(ICAS-4) - XKISOI(ICAS-4)
  630.      &            -(EPSI(ICAS-4) - EPSII(ICAS-4))*XHH2I)
  631.           IF (ABS(EPSI(ICAS-4) - EPSII(ICAS-4)).LT.EPSILO2) THEN
  632.             XHH2I=XHH2
  633.           ELSE
  634.             XHH2I = (XKISO(ICAS-4) - XKISOI(ICAS-4))/
  635.      &              (EPSI(ICAS-4) - EPSII(ICAS-4))
  636.           ENDIF
  637.           XKISOI(ICAS-4) = XKISO(ICAS-4)
  638.           EPSII(ICAS-4) = EPSI(ICAS-4)
  639.         ENDIF
  640. C
  641. C  Test de convergence
  642. C
  643. C        IF (ABS(FFF(ICAS)).LT.EPSIL(ICAS)) GOTO 1020
  644.         IF ((ABS(FFF(ICAS)).LT.EPSIL(ICAS))
  645.      &      .AND.(ABS(R(1)).LT.EPSIL(ICAS))) GOTO 1020
  646.       ENDDO
  647. C
  648. C  Pas de convergence
  649. C
  650.       SEGSUP,QUASIN
  651.       GOTO 6000
  652. C
  653.  1020 CONTINUE
  654.       SEGSUP,QUASIN
  655. C
  656.         IF (GAMBDA(ICAS).LT.EPSILO2) THEN
  657.            WRITE (IOIMP,*)'BETOCY: Multiplicateurs plastiques < 0'
  658.            KERRE = 2
  659.            RETURN
  660.          ELSE
  661.            CONTINUE
  662.          ENDIF
  663. C
  664. C  On verifie qu'on ne viole pas les autres criteres
  665. C
  666.         IF (IHARD.EQ.2) THEN
  667.           IHARDC=1
  668.         ELSE
  669.           IHARDC=2
  670.         ENDIF
  671. C
  672.         IF (ICAS.LE.2) THEN
  673.           FG(1) = XKISO(1)*RESIST(5)
  674.           FG(2) = XKISO(2)*RESIST(6)
  675.         ELSE
  676.           IF (ICAS.LE.4) THEN
  677.             FG(1) = XKISO(2)*RESIST(6)
  678.             FG(2) = XKISO(1)*RESIST(5)
  679.           ELSE
  680.             IF (ICAS.EQ.5) THEN
  681.               FG(1) = XKISO(1)*RESIST(5)
  682.               FG(2) = XKISO(2)*RESIST(6)
  683.             ELSE
  684.               FG(1) = XKISO(2)*RESIST(6)
  685.               FG(2) = XKISO(1)*RESIST(5)
  686.             ENDIF
  687.           ENDIF
  688.         ENDIF
  689. C
  690. C  Test GS
  691. C
  692.         IF (ICAS.GE.5) THEN
  693. C Regles de suivi si on est sur la grande surface
  694. C 1- Direction de l'ecoulement
  695.           SIGHP(IHARD) = FG(1) -RESIST(2*ICAS-10+IHARD)
  696. C 2- Direction perpendiculaire a l'ecoulement
  697. C (si on se trouve dans le coin)
  698.             XLIMIT = FG(1) -RESIST(2*ICAS-10+IHARDC)
  699.             IF (ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  700.               SIGHP(IHARDC) = XLIMIT
  701.               KONTAC(2*ICAS-10+IHARDC)=1
  702.             ENDIF
  703.        ENDIF
  704. C
  705. C
  706. C  On teste les autres mecanismes
  707. C
  708. C Test des contacts lateraux
  709. C     Meme signe signtc
  710.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  711.           ICASB = ICAS + 1
  712.         ELSE
  713.           IF ((ICAS.EQ.2).OR.(ICAS.EQ.4)) THEN
  714.             ICASB = ICAS -1
  715.           ELSE
  716.             ICASB = 2*ICAS-10+IHARDC
  717.           ENDIF
  718.         ENDIF
  719.         IF (KONTAC(ICASB).EQ.0) THEN
  720.           CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC), RESIST(ICASB),
  721.      &ISGNTC,XFFFS)
  722.           IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  723.      &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  724.             IHARDB = IHARDC
  725.             ISGNTCB = ISGNTC
  726.             GOTO 1100
  727.           ELSE
  728.             IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  729.      &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  730.               KONTAC(ICASB)=1
  731.             ENDIF
  732.           ENDIF
  733.         ELSE
  734.           CALL CRITPE(SIGTP(IHARDC), XZER, FG(1), ISGNTC,XFFFS)
  735.           IF (ICASB.LE.2) THEN
  736.             ICASB = 5
  737.           ELSE
  738.             ICASB = 6
  739.           ENDIF
  740.           IF (XFFFS.GE.EPSIL(ICASB)) THEN
  741.             IHARDB = IHARDC
  742.             ISGNTCB = ISGNTC
  743.             GOTO 1100
  744.           ENDIF
  745.         ENDIF
  746. C     Autre signe signtc
  747.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  748.           ICASB = ICAS + 1
  749.         ELSE
  750.           IF ((ICAS.EQ.2).OR.(ICAS.EQ.4)) THEN
  751.             ICASB = ICAS - 1
  752.           ELSE
  753.             ICASB = 2*ICAS-10+IHARDC
  754.           ENDIF
  755.         ENDIF
  756.         IF (ICASB.LE.2) THEN
  757.            ICASB = ICASB +2
  758.         ELSE
  759.            ICASB = ICASB -2
  760.         ENDIF
  761.         IF (KONTAC(ICASB).EQ.0) THEN
  762.           CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC), RESIST(ICASB),
  763.      &-ISGNTC,XFFFS)
  764.           IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  765.      &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  766.             IHARDB = IHARDC
  767.             ISGNTCB = -ISGNTC
  768.             GOTO 1100
  769.           ELSE
  770.             IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  771.      &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  772.               KONTAC(ICASB)=1
  773.             ENDIF
  774.           ENDIF
  775.         ELSE
  776.           CALL CRITPE(SIGTP(IHARDC), XZER, FG(2), -ISGNTC,XFFFS)
  777.           IF (ICASB.LE.2) THEN
  778.             ICASB = 5
  779.           ELSE
  780.             ICASB = 6
  781.           ENDIF
  782.           IF (XFFFS.GE.EPSIL(ICASB)) THEN
  783.             IHARDB = IHARDC
  784.             ISGNTCB = -ISGNTC
  785.             GOTO 1100
  786.           ENDIF
  787.         ENDIF
  788. C Test du contact frontal
  789.         IF (ICAS.LE.4) THEN
  790.           CALL CRITPE(SIGTP(IHARD), XZER, FG(1), ISGNTC,XFFFS)
  791.           IF (((XFFFS.GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  792.      &((XFFFS.GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  793.             GOTO 1040
  794.           ELSE
  795.             IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  796.      &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  797.               KONTAC(ICAS) = 1
  798.             ENDIF
  799.           ENDIF
  800.         ENDIF
  801. C
  802. C Il faut faire suivre la petite surface si la grande se retrecit trop vite
  803. C
  804.         IF (ICAS.LE.4) THEN
  805.             IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  806.               ICASB = ICAS +1
  807.             ELSE
  808.               ICASB = ICAS -1
  809.             ENDIF
  810.             XLIMIT = FG(1) -RESIST(ICASB)
  811.             IF (ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  812.               SIGHP(IHARDC) = XLIMIT
  813.               KONTAC(ICASB) = 1
  814.             ENDIF
  815.             IF (ICASB.LE.2) THEN
  816.               ICASB = ICASB +2
  817.             ELSE
  818.               ICASB = ICASB -2
  819.             ENDIF
  820.             XLIMIT = FG(2) -RESIST(ICASB)
  821.             IF (-ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  822.               SIGHP(IHARDC) = XLIMIT
  823.               KONTAC(ICASB) = 1
  824.             ENDIF
  825.         ENDIF
  826. C
  827. C  Mis a jour des indicateurs de contact
  828. C
  829.         IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  830.      &(RESIST(1)-RESIST(3))) THEN
  831.           IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  832.           IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  833.         ELSE
  834.           KONTAC(1) = 1
  835.           KONTAC(2) = 1
  836.           KONTAC(3) = 1
  837.           KONTAC(4) = 1
  838.         ENDIF
  839. C
  840.         ICASB = ICAS
  841.         GOTO 9999
  842. C
  843. C Cas ou on trouve le contact frontal
  844. C
  845.  1040   CONTINUE
  846.           SIGHP(1) = SIGHP0(1)
  847.           SIGHP(2) = SIGHP0(2)
  848.           SIGT(1) = SIGI(1) + DSIGI(1)
  849.           SIGT(2) = SIGI(2) + DSIGI(2)
  850.           SIGT(4) = SIGI(4) + DSIGI(4)
  851.           GOTO 3000
  852. C
  853. C Cas ou on viole un critere lateral
  854. C
  855.  1100   CONTINUE
  856.           SIGHP(1) = SIGHP0(1)
  857.           SIGHP(2) = SIGHP0(2)
  858.           SIGT(1) = SIGI(1) + DSIGI(1)
  859.           SIGT(2) = SIGI(2) + DSIGI(2)
  860.           SIGT(4) = SIGI(4) + DSIGI(4)
  861.           IF (ICAS.EQ.ICASB) THEN
  862.             IF (ICAS.GE.5) THEN
  863.               GOTO 1000
  864.             ELSE
  865.               KERRE = 2
  866.               RETURN
  867.             ENDIF
  868.           ELSE
  869.             GOTO 4000
  870.           ENDIF
  871. C================================================================
  872. C                 CAS 2
  873. C    Un seul mecanisme de chaque surface
  874. C    mais entree en contact des 2 mecanismes
  875. C================================================================
  876.  3000  CONTINUE
  877. C
  878.         CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  879. C
  880. C
  881. C Calcul de la normale a l'estimation elastique et du jacobien initial
  882. C
  883.         IF (IHARD.EQ.1) THEN
  884.           IF (SIGTP(1).GE.SIGTP(2)) THEN
  885.             ISGN2 = 1
  886.           ELSE
  887.             ISGN2 = -1
  888.           ENDIF
  889.         ELSE
  890.           IF (IHARD.EQ.2) THEN
  891.             IF (SIGTP(2).GE.SIGTP(1)) THEN
  892.               ISGN2 = 1
  893.             ELSE
  894.               ISGN2 = -1
  895.             ENDIF
  896.           ENDIF
  897.         ENDIF
  898. C
  899.  3010 CONTINUE
  900. C
  901. C  Definition du segment
  902. C
  903.       ITM = 30
  904.       MN = 2
  905.       SEGINI,QUASIN
  906.       IT = -1
  907.       CALL ZERO(D,1,MN)
  908. C
  909. C  Calcul de xhh2 et xkiso
  910. C
  911.       IF (ICAS.LE.2) THEN
  912.         CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  913.      &XKISO(1),XHH2,KERRE)
  914.         CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  915.      &XKISO(2),DUMM,KERRE)
  916.       ELSE
  917.         CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  918.      &XKISO(1),DUMM,KERRE)
  919.         CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  920.      &XKISO(2),XHH2,KERRE)
  921.       ENDIF
  922.       EPSI(1) = EPSI0(1)
  923.       XKISO0(1) = XKISO(1)
  924.       EPSI(2) = EPSI0(2)
  925.       XKISO0(2) = XKISO(2)
  926.       XHH2I = XHH2
  927.       EPSII(1) = EPSI0(1)
  928.       XKISOI(1) = XKISO(1)
  929.       EPSII(2) = EPSI0(2)
  930.       XKISOI(2) = XKISO(2)
  931. C
  932.        DO 3011,IE1=1,6
  933.          FFF(IE1) = XZER
  934.  3011  CONTINUE
  935. C
  936.         IF (ISGNTC.GT.0) THEN
  937.           FG(1) = XKISO(1)*RESIST(5)
  938.         ELSE
  939.           FG(1) = XKISO(2)*RESIST(6)
  940.         ENDIF
  941.         CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD)
  942.      &              ,RESIST(ICAS),ISGNTC,FFF(ICAS))
  943.         CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC, XFFFS)
  944. C
  945.       CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  946.       XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)-XHH1
  947.       XJ00(2,2) = -XXYY(XNOR1,DSIGI)
  948.       XJ00(1,2) = -XXYY(XNOR1,DSIGI)
  949.       XJ00(2,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)
  950.       IF (ISGNTC.LT.0) THEN
  951.         XJ00(2,1)=XJ00(2,1)+RESIST(6)*XHH2*DEGRAD*ISGNTC*RESIST(ICAS)
  952.       ENDIF
  953.       R(1) = FFF(ICAS)
  954.       R(2) = XFFFS
  955.       GAMBDA(1)=XZER
  956.       GAMBDA(2)=XZER
  957.       GAMBDA(3)=XZER
  958.       GAMBDA(4)=XZER
  959.       GAMBDA(5)=XZER
  960.       GAMBDA(6)=XZER
  961. C
  962. C  Inverse du jacobien
  963. C
  964.       XH0(1,1) = XJ00(1,1)
  965.       XH0(1,2) = XJ00(1,2)
  966.       XH0(2,1) = XJ00(2,1)
  967.       XH0(2,2) = XJ00(2,2)
  968.       CALL INVALM (XH0, MN , MN , IRD, 1.D-12)
  969.       IF (IRD.NE.0) THEN
  970.         SEGSUP,QUASIN
  971.         GOTO 6000
  972.       ENDIF
  973. C
  974. C  On rentre dans les iterations internes
  975. C
  976.       DO I=0,ITM
  977. C
  978. C  Appel de Broyden
  979. C
  980.         CALL BROYDE (QUASIN)
  981. C
  982.         GAMBDA(ICAS) = D(1)
  983.         SPLITT = D(2)
  984. C
  985. C  Calcul du phi et des termes de la matrice A et A-1
  986. C
  987.         SITPRO(1) = SIGI(1) + ((UN-SPLITT)*DSIGI(1))
  988.      &            - UNDEMI*(GAMBDA(1) + GAMBDA(2))*CPI
  989.      &            + UNDEMI*(GAMBDA(3) + GAMBDA(4))*CPI
  990.         SITPRO(2) = SIGI(2) + ((UN-SPLITT)*DSIGI(2))
  991.      &            - UNDEMI*(GAMBDA(1) + GAMBDA(2))*CPI
  992.      &            + UNDEMI*(GAMBDA(3) + GAMBDA(4))*CPI
  993.         SITPRO(4) = SIGI(4) + (UN-SPLITT)*DSIGI(4)
  994. C
  995.         SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  996.         SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  997.         EPSI(1)  = EPSI0(1)
  998.         WORK =  GAMBDA(ICAS)*ISGNTC*RESIST(ICAS)
  999.         EPSI(2)  = EPSI0(2) + DEGRAD*WORK
  1000.         CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  1001.      &XKISO(2),XHH2,KERRE)
  1002. C
  1003.         PHI = ISGNTC*(-UNDEMI*(SITPRO(1)+SITPRO(2))+SIGHP(IHARD)
  1004.      &    + RESIST(ICAS))
  1005. C
  1006. C
  1007.         SOMLAM = GAMBDA(1) + GAMBDA(2) + GAMBDA(3) + GAMBDA(4)
  1008.         AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  1009.         AINV(2) = GGG*SOMLAM/(DEUX*PHI + DEUX*GGG*SOMLAM)
  1010.         AINV(3) = PHI/(PHI + GGG*SOMLAM)
  1011. C
  1012. C   Calcul des contraintes a l'aide de AINV
  1013. C
  1014.         SIGT(1) = AINV(1)*SITPRO(1)+AINV(2)*SITPRO(2)
  1015.         SIGT(2) = AINV(2)*SITPRO(1)+AINV(1)*SITPRO(2)
  1016.         SIGT(4) = AINV(3)*SITPRO(4)
  1017. C
  1018. C   Nouvelles contraintes principales
  1019. C
  1020.         CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1021. C
  1022.         NMEC = 0
  1023. C
  1024.         IF (ISGNTC.GT.0) THEN
  1025.           FG(1) = XKISO(1)*RESIST(5)
  1026.         ELSE
  1027.           FG(1) = (XKISOI(2) + XHH2I*(EPSI(2)-EPSII(2)))*RESIST(6)
  1028.         ENDIF
  1029.         CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD)
  1030.      &              ,RESIST(ICAS),ISGNTC,FFF(ICAS))
  1031.         R(1) = + FFF(ICAS)
  1032.         IF (ABS(FFF(ICAS)).GE.EPSIL(ICAS)) THEN
  1033.           NMEC = NMEC + 1
  1034.         ENDIF
  1035.         CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC, XFFFS)
  1036.         R(2) =  XFFFS
  1037.         IF (ICAS.GE.3) THEN
  1038.            R(2) = R(2) + RESIST(6)*(XKISO(2) - XKISOI(2)
  1039.      &              - (EPSI(2) - EPSII(2))*XHH2I)
  1040.           IF (ABS(EPSI(2) - EPSII(2)).LT.EPSILO2) THEN
  1041.             XHH2I=XHH2
  1042.           ELSE
  1043.             XHH2I = (XKISO(2) - XKISOI(2))/
  1044.      &              (EPSI(2) - EPSII(2))
  1045.           ENDIF
  1046.           EPSII(2) = EPSI(2)
  1047.           XKISOI(2) = XKISO(2)
  1048.         ENDIF
  1049.         IF (((ABS(XFFFS).GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  1050.      &((ABS(XFFFS).GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  1051.             NMEC = NMEC + 1
  1052.         ENDIF
  1053. C
  1054. C  Test de convergence
  1055. C
  1056.         IF (NMEC.EQ.0) GOTO 3020
  1057.       ENDDO
  1058. C
  1059. C  Pas de convergence
  1060. C
  1061.       SEGSUP,QUASIN
  1062.       GOTO 6000
  1063. C
  1064.  3020 CONTINUE
  1065. C
  1066. C  On verifie que la convergence ne s'est pas faite sur un point debile
  1067. C
  1068.       SEGSUP,QUASIN
  1069.       IF (((SPLITT-1.).LE.EPSILO).AND.(SPLITT.GE.(-1.*EPSILO))) THEN
  1070. C
  1071.         IF (IHARD.EQ.2) THEN
  1072.           IHARDC=1
  1073.         ELSE
  1074.           IHARDC=2
  1075.         ENDIF
  1076.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) ICASC = ICAS + 1
  1077.         IF ((ICAS.EQ.2).OR.(ICAS.EQ.4)) ICASC = ICAS - 1
  1078. C
  1079.         IF (ICAS.LE.2) THEN
  1080.           FG(1) = XKISO(1)*RESIST(5)
  1081.           FG(2) = XKISO(2)*RESIST(6)
  1082.         ELSE
  1083.           FG(1) = XKISO(2)*RESIST(6)
  1084.           FG(2) = XKISO(1)*RESIST(5)
  1085.         ENDIF
  1086. C
  1087. C Test des contacts lateraux
  1088. C     Meme signe signtc
  1089.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  1090.           ICASB = ICAS + 1
  1091.         ELSE
  1092.           ICASB = ICAS -1
  1093.         ENDIF
  1094.         IF (KONTAC(ICASB).EQ.0) THEN
  1095.           CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC),
  1096.      & RESIST(ICASB), ISGNTC,XFFFS)
  1097.           IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1098.      &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1099.             IHARDB = IHARDC
  1100.             ISGNTCB = ISGNTC
  1101.             GOTO 3100
  1102.           ELSE
  1103.             IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1104.      &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1105.               KONTAC(ICASB)=1
  1106.             ENDIF
  1107.           ENDIF
  1108.         ELSE
  1109.           CALL CRITPE(SIGTP(IHARDC), XZER, FG(1), ISGNTC,XFFFS)
  1110.           IF (ICASB.LE.2) THEN
  1111.              ICASB = 5
  1112.           ELSE
  1113.              ICASB = 6
  1114.           ENDIF
  1115.           IF (XFFFS.GE.EPSIL(ICASB)) THEN
  1116.             IHARDB = IHARDC
  1117.             ISGNTCB = -ISGNTC
  1118.             GOTO 3100
  1119.           ENDIF
  1120.         ENDIF
  1121. C     Autre signe signtc
  1122.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  1123.           ICASB = ICAS + 1
  1124.         ELSE
  1125.           ICASB = ICAS -1
  1126.         ENDIF
  1127.         IF (ICASB.LE.2) THEN
  1128.            ICASB = ICASB +2
  1129.         ELSE
  1130.            ICASB = ICASB -2
  1131.         ENDIF
  1132.         IF (KONTAC(ICASB).EQ.0) THEN
  1133.           CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC),RESIST(ICASB),
  1134.      & -ISGNTC,XFFFS)
  1135.           IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1136.      &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1137.             IHARDB = IHARDC
  1138.             ISGNTCB = ISGNTC
  1139.             GOTO 3100
  1140.           ELSE
  1141.             IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1142.      &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1143.               KONTAC(ICASB)=1
  1144.             ENDIF
  1145.           ENDIF
  1146.         ELSE
  1147.           CALL CRITPE(SIGTP(IHARDC), XZER, FG(2), -ISGNTC,XFFFS)
  1148.           IF (ICASB.LE.2) THEN
  1149.              ICASB = 5
  1150.           ELSE
  1151.              ICASB = 6
  1152.           ENDIF
  1153.           IF (XFFFS.GE.EPSIL(ICASB)) THEN
  1154.             IHARDB = IHARDC
  1155.             ISGNTCB = -ISGNTC
  1156.             GOTO 3100
  1157.           ENDIF
  1158.         ENDIF
  1159. C
  1160. C Il faut faire suivre la petite surface si la grande se retrecit trop vite
  1161. C
  1162.         IF (ICAS.LE.4) THEN
  1163.             IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  1164.               ICASB = ICAS +1
  1165.             ELSE
  1166.               ICASB = ICAS -1
  1167.             ENDIF
  1168.             XLIMIT = FG(1) -RESIST(ICASB)
  1169.             IF (ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  1170.               SIGHP(IHARDC) = XLIMIT
  1171.               KONTAC(ICASB) = 1
  1172.             ENDIF
  1173.             IF (ICASB.LE.2) THEN
  1174.               ICASB = ICASB +2
  1175.             ELSE
  1176.               ICASB = ICASB -2
  1177.             ENDIF
  1178.             XLIMIT = FG(2) -RESIST(ICASB)
  1179.             IF (-ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  1180.               SIGHP(IHARDC) = XLIMIT
  1181.               KONTAC(ICASB) = 1
  1182.             ENDIF
  1183.         ENDIF
  1184. C
  1185. C
  1186. C  Mis a jour des indicateurs de contact
  1187. C
  1188.         KONTAC(ICAS) = 1
  1189.         IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  1190.      &(RESIST(1)-RESIST(3))) THEN
  1191.           IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  1192.           IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  1193.         ELSE
  1194.           KONTAC(1) = 1
  1195.           KONTAC(2) = 1
  1196.           KONTAC(3) = 1
  1197.           KONTAC(4) = 1
  1198.         ENDIF
  1199.         SIGI(1) = SIGT(1)
  1200.         SIGI(2) = SIGT(2)
  1201.         SIGI(4) = SIGT(4)
  1202.         DSIGI(1) = DSIGI(1)*SPLITT
  1203.         DSIGI(2) = DSIGI(2)*SPLITT
  1204.         DSIGI(4) = DSIGI(4)*SPLITT
  1205.         SIGT(1) = SIGI(1) + DSIGI(1)
  1206.         SIGT(2) = SIGI(2) + DSIGI(2)
  1207.         SIGT(4) = SIGI(4) + DSIGI(4)
  1208.         SIGHP0(1) = SIGHP(1)
  1209.         SIGHP0(2) = SIGHP(2)
  1210.         EPSI0(1) = EPSI(1)
  1211.         EPSI0(2) = EPSI(2)
  1212.         CPHI0 = CPHI
  1213.         SPHI0 = SPHI
  1214.         NPASS2 = 0
  1215.         IF (ICAS.LE.2) THEN
  1216.            ICAS = 5
  1217.         ELSE
  1218.            ICAS = 6
  1219.         ENDIF
  1220.         GOTO 1000
  1221. C
  1222. C
  1223.       ELSE
  1224.         GOTO 6000
  1225.       ENDIF
  1226. C
  1227. C Cas ou on viole un critere lateral
  1228. C
  1229.  3100   CONTINUE
  1230.           SIGHP(1) = SIGHP0(1)
  1231.           SIGHP(2) = SIGHP0(2)
  1232.           SIGT(1) = SIGI(1) + DSIGI(1)
  1233.           SIGT(2) = SIGI(2) + DSIGI(2)
  1234.           SIGT(4) = SIGI(4) + DSIGI(4)
  1235.           CALL CONPRI ( SIGT ,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1236.           CALL VALEUR(EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1),
  1237.      &DUMM,KERRE)
  1238.           CALL VALEUR(EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2),
  1239.      &DUMM,KERRE)
  1240.             CALL CRITPE(SIGTP(IHARD), SIGHP(IHARD), RESIST(ICAS),
  1241.      & ISGNTC, FFF(ICAS))
  1242.           IF (ICASB.LE.4) THEN
  1243.             CALL CRITPE(SIGTP(IHARDB), SIGHP(IHARDB), RESIST(ICASB),
  1244.      & ISGNTCB, FFF(ICASB))
  1245.           ELSE
  1246.             FG(1) = XKISO(ICASB-4)*RESIST(ICASB)
  1247.             CALL CRITPE (SIGTP(IHARDB), XZER, FG(1),
  1248.      & ISGNTCB,FFF(ICASB))
  1249.           ENDIF
  1250.           IF (ICAS.LE.2) THEN
  1251.             FG(1) = XKISO(1)*RESIST(5)
  1252.           ELSE
  1253.             FG(1) = XKISO(2)*RESIST(6)
  1254.           ENDIF
  1255.           CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC, XFFFS)
  1256.           GOTO 5000
  1257. C================================================================
  1258. C                 CAS 4
  1259. C           Deux mecanismes
  1260. C       Petite et Grande surfaces (18/07)
  1261. C================================================================
  1262.  4000  CONTINUE
  1263. C
  1264.         CALL CONPRI ( SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1265. C
  1266. C
  1267. C Calcul de la normale a l'estimation elastique et du jacobien initial
  1268. C
  1269.         IF (IHARD.EQ.1) THEN
  1270.           IF (SIGTP(1).GE.SIGTP(2)) THEN
  1271.             ISGN2 = 1
  1272.           ELSE
  1273.             ISGN2 = -1
  1274.           ENDIF
  1275.         ELSE
  1276.           IF (IHARD.EQ.2) THEN
  1277.             IF (SIGTP(2).GE.SIGTP(1)) THEN
  1278.               ISGN2 = 1
  1279.             ELSE
  1280.               ISGN2 = -1
  1281.             ENDIF
  1282.           ENDIF
  1283.         ENDIF
  1284. C
  1285.         IF (IHARDB.EQ.1) THEN
  1286.           IF (SIGTP(1).GE.SIGTP(2)) THEN
  1287.             ISGN2B = 1
  1288.           ELSE
  1289.             ISGN2B = -1
  1290.           ENDIF
  1291.         ELSE
  1292.           IF (IHARDB.EQ.2) THEN
  1293.             IF (SIGTP(1).GE.SIGTP(2)) THEN
  1294.               ISGN2B = -1
  1295.             ELSE
  1296.               ISGN2B = 1
  1297.             ENDIF
  1298.           ENDIF
  1299.         ENDIF
  1300. C
  1301. C
  1302. C
  1303.  4010 CONTINUE
  1304. C
  1305. C  Definition du segment
  1306. C
  1307.       ITM = 30
  1308.       MN = 2
  1309.       SEGINI,QUASIN
  1310.       IT = -1
  1311.       CALL ZERO(D,1,MN)
  1312. C
  1313. C  XKISO et XHH2
  1314. C
  1315.         IF (ICAS.LE.4) THEN
  1316.           IF (ICASB.LE.4) THEN
  1317.            CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  1318.      &XKISO(1),DUMM,KERRE)
  1319.            CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  1320.      &XKISO(2),DUMM,KERRE)
  1321.           ELSE
  1322.             IF (ICASB.EQ.5) THEN
  1323.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1324.      &         ,XHH2B,KERRE)
  1325.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1326.      &         ,DUMM,KERRE)
  1327.             ELSE
  1328.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1329.      &         ,DUMM,KERRE)
  1330.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1331.      &         ,XHH2B,KERRE)
  1332.             ENDIF
  1333.           ENDIF
  1334.         ELSE
  1335.           IF (ICASB.LE.4) THEN
  1336.             IF (ICAS.EQ.5) THEN
  1337.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1338.      &         ,XHH2,KERRE)
  1339.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1340.      &         ,DUMM,KERRE)
  1341.             ELSE
  1342.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1343.      &         ,DUMM,KERRE)
  1344.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1345.      &         ,XHH2,KERRE)
  1346.             ENDIF
  1347.           ELSE
  1348.             IF (ICASB.EQ.5) THEN
  1349.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1350.      &         ,XHH2B,KERRE)
  1351.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1352.      &         ,XHH2,KERRE)
  1353.             ELSE
  1354.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1355.      &         ,XHH2,KERRE)
  1356.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1357.      &         ,XHH2B,KERRE)
  1358.             ENDIF
  1359.           ENDIF
  1360.         ENDIF
  1361.       EPSI(1) = EPSI0(1)
  1362.       EPSI(2) = EPSI0(2)
  1363.       XHH2I = XHH2
  1364.       XHH2IB = XHH2B
  1365.       XKISO0(1) = XKISO(1)
  1366.       XKISO0(2) = XKISO(2)
  1367.       EPSII(1) = EPSI0(1)
  1368.       EPSII(2) = EPSI0(2)
  1369.       XKISOI(1) = XKISO(1)
  1370.       XKISOI(2) = XKISO(2)
  1371. C
  1372.        DO 4011,IE1=1,6
  1373.          FFF(IE1) = XZER
  1374.  4011  CONTINUE
  1375. C
  1376.         IF (ICAS.LE.4) THEN
  1377.           CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD),RESIST(ICAS),
  1378.      &               ISGNTC,FFF(ICAS))
  1379.         ELSE
  1380.           FG(1) = XKISOI(ICAS-4)*RESIST(ICAS)
  1381.           CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FFF(ICAS))
  1382.         ENDIF
  1383. C
  1384.         IF (ICASB.LE.4) THEN
  1385.           CALL CRITPE (SIGTP(IHARDB),SIGHP0(IHARDB),RESIST(ICASB),
  1386.      &               ISGNTCB,FFF(ICASB))
  1387.         ELSE
  1388.           FG(1) = XKISOI(ICASB-4)*RESIST(ICASB)
  1389.           CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,FFF(ICASB))
  1390.         ENDIF
  1391. C
  1392. C  Jcobien initial
  1393. C
  1394.       CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  1395.       CALL NORMBL(SIGT,RESIST(ICASB),ISGNTCB,ISGN2B,1,XNOR2)
  1396. C
  1397.       XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)
  1398.       XJ00(1,2) = -XXCCYY(XNOR1, XNOR2, YOUN, XNU)
  1399.       XJ00(2,1) = -XXCCYY(XNOR2, XNOR1, YOUN, XNU)
  1400.       XJ00(2,2) = -XXCCYY(XNOR2, XNOR2, YOUN, XNU)
  1401. C
  1402. C  Premier Mecanisme
  1403. C
  1404.       IF (ICAS.LE.4) THEN
  1405.         XJ00(1,1) = XJ00(1,1) - XHH1
  1406.       ELSE
  1407.         XJ00(1,1) = XJ00(1,1)-ISGNTC*RESIST(ICAS)*XHH2
  1408.         IF ((ICASB.LE.4).AND.(ICAS.EQ.6)) THEN
  1409.           XJ00(1,2) = XJ00(1,2)
  1410.      & +RESIST(6)*XHH2*DEGRAD*ISGNTCB*RESIST(ICASB)
  1411.         ENDIF
  1412.       ENDIF
  1413. C
  1414. C  Second Mecanisme
  1415. C
  1416.       IF (ICASB.LE.4) THEN
  1417.           XJ00(2,2) = XJ00(2,2) - XHH1
  1418.       ELSE
  1419.         XJ00(2,2) = XJ00(2,2)-ISGNTCB*RESIST(ICASB)*XHH2B
  1420.         IF ((ICAS.LE.4).AND.(ICASB.EQ.6)) THEN
  1421.            XJ00(2,1) = XJ00(2,1)
  1422.      & +RESIST(6)*XHH2B*DEGRAD*ISGNTC*RESIST(ICAS)
  1423.         ENDIF
  1424.       ENDIF
  1425. C
  1426. C
  1427.       R(1) = + FFF(ICAS)
  1428.       R(2) = + FFF(ICASB)
  1429.       GAMBDA(1) = XZER
  1430.       GAMBDA(2) = XZER
  1431.       GAMBDA(3) = XZER
  1432.       GAMBDA(4) = XZER
  1433.       GAMBDA(5) = XZER
  1434.       GAMBDA(6) = XZER
  1435. C
  1436. C  Inverse du jacobien
  1437. C
  1438.       XH0(1,1) = XJ00(1,1)
  1439.       XH0(1,2) = XJ00(1,2)
  1440.       XH0(2,1) = XJ00(2,1)
  1441.       XH0(2,2) = XJ00(2,2)
  1442.       CALL INVALM (XH0, MN , MN , IRD, 1.D-12)
  1443.       IF (IRD.NE.0) THEN
  1444.           SEGSUP,QUASIN
  1445.           GOTO 6000
  1446.       ENDIF
  1447. C
  1448. C  On rentre dans les iterations internes
  1449. C
  1450.       DO I=0,ITM
  1451. C
  1452. C  Appel de Broyden
  1453. C
  1454.         CALL BROYDE (QUASIN)
  1455. C
  1456.         GAMBDA(ICAS) = D(1)
  1457.         GAMBDA(ICASB) = D(2)
  1458. C
  1459. C  Calcul du phi et des termes de la matrice A et A-1
  1460. C
  1461.         SITPRO(1) = SIGI(1) + DSIGI(1)
  1462.      &       - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  1463.      &       + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  1464.         SITPRO(2) = SIGI(2) + DSIGI(2)
  1465.      &       - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  1466.      &       + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  1467.         SITPRO(4) = SIGI(4) + DSIGI(4)
  1468. C
  1469. C  PHI et AINV
  1470. C
  1471.         SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  1472.         SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  1473.         EPSI(1) = EPSI0(1) + GAMBDA(5)
  1474.         IF (ICAS.LE.4) THEN
  1475.           WORK = GAMBDA(ICAS)*RESIST(ICAS)*ISGNTC
  1476.         ELSE
  1477.           WORK = XZER
  1478.         ENDIF
  1479.         IF (ICASB.LE.4) THEN
  1480.            WORK = WORK + GAMBDA(ICASB)*RESIST(ICASB)*ISGNTCB
  1481.         ENDIF
  1482.         EPSI(2) = EPSI0(2) + GAMBDA(6) + DEGRAD*WORK
  1483. C
  1484.         IF (ICAS.LE.4) THEN
  1485.           IF (ICASB.LE.4) THEN
  1486.            CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  1487.      &XKISO(1),DUMM,KERRE)
  1488.            CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  1489.      &XKISO(2),DUMM,KERRE)
  1490.           ELSE
  1491.             IF (ICASB.EQ.5) THEN
  1492.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1493.      &         ,XHH2B,KERRE)
  1494.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1495.      &         ,DUMM,KERRE)
  1496.             ELSE
  1497.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1498.      &         ,DUMM,KERRE)
  1499.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1500.      &         ,XHH2B,KERRE)
  1501.             ENDIF
  1502.           ENDIF
  1503.         ELSE
  1504.           IF (ICASB.LE.4) THEN
  1505.             IF (ICAS.EQ.5) THEN
  1506.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1507.      &         ,XHH2,KERRE)
  1508.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1509.      &         ,DUMM,KERRE)
  1510.             ELSE
  1511.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1512.      &         ,DUMM,KERRE)
  1513.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1514.      &         ,XHH2,KERRE)
  1515.             ENDIF
  1516.           ELSE
  1517.             IF (ICASB.EQ.5) THEN
  1518.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1519.      &         ,XHH2B,KERRE)
  1520.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1521.      &         ,XHH2,KERRE)
  1522.             ELSE
  1523.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1524.      &         ,XHH2,KERRE)
  1525.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1526.      &         ,XHH2B,KERRE)
  1527.             ENDIF
  1528.           ENDIF
  1529.         ENDIF
  1530. C
  1531. C  Calcul des phi
  1532. C
  1533.         IF (ICAS.LE.4) THEN
  1534.            PHI1=ISGNTC*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  1535.      &      + SIGHP(IHARD) + RESIST(ICAS))
  1536.         ELSE
  1537.            PHI1= ISGNTC*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  1538.      &       +(XKISOI(ICAS-4)+XHH2I*(EPSI(ICAS-4)-EPSII(ICAS-4)))
  1539.      &        *RESIST(ICAS))
  1540.         ENDIF
  1541. C
  1542.         IF (ICASB.LE.4) THEN
  1543.            PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  1544.      &   + SIGHP(IHARDB)+  RESIST(ICASB))
  1545.         ELSE
  1546.            PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  1547.      &     +(XKISOI(ICASB-4)+XHH2IB*(EPSI(ICASB-4)-EPSII(ICASB-4)))
  1548.      &       *RESIST(ICASB))
  1549.         ENDIF
  1550. C
  1551. C   Test sur phi1 et phi2
  1552. C
  1553.         IF ((ABS(PHI2)).GT.(ABS(EPSILO*RESIST(ICASB)))) THEN
  1554.            SOMLAM = GAMBDA(ICAS) + GAMBDA(ICASB)*PHI1/PHI2
  1555.            PHI = PHI1
  1556.         ELSE
  1557.           IF ((ABS(PHI1)).GT.(ABS(EPSILO*RESIST(ICAS)))) THEN
  1558.             SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)*PHI2/PHI1
  1559.             PHI = PHI2
  1560.           ELSE
  1561.             SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)
  1562.             PHI = PHI1
  1563.           ENDIF
  1564.         ENDIF
  1565. C
  1566.         AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  1567.         AINV(2) = (GGG*SOMLAM)/(DEUX*PHI + DEUX*GGG*SOMLAM)
  1568.         AINV(3) = PHI/(PHI + GGG*SOMLAM)
  1569. C
  1570. C   Calcul des contraintes a l'aide de AINV
  1571. C
  1572.         SIGT(1) = AINV(1)*SITPRO(1)+AINV(2)*SITPRO(2)
  1573.         SIGT(2) = AINV(2)*SITPRO(1)+AINV(1)*SITPRO(2)
  1574.         SIGT(4) = AINV(3)*SITPRO(4)
  1575. C
  1576. C   Nouvelles contraintes principales
  1577. C
  1578.         CALL CONPRI ( SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1579. C
  1580. C
  1581.         IF (ICAS.LE.4) THEN
  1582.           CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD),RESIST(ICAS),
  1583.      &               ISGNTC,FFF(ICAS))
  1584.           R(1)= + FFF(ICAS)
  1585.         ELSE
  1586.           FG(1) = (XKISOI(ICAS-4)+XHH2I*(EPSI(ICAS-4)-EPSII(ICAS-4)))
  1587.      &           *RESIST(ICAS)
  1588.           CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FFF(ICAS))
  1589.           R(1) = + FFF(ICAS)
  1590.      &         - ISGNTC*RESIST(ICAS)
  1591.      &         * (XKISO(ICAS-4) - XKISOI(ICAS-4)
  1592.      &            -(EPSI(ICAS-4) - EPSII(ICAS-4))*XHH2I)
  1593.           IF (ABS(EPSI(ICAS-4) - EPSII(ICAS-4)).LT.EPSILO2) THEN
  1594.             XHH2I=XHH2
  1595.           ELSE
  1596.             XHH2I = (XKISO(ICAS-4) - XKISOI(ICAS-4))/
  1597.      &              (EPSI(ICAS-4) - EPSII(ICAS-4))
  1598.           ENDIF
  1599.           EPSII(ICAS-4) = EPSI(ICAS-4)
  1600.           XKISOI(ICAS-4) = XKISO(ICAS-4)
  1601.         ENDIF
  1602. C
  1603.         IF (ICASB.LE.4) THEN
  1604.           CALL CRITPE (SIGTP(IHARDB),SIGHP(IHARDB),RESIST(ICASB),
  1605.      &               ISGNTCB,FFF(ICASB))
  1606.           R(2)= + FFF(ICASB)
  1607.         ELSE
  1608.           FG(1) = (XKISOI(ICASB-4)+XHH2IB*(EPSI(ICASB-4) -
  1609.      &     EPSII(ICASB-4)))*RESIST(ICASB)
  1610.           CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,FFF(ICASB))
  1611.           R(2) = + FFF(ICASB)
  1612.      &         - ISGNTCB*RESIST(ICASB)
  1613.      &         * (XKISO(ICASB-4) - XKISOI(ICASB-4)
  1614.      &            -(EPSI(ICASB-4) - EPSII(ICASB-4))*XHH2IB)
  1615.           IF (ABS(EPSI(ICASB-4) - EPSII(ICASB-4)).LT.EPSILO2) THEN
  1616.             XHH2IB=XHH2B
  1617.           ELSE
  1618.             XHH2IB = (XKISO(ICASB-4) - XKISOI(ICASB-4))/
  1619.      &              (EPSI(ICASB-4) - EPSII(ICASB-4))
  1620.           ENDIF
  1621.           XKISOI(ICASB-4) = XKISO(ICASB-4)
  1622.           EPSII(ICASB-4) = EPSI(ICASB-4)
  1623.         ENDIF
  1624. C
  1625. C  Test de convergence
  1626. C
  1627.         IF ((ABS(FFF(ICAS)).LE.EPSIL(ICAS))
  1628.      &.AND.(ABS(FFF(ICASB)).LE.EPSIL(ICASB)).AND.
  1629.      &(ABS(R(1)).LE.EPSIL(ICAS)).AND.(ABS(R(2)).LE.EPSIL(ICASB))) THEN
  1630.            GOTO 4020
  1631.         ENDIF
  1632.       ENDDO
  1633. C
  1634. C  Pas de convergence
  1635. C
  1636.       SEGSUP,QUASIN
  1637.       GOTO 6000
  1638. C
  1639.  4020 CONTINUE
  1640.       SEGSUP,QUASIN
  1641. C
  1642.          IF ((GAMBDA(ICAS).LT.EPSILO2).OR.
  1643.      &(GAMBDA(ICASB).LT.EPSILO2)) THEN
  1644.  4030 CONTINUE
  1645.              IF (NPASS.LE.5) THEN
  1646.                 NPASS = NPASS + 1
  1647.              ELSE
  1648.                GOTO 6000
  1649.              ENDIF
  1650.            IF ((GAMBDA(ICAS).LT.EPSILO2).AND.
  1651.      &(GAMBDA(ICASB).LT.EPSILO2)) THEN
  1652. C  Ce cas n'arrive que si les normales ont ete mal estimes:
  1653. C   (increment elastique trop grand) --> on bissecte
  1654.               GOTO 6000
  1655.            ELSE
  1656.              IF (GAMBDA(ICAS).LT.EPSILO2) THEN
  1657.                ICAS = ICASB
  1658.                IHARD = IHARDB
  1659.                ISGNTC = ISGNTCB
  1660.              ELSE
  1661.                CONTINUE
  1662.              ENDIF
  1663.            ENDIF
  1664.            SIGHP(1) = SIGHP0(1)
  1665.            SIGHP(2) = SIGHP0(2)
  1666.            SIGT(1) = SIGI(1) + DSIGI(1)
  1667.            SIGT(2) = SIGI(2) + DSIGI(2)
  1668.            SIGT(4) = SIGI(4) + DSIGI(4)
  1669.            GOTO 1000
  1670.          ENDIF
  1671. C  On teste les mecanismes de la grande surface
  1672. C  Si un des criteres est viole (et donc s'il n'y a pas contact)
  1673. C  il faut trouver le point de contact
  1674. C
  1675.       CALL CONPRI ( SIGT ,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1676. C
  1677.       NCON = 0
  1678.       FINT1 = 0.
  1679.       FINT2 = 0.
  1680. C
  1681.       XFFFS = XZER
  1682.       XFFFSB = XZER
  1683.       IF (ICAS.LE.4) THEN
  1684.         IF (ICAS.LE.2) THEN
  1685.           FG(1) = XKISO(1)*RESIST(5)
  1686.         ELSE
  1687.           FG(1) = XKISO(2)*RESIST(6)
  1688.         ENDIF
  1689.         CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,XFFFS)
  1690.         IF (((XFFFS.GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  1691.      &       ((XFFFS.GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  1692.           NCON = NCON +1
  1693.         ELSE
  1694.           IF (((ICAS.LE.2).AND.(ABS(XFFFS).LT.EPSIL(5))).OR.
  1695.      &((ICAS.GE.3).AND.(ABS(XFFFS).LT.EPSIL(6)))) THEN
  1696.             KONTAC(ICAS) = 1
  1697.           ELSE
  1698.             CONTINUE
  1699.           ENDIF
  1700.         ENDIF
  1701.       ENDIF
  1702. C
  1703.       IF (ICASB.LE.4) THEN
  1704.         IF (ICASB.LE.2) THEN
  1705.           FG(1) = XKISO(1)*RESIST(5)
  1706.         ELSE
  1707.           FG(1) = XKISO(2)*RESIST(6)
  1708.         ENDIF
  1709.         CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,XFFFSB)
  1710.         IF (((XFFFSB.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1711.      &       ((XFFFSB.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1712.           NCON = NCON +1
  1713.         ELSE
  1714.           IF (((ICASB.LE.2).AND.(ABS(XFFFSB).LT.EPSIL(5))).OR.
  1715.      &((ICASB.GT.3).AND.(ABS(XFFFSB).LT.EPSIL(6)))) THEN
  1716.             KONTAC(ICASB) = 1
  1717.           ELSE
  1718.             CONTINUE
  1719.           ENDIF
  1720.         ENDIF
  1721.       ENDIF
  1722. C
  1723.       IF (NCON.EQ.1) THEN
  1724.         IF (((XFFFS.GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  1725.      &       ((XFFFS.GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  1726.           GOTO 4040
  1727.         ELSE
  1728.           ICASI = ICAS
  1729.           ISGNTCI = ISGNTC
  1730.           IHARDI = IHARD
  1731.           ICAS = ICASB
  1732.           ISGNTC = ISGNTCB
  1733.           IHARD = IHARDB
  1734.           ICASB = ICASI
  1735.           ISGNTCB = ISGNTCI
  1736.           IHARDB = IHARDI
  1737.           GOTO 4040
  1738.         ENDIF
  1739.       ELSE
  1740.         IF (NCON.EQ.2) THEN
  1741.           GOTO 4040
  1742.         ELSE
  1743.           CONTINUE
  1744.         ENDIF
  1745.         CONTINUE
  1746.       ENDIF
  1747. C
  1748. C Regles de suivi avant de sortir normalement
  1749. C   Direction de l'ecoulement
  1750.       IF (ICAS.GE.5) THEN
  1751.         SIGHP(IHARD) = XKISO(ICAS-4)*RESIST(ICAS)
  1752.      &              -RESIST(2*ICAS-10+IHARD)
  1753.       ENDIF
  1754.       IF (ICASB.GE.5) THEN
  1755.         SIGHP(IHARDB) = XKISO(ICASB-4)*RESIST(ICASB)
  1756.      &              -RESIST(2*ICASB-10+IHARDB)
  1757.       ENDIF
  1758. C
  1759. C  Mis a jour des indicateurs de contact
  1760. C
  1761.       IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  1762.      &(RESIST(1)-RESIST(3))) THEN
  1763.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  1764.         IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  1765.         IF ((ICASB.EQ.1).OR.(ICASB.EQ.2)) KONTAC(ICASB+2)=0
  1766.         IF ((ICASB.EQ.3).OR.(ICASB.EQ.4)) KONTAC(ICASB-2)=0
  1767.       ELSE
  1768.         KONTAC(1) = 1
  1769.         KONTAC(2) = 1
  1770.         KONTAC(3) = 1
  1771.         KONTAC(4) = 1
  1772.       ENDIF
  1773. C
  1774.       GOTO 9999
  1775. C
  1776. C  On va au contact
  1777.  4040 CONTINUE
  1778.       SIGHP(1) = SIGHP0(1)
  1779.       SIGHP(2) = SIGHP0(2)
  1780.       SIGT(1) = SIGI(1) + DSIGI(1)
  1781.       SIGT(2) = SIGI(2) + DSIGI(2)
  1782.       SIGT(4) = SIGI(4) + DSIGI(4)
  1783.       EPSI(1) = EPSI0(1)
  1784.       EPSI(2) = EPSI0(2)
  1785.       GOTO 5000
  1786. C================================================================
  1787. C                 CAS 5
  1788. C    Deux  mecanismes avec
  1789. C    Entree en contact de 2 mecanismes
  1790. C================================================================
  1791.  5000  CONTINUE
  1792. C
  1793.         CALL CONPRI ( SIGT , CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1794. C
  1795. C  ISGN2 et ISGN2B
  1796.       IF (IHARD.EQ.1) THEN
  1797.         IF (SIGTP(1).GE.SIGTP(2)) THEN
  1798.           ISGN2 = 1
  1799.         ELSE
  1800.           ISGN2 = -1
  1801.         ENDIF
  1802.       ELSE
  1803.         IF (IHARD.EQ.2) THEN
  1804.           IF (SIGTP(1).GE.SIGTP(2)) THEN
  1805.             ISGN2 = -1
  1806.           ELSE
  1807.             ISGN2 = 1
  1808.           ENDIF
  1809.         ENDIF
  1810.       ENDIF
  1811. C
  1812.         IF (IHARDB.EQ.1) THEN
  1813.           IF (SIGTP(1).GE.SIGTP(2)) THEN
  1814.             ISGN2B = 1
  1815.           ELSE
  1816.             ISGN2B = -1
  1817.           ENDIF
  1818.         ELSE
  1819.           IF (IHARDB.EQ.2) THEN
  1820.             IF (SIGTP(1).GE.SIGTP(2)) THEN
  1821.               ISGN2B = -1
  1822.             ELSE
  1823.               ISGN2B = 1
  1824.             ENDIF
  1825.           ENDIF
  1826.         ENDIF
  1827. C
  1828.  5010 CONTINUE
  1829. C
  1830. C  Definition du segment
  1831. C
  1832.       ITM = 30
  1833.       MN = 3
  1834.       SEGINI,QUASIN
  1835.       IT = -1
  1836.       CALL ZERO(D,1,MN)
  1837. C  XKISO et XHH2
  1838.         IF (ICASB.LE.4) THEN
  1839.           IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1840.             CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  1841.      &XKISO(1),XHH2,KERRE)
  1842.             CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  1843.      &XKISO(2),DUMM,KERRE)
  1844.           ELSE
  1845.             CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  1846.      &XKISO(1),DUMM,KERRE)
  1847.             CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  1848.      &XKISO(2),XHH2,KERRE)
  1849.           ENDIF
  1850.         ELSE
  1851.           IF (ICASB.EQ.5) THEN
  1852.             IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1853.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1854.      &         ,XHH2,KERRE)
  1855.               XHH2B=XHH2
  1856.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1857.      &         ,DUMM,KERRE)
  1858.             ELSE
  1859.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1860.      &         ,XHH2B,KERRE)
  1861.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1862.      &         ,XHH2,KERRE)
  1863.             ENDIF
  1864.           ELSE
  1865.             IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1866.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1867.      &         ,XHH2,KERRE)
  1868.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1869.      &         ,XHH2B,KERRE)
  1870.             ELSE
  1871.               CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1872.      &         ,DUMM,KERRE)
  1873.               CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1874.      &         ,XHH2,KERRE)
  1875.               XHH2B = XHH2
  1876.             ENDIF
  1877.           ENDIF
  1878.         ENDIF
  1879. C
  1880.         XKISO0(1)=XKISO(1)
  1881.         EPSI(1)=EPSI0(1)
  1882.         XKISO0(2)=XKISO(2)
  1883.         EPSI(2)=EPSI0(2)
  1884.         XHH2I = XHH2
  1885.         XHH2IB = XHH2B
  1886.         XKISOI(1)=XKISO(1)
  1887.         XKISOI(2)=XKISO(2)
  1888.         EPSII(1)=EPSI0(1)
  1889.         EPSII(2)=EPSI0(2)
  1890. C
  1891.        DO 5011,IE1=1,6
  1892.          FFF(IE1) = XZER
  1893.  5011  CONTINUE
  1894. C
  1895.           CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD),RESIST(ICAS),
  1896.      &               ISGNTC,FFF(ICAS))
  1897. C
  1898.         IF (ICASB.LE.4) THEN
  1899.           CALL CRITPE (SIGTP(IHARDB),SIGHP0(IHARDB),RESIST(ICASB),
  1900.      &               ISGNTCB,FFF(ICASB))
  1901.         ELSE
  1902.           FG(1) = XKISOI(ICASB -4)*RESIST(ICASB)
  1903.           CALL CRITPE (SIGTP(IHARDB), XZER, FG(1),ISGNTCB,FFF(ICASB))
  1904.         ENDIF
  1905. C
  1906.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1907.           FG(1) = XKISO(1)*RESIST(5)
  1908.           CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC,XFFFS)
  1909.         ELSE
  1910.           FG(1) = XKISO(2)*RESIST(6)
  1911.           CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,XFFFS)
  1912.         ENDIF
  1913. C
  1914.         CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  1915.         CALL NORMBL(SIGT,RESIST(ICASB),ISGNTCB,ISGN2B,1,XNOR2)
  1916. C
  1917. C   Jacobien
  1918. C
  1919.       XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU) - XHH1
  1920.       XJ00(1,2) = -XXCCYY(XNOR1, XNOR2, YOUN, XNU)
  1921.       XJ00(1,3) = -XXYY(XNOR1 , DSIGI)
  1922.       XJ00(2,1) = -XXCCYY(XNOR2, XNOR1, YOUN, XNU)
  1923.       XJ00(2,2) = -XXCCYY(XNOR2, XNOR2, YOUN, XNU)
  1924.       XJ00(2,3) = -XXYY(XNOR2,DSIGI)
  1925.       XJ00(3,1) = -XXCCYY(XNOR1,XNOR1,YOUN,XNU)
  1926.       XJ00(3,2) = -XXCCYY(XNOR1,XNOR2,YOUN,XNU)
  1927.       XJ00(3,3) = -XXYY(XNOR1,DSIGI)
  1928. C
  1929.       IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) THEN
  1930.         XJ00(3,1)=XJ00(3,1)+RESIST(6)*XHH2*DEGRAD*ISGNTC*RESIST(ICAS)
  1931.       ENDIF
  1932. C
  1933.       IF (ICASB.LE.4) THEN
  1934.         XJ00(2,2) = XJ00(2,2) - XHH1
  1935.         IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) THEN
  1936.           XJ00(3,2)=XJ00(3,2)+RESIST(6)*XHH2
  1937.      &*DEGRAD*ISGNTCB*RESIST(ICASB)
  1938.         ENDIF
  1939.       ELSE
  1940.         XJ00(2,2) = XJ00(2,2) - RESIST(ICASB)*XHH2B*ISGNTCB
  1941.       ENDIF
  1942. C
  1943. C
  1944.       R(1) = FFF(ICAS)
  1945.       R(2) = FFF(ICASB)
  1946.       R(3) = XFFFS
  1947. C
  1948.       SPLITT = XZER
  1949.       GAMBDA(1) = XZER
  1950.       GAMBDA(2) = XZER
  1951.       GAMBDA(3) = XZER
  1952.       GAMBDA(4) = XZER
  1953.       GAMBDA(5) = XZER
  1954.       GAMBDA(6) = XZER
  1955. C
  1956. C  Inverse du jacobien
  1957. C
  1958.       XH0(1,1) = XJ00(1,1)
  1959.       XH0(1,2) = XJ00(1,2)
  1960.       XH0(1,3) = XJ00(1,3)
  1961.       XH0(2,1) = XJ00(2,1)
  1962.       XH0(2,2) = XJ00(2,2)
  1963.       XH0(2,3) = XJ00(2,3)
  1964.       XH0(3,1) = XJ00(3,1)
  1965.       XH0(3,2) = XJ00(3,2)
  1966.       XH0(3,3) = XJ00(3,3)
  1967.       CALL INVALM (XH0, MN , MN , IRD, 1.D-12)
  1968.       IF (IRD.NE.0) THEN
  1969.         SEGSUP,QUASIN
  1970.         GOTO 6000
  1971.       ENDIF
  1972. C
  1973. C  On rentre dans les iterations internes
  1974. C
  1975.       DO I=0,ITM
  1976. C
  1977. C  Appel de Broyden
  1978. C
  1979.         CALL BROYDE (QUASIN)
  1980. C
  1981.         GAMBDA(ICAS) = D(1)
  1982.         GAMBDA(ICASB) = D(2)
  1983.         SPLITT = D(3)
  1984. C
  1985. C  Calcul du phi et des termes de la matrice A et A-1
  1986. C
  1987.         SITPRO(1) = SIGI(1) + (UN-SPLITT)*DSIGI(1)
  1988.      &       - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  1989.      &       + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  1990.         SITPRO(2) = SIGI(2) + (UN-SPLITT)*DSIGI(2)
  1991.      &       - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  1992.      &       + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  1993.         SITPRO(4) = SIGI(4) + (UN-SPLITT)*DSIGI(4)
  1994. C
  1995. C  PHI et AINV
  1996. C
  1997.         SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  1998.         SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  1999.         EPSI(1) = EPSI0(1) + GAMBDA(5)
  2000.         WORK = ISGNTC*RESIST(ICAS)*GAMBDA(ICAS)
  2001.         IF (ICASB.LE.4) THEN
  2002.           WORK = WORK + ISGNTCB*RESIST(ICASB)*GAMBDA(ICASB)
  2003.         ENDIF
  2004.         EPSI(2)=EPSI0(2)+ GAMBDA(6) + DEGRAD*WORK
  2005. C
  2006.         IF (ICASB.LE.4) THEN
  2007.           IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2008.             CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  2009.      &XKISO(1),XHH2,KERRE)
  2010.             CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  2011.      &XKISO(2),DUMM,KERRE)
  2012.           ELSE
  2013.             CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  2014.      &XKISO(1),DUMM,KERRE)
  2015.             CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  2016.      &XKISO(2),XHH2,KERRE)
  2017.           ENDIF
  2018.         ELSE
  2019.           IF (ICASB.EQ.5) THEN
  2020.             IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2021.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2022.      &         ,XHH2,KERRE)
  2023.               XHH2B=XHH2
  2024.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2025.      &         ,DUMM,KERRE)
  2026.             ELSE
  2027.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2028.      &         ,XHH2B,KERRE)
  2029.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2030.      &         ,XHH2,KERRE)
  2031.             ENDIF
  2032.           ELSE
  2033.             IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2034.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2035.      &         ,XHH2,KERRE)
  2036.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2037.      &         ,XHH2B,KERRE)
  2038.             ELSE
  2039.               CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2040.      &         ,DUMM,KERRE)
  2041.               CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2042.      &         ,XHH2,KERRE)
  2043.               XHH2B = XHH2
  2044.             ENDIF
  2045.           ENDIF
  2046.         ENDIF
  2047. C
  2048. C  Calcul des phi
  2049. C
  2050.         PHI1=ISGNTC*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  2051.      &      + SIGHP(IHARD) + RESIST(ICAS))
  2052. C
  2053.         IF (ICASB.LE.4) THEN
  2054.            PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  2055.      &   + SIGHP(IHARDB)+  RESIST(ICASB))
  2056.         ELSE
  2057.            PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  2058.      &    + (XKISOI(ICASB-4)+XHH2IB*(EPSI(ICASB-4)-EPSII(ICASB-4)))
  2059.      &      *RESIST(ICASB))
  2060.         ENDIF
  2061. C
  2062. C   Test sur phi1 et phi2
  2063. C
  2064.         IF ((ABS(PHI2)).GT.(ABS(EPSILO*RESIST(ICASB)))) THEN
  2065.            SOMLAM = GAMBDA(ICAS) + GAMBDA(ICASB)*PHI1/PHI2
  2066.            PHI = PHI1
  2067.         ELSE
  2068.           IF ((ABS(PHI1)).GT.(ABS(EPSILO*RESIST(ICAS)))) THEN
  2069.             SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)*PHI2/PHI1
  2070.             PHI = PHI2
  2071.           ELSE
  2072.             SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)
  2073.             PHI = PHI1
  2074.           ENDIF
  2075.         ENDIF
  2076. C
  2077.         AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  2078.         AINV(2) = (GGG*SOMLAM)/(DEUX*PHI + DEUX*GGG*SOMLAM)
  2079.         AINV(3) = PHI/(PHI + GGG*SOMLAM)
  2080. C
  2081. C   Calcul des contraintes a l'aide de AINV
  2082. C
  2083.         SIGT(1) = AINV(1)*SITPRO(1)+AINV(2)*SITPRO(2)
  2084.         SIGT(2) = AINV(2)*SITPRO(1)+AINV(1)*SITPRO(2)
  2085.         SIGT(4) = AINV(3)*SITPRO(4)
  2086. C
  2087. C   Nouvelles contraintes principales
  2088. C
  2089.         CALL CONPRI ( SIGT , CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  2090. C
  2091.           CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD),RESIST(ICAS),
  2092.      &               ISGNTC,FFF(ICAS))
  2093.           R(1)= + FFF(ICAS)
  2094. C
  2095.         IF (ICASB.LE.4) THEN
  2096.           CALL CRITPE (SIGTP(IHARDB),SIGHP(IHARDB),RESIST(ICASB),
  2097.      &               ISGNTCB,FFF(ICASB))
  2098.           R(2)= + FFF(ICASB)
  2099.         ELSE
  2100.           FG(1) = (XKISOI(ICASB -4)+XHH2IB*(EPSI(ICASB-4)
  2101.      & - EPSII(ICASB-4)))*RESIST(ICASB)
  2102.           CALL CRITPE (SIGTP(IHARDB), XZER, FG(1),ISGNTCB,FFF(ICASB))
  2103.           R(2) = + FFF(ICASB)
  2104.      &         - ISGNTCB*RESIST(ICASB)
  2105.      &         * (XKISO(ICASB-4) - XKISOI(ICASB-4)
  2106.      &            -(EPSI(ICASB-4) - EPSII(ICASB-4))*XHH2IB)
  2107.           IF (ABS(EPSI(ICASB-4) - EPSII(ICASB-4)).LT.EPSILO2) THEN
  2108.             XHH2IB=XHH2B
  2109.           ELSE
  2110.             XHH2IB = (XKISO(ICASB-4) - XKISOI(ICASB-4))/
  2111.      &              (EPSI(ICASB-4) - EPSII(ICASB-4))
  2112.           ENDIF
  2113.         ENDIF
  2114. C
  2115.         IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2116.           FG(1) = XKISO(1)*RESIST(5)
  2117.           CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC,XFFFS)
  2118.           R(3) = XFFFS
  2119.         ELSE
  2120.           FG(1) = (XKISOI(2)+XHH2I*(EPSI(2)-EPSII(2)))*RESIST(6)
  2121.           CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,XFFFS)
  2122.           R(3) = XFFFS + RESIST(6)*
  2123.      &     (XKISO(2) - XKISOI(2)-(EPSI(2) - EPSII(2))*XHH2I)
  2124.           IF (ABS(EPSI(2) - EPSII(2)).LT.EPSILO2) THEN
  2125.             XHH2I=XHH2
  2126.           ELSE
  2127.             XHH2I = (XKISO(2) - XKISOI(2))/
  2128.      &              (EPSI(2) - EPSII(2))
  2129.           ENDIF
  2130.         ENDIF
  2131.           XKISOI(1) = XKISO(1)
  2132.           EPSII(1) = EPSI(1)
  2133.           XKISOI(2) = XKISO(2)
  2134.           EPSII(2) = EPSI(2)
  2135. C
  2136. C  Test de convergence
  2137. C
  2138.         IF (ICAS.LE.2) THEN
  2139.           IF ((ABS(FFF(ICAS)).LT.EPSIL(ICAS)).AND.
  2140.      &      (ABS(R(1)).LT.EPSIL(ICAS)).AND.
  2141.      &      (ABS(FFF(ICASB)).LT.EPSIL(ICASB)).AND.
  2142.      &      (ABS(R(2)).LT.EPSIL(ICASB)).AND.
  2143.      &      (ABS(XFFFS).LT.EPSIL(5)).AND.
  2144.      &      (ABS(R(3)).LT.EPSIL(5))) THEN
  2145.              GOTO 5020
  2146.           ENDIF
  2147.         ELSE
  2148.           IF ((ABS(FFF(ICAS)).LT.EPSIL(ICAS)).AND.
  2149.      &      (ABS(R(1)).LT.EPSIL(ICAS)).AND.
  2150.      &      (ABS(FFF(ICASB)).LT.EPSIL(ICASB)).AND.
  2151.      &      (ABS(R(2)).LT.EPSIL(ICASB)).AND.
  2152.      &      (ABS(XFFFS).LT.EPSIL(6)).AND.
  2153.      &      (ABS(R(3)).LT.EPSIL(6))) THEN
  2154.              GOTO 5020
  2155.           ENDIF
  2156.         ENDIF
  2157.       ENDDO
  2158. C
  2159. C  Pas de convergence
  2160. C
  2161.         SEGSUP,QUASIN
  2162.         GOTO 6000
  2163. C
  2164.  5020 CONTINUE
  2165.         SEGSUP,QUASIN
  2166. C
  2167. C
  2168.        IF (((SPLITT-1.).LE.EPSILO).AND.(SPLITT.GE.(-1.*EPSILO)).AND.
  2169.      &(GAMBDA(ICAS).GE.EPSILO2)
  2170.      &.AND.(GAMBDA(ICASB).GE.EPSILO2)) THEN
  2171. C
  2172. C  On a les contraintes, les variables d'ecrouissage cinemat. et tgphi
  2173. C  Mis a jour des variables d'ecrouissage isotrope ou regles de suivi
  2174. C
  2175.  
  2176.       IF (NCON.EQ.2) THEN
  2177.         IF (ICASB.LE.4) THEN
  2178. C
  2179.            IF (ICASB.LE.2) THEN
  2180.              FG(1) = XKISO(1)*RESIST(5)
  2181.            ELSE
  2182.              FG(1) = XKISO(2)*RESIST(6)
  2183.            ENDIF
  2184.            CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,XFFFSB)
  2185. C
  2186. C  On regarde si on ne s'est pas trompe d'ordre
  2187. C
  2188.           IF (((XFFFSB.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  2189.      &         ((XFFFSB.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  2190.             ICASI = ICAS
  2191.             ISGNTCI = ISGNTC
  2192.             IHARDI = IHARD
  2193.             ICAS = ICASB
  2194.             ISGNTC = ISGNTCB
  2195.             IHARD = IHARDB
  2196.             ICASB = ICASI
  2197.             ISGNTCB = ISGNTCI
  2198.             IHARDB = IHARDI
  2199.             GOTO 4040
  2200.           ELSE
  2201.             IF ((ICASB.LE.2).AND.(ABS(XFFFSB).LT.EPSIL(5))) THEN
  2202.                KONTAC(ICASB) = 1
  2203.                ICASB=5
  2204.             ENDIF
  2205.             IF (((ICASB.EQ.3).OR.(ICASB.EQ.4))
  2206.      &          .AND.(ABS(XFFFSB).LT.EPSIL(6))) THEN
  2207.               KONTAC(ICASB) = 1
  2208.               ICASB=6
  2209.             ENDIF
  2210.           ENDIF
  2211. C
  2212.         ELSE
  2213.           CONTINUE
  2214.         ENDIF
  2215.       ELSE
  2216.         CONTINUE
  2217.       ENDIF
  2218. C
  2219. C Le contact est realise
  2220. C
  2221.        KONTAC(ICAS) = 1
  2222.        IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  2223.      &(RESIST(1)-RESIST(3))) THEN
  2224.          IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  2225.          IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  2226.          IF ((ICASB.EQ.1).OR.(ICASB.EQ.2)) KONTAC(ICASB+2)=0
  2227.          IF ((ICASB.EQ.3).OR.(ICASB.EQ.4)) KONTAC(ICASB-2)=0
  2228.        ELSE
  2229.          KONTAC(1) = 1
  2230.          KONTAC(2) = 1
  2231.          KONTAC(3) = 1
  2232.          KONTAC(4) = 1
  2233.        ENDIF
  2234. C
  2235.        SIGI(1) = SIGT(1)
  2236.        SIGI(2) = SIGT(2)
  2237.        SIGI(4) = SIGT(4)
  2238.        DSIGI(1) = DSIGI(1)*SPLITT
  2239.        DSIGI(2) = DSIGI(2)*SPLITT
  2240.        DSIGI(4) = DSIGI(4)*SPLITT
  2241.        SIGT(1) = SIGI(1) + DSIGI(1)
  2242.        SIGT(2) = SIGI(2) + DSIGI(2)
  2243.        SIGT(4) = SIGI(4) + DSIGI(4)
  2244.        SIGHP0(1) = SIGHP(1)
  2245.        SIGHP0(2) = SIGHP(2)
  2246.        EPSI0(1) = EPSI(1)
  2247.        EPSI0(2) = EPSI(2)
  2248.        NPASS2 = 0
  2249.        CPHI0 = CPHI
  2250.        SPHI0 = SPHI
  2251.        IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2252.           ICAS = 5
  2253.        ELSE
  2254.           ICAS = 6
  2255.        ENDIF
  2256.         IF (ICASB.LE.4) THEN
  2257.           GOTO 4000
  2258.         ELSE
  2259.           IF (ICASB.EQ.ICAS) THEN
  2260. C On passe de 2 mecanismes a 1 mecanisme
  2261.             GOTO 1000
  2262.           ELSE
  2263.             GOTO 4000
  2264.           ENDIF
  2265.         ENDIF
  2266. C
  2267. C
  2268.        ELSE
  2269.  5040     CONTINUE
  2270.           IF ((GAMBDA(ICASB).LT.EPSILO2).AND.(NPASS.LE.3)) THEN
  2271.             SIGHP(1) = SIGHP0(1)
  2272.             SIGHP(2) = SIGHP0(2)
  2273.             SIGT(1) = SIGI(1) + DSIGI(1)
  2274.             SIGT(2) = SIGI(2) + DSIGI(2)
  2275.             SIGT(4) = SIGI(4) + DSIGI(4)
  2276.             IF (NPASS.LT.5) THEN
  2277.               NPASS = NPASS + 1
  2278.               GOTO 3000
  2279.             ELSE
  2280.               GOTO 6000
  2281.             ENDIF
  2282.           ELSE
  2283.             IF ((GAMBDA(ICAS).LT.EPSILO2).AND.(NPASS2.EQ.0)) THEN
  2284.               NPASS2 = 1
  2285.               SIGHP(1) = SIGHP0(1)
  2286.               SIGHP(2) = SIGHP0(2)
  2287.               SIGT(1) = SIGI(1) + DSIGI(1)
  2288.               SIGT(2) = SIGI(2) + DSIGI(2)
  2289.               SIGT(4) = SIGI(4) + DSIGI(4)
  2290.               ICAS = ICASB
  2291.               IHARD = IHARDB
  2292.               ISGNTC = ISGNTCB
  2293.               IF (NPASS.LT.5) THEN
  2294.                 NPASS = NPASS + 1
  2295.                 GOTO 1000
  2296.               ELSE
  2297.                 GOTO 6000
  2298.               ENDIF
  2299.             ELSE
  2300.               GOTO 6000
  2301.             ENDIF
  2302.           ENDIF
  2303.         ENDIF
  2304. C
  2305. C On bissecte
  2306. C
  2307.  6000 CONTINUE
  2308.       NDICHO = NDICHO + 1
  2309.       IF (NDICHO.LE.90) THEN
  2310.         IDICHO = 1
  2311.         DSIGIDI(1) = DSIGIDI(1) + UNDEMI*DSIGI(1)
  2312.         DSIGIDI(2) = DSIGIDI(2) + UNDEMI*DSIGI(2)
  2313.         DSIGIDI(4) = DSIGIDI(4) + UNDEMI*DSIGI(4)
  2314.         DSIGI(1) = DSIGI(1)*UNDEMI
  2315.         DSIGI(2) = DSIGI(2)*UNDEMI
  2316.         DSIGI(4) = DSIGI(4)*UNDEMI
  2317.         CPHI0 = CPHI00
  2318.         SPHI0 = SPHI00
  2319.         GOTO 10
  2320.       ELSE
  2321.         KERRE = 2
  2322.         RETURN
  2323.       ENDIF
  2324. C
  2325. C=================================================================
  2326. C      MIS A JOUR DE TOUTES LES VARIABLES ET CONTRAINTES AVANT
  2327. C      DE SORTIR
  2328. C=================================================================
  2329.  9999 CONTINUE
  2330. C  Cas de la bissection
  2331.       IF (IDICHO.EQ.1) THEN
  2332.         IDICHO = 0
  2333.         NDICHO = 0
  2334.         NPASS = 0
  2335.         DSIGI(1) = DSIGIDI(1)
  2336.         DSIGI(2) = DSIGIDI(2)
  2337.         DSIGI(4) = DSIGIDI(4)
  2338.         SIGI(1) = SIGT(1)
  2339.         SIGI(2) = SIGT(2)
  2340.         SIGI(4) = SIGT(4)
  2341.         DSIGIDI(1) = XZER
  2342.         DSIGIDI(2) = XZER
  2343.         DSIGIDI(4) = XZER
  2344.         SIGHP0(1) = SIGHP(1)
  2345.         SIGHP0(2) = SIGHP(2)
  2346.         EPSI0(1) = EPSI(1)
  2347.         EPSI0(2) = EPSI(2)
  2348.         CPHI0 = CPHI
  2349.         CPHI00 = CPHI0
  2350.         SPHI0 = SPHI
  2351.         SPHI00 = SPHI0
  2352.         KONTA0(1) = KONTAC(1)
  2353.         KONTA0(2) = KONTAC(2)
  2354.         KONTA0(3) = KONTAC(3)
  2355.         KONTA0(4) = KONTAC(4)
  2356.         GOTO 10
  2357.       ENDIF
  2358. C  STRESS
  2359.       SIGF(1) = SIGT(1)
  2360.       SIGF(2) = SIGT(2)
  2361.       SIGF(4) = SIGT(4)
  2362.       SIGF(3) = XZER
  2363. C  PLASTIC STRAIN
  2364.       CINV(1) = UN/YOUN
  2365.       CINV(2) = -XNU/YOUN
  2366.       CINV(3) = UN/GGG
  2367.       DEFP(1) =CINV(1)*(SIG0(1) + DSIGT(1) - SIGT(1))
  2368.      &        +CINV(2)*(SIG0(2) + DSIGT(2) - SIGT(2))
  2369.       DEFP(2) =CINV(2)*(SIG0(1) + DSIGT(1) - SIGT(1))
  2370.      &        +CINV(1)*(SIG0(2) + DSIGT(2) - SIGT(2))
  2371.       DEFP(4) =CINV(3)*(SIG0(4) + DSIGT(4) - SIGT(4))
  2372.       DEFP(3) = XZER
  2373. C
  2374. C Cas ou la grande surface est a l'interieur de la petite
  2375. C
  2376.       IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).LE.
  2377.      &(RESIST(1)-RESIST(3))) THEN
  2378.         KONTAC(1) = 1
  2379.         KONTAC(2) = 1
  2380.         KONTAC(3) = 1
  2381.         KONTAC(4) = 1
  2382.       ENDIF
  2383.       VARF(1) = SIGHP(1)
  2384.       VARF(2) = SIGHP(2)
  2385.       VARF(3) = EPSI(1)
  2386.       VARF(4) = EPSI(2)
  2387.       VARF(5) = CPHI - 1.D0
  2388.       VARF(6) = SPHI
  2389.       VARF(7) = KONTAC(1)
  2390.       VARF(8) = KONTAC(2)
  2391.       VARF(9) = KONTAC(3)
  2392.       VARF(10) = KONTAC(4)
  2393.       VARF(11) = ICAS
  2394.       VARF(12) = ICASB
  2395.       VARF(13) = GAMBDA(ICAS)
  2396.       VARF(14) = GAMBDA(ICASB)
  2397. C=====================================================
  2398. C            FIN DE LA ROUTINE DU MODELE
  2399. C======================================================
  2400.       RETURN
  2401.       END
  2402.  
  2403.  
  2404.  
  2405.  
  2406.  
  2407.  
  2408.  
  2409.  
  2410.  
  2411.  
  2412.  
  2413.  

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