Télécharger la plaquette Cast3M


fr en

Télécharger mistrl.eso

Retour à la liste

Numérotation des lignes : 

mistrl

  1. C MISTRL    SOURCE    PV        11/03/07    21:17:30     6885
  2.  
  3. C-----------------------------------------------------------------------
  4. C Sous-programme d'interface entre CAST3M et module MISTRAL version 2.0
  5. C-----------------------------------------------------------------------
  6.       SUBROUTINE MISTRL (TEMP0,T0,FI0, SIG0, VAR0,
  7.      &                   IFOURB, NSTRS,
  8.      &                   DT, TF,FIF, DEPST, XMAT,TXR,IDIM,
  9.      &                   PDILT, NDPI,NDVP,NXX,NPSI,
  10.      &                   PCOHI,PECOU,PEDIR,PRVCE,PECRX,PDVDI, PCROI,
  11.      &                   NPINCR,PINCR,
  12.      &                   SIGF,VARF,EPINF)
  13. C-----------------------------------------------------------------------
  14. C     Entrees :
  15. C     ---------
  16. C     TEMP0    : instant de debut du pas de temps
  17. C     T0       : temperature au debut du pas de temps
  18. C     FI0      : flux de neutrons rapides au debut du pas de temps
  19. C     SIG0()   : contraintes au debut du pas de temps *
  20. C     VAR0()   : variables internes materiau au debut du pas de temps **
  21. C     IFOURB   : type de calcul CAST3M :
  22. C                =  2 : tridimensionnel (3D)
  23. C                =  0 : axisymetrie (2D)
  24. C                = -1 : deformations planes (2D)
  25. C                = -2 : contraintes planes (2D)
  26. C                = -3 : deformations planes generalisees (2D)
  27. C                =  3 a 11 : unidimensionnel plan (massif 1D)
  28. C                = 12 a 14 : unidimensionnel axisymetrique (massif 1D)
  29. C     NSTRS    : nombre de composantes de contraintes et de deformations
  30. C                pour CAST3M
  31. C     DT       : pas de temps
  32. C     TF       : temperature a la fin du pas de temps
  33. C     FIF      : flux de neutrons rapides a la fin du pas de temps
  34. C     DEPST()  : increment impose de deformation totale *
  35. C     XMAT()   : caracteristiques materiau selon CAST3M
  36. C     TXR(l,c) : composante l du vecteur c de la base locale de CAST3M (element
  37. C                fini) dans la base generale de CAST3M (ensemble du modele)
  38. C     IDIM     : dimension de l'espace dans CAST3M : 1 a 3 (pour DIMENSION TXR)
  39. C     PDILT    : tableau des parametres relatifs a la dilatation thermique
  40. C     NDPI     : nombre de deformations plastiques instantanees (0 ou 1)
  41. C     NDVP     : nombre de deformations viscoplastiques (0 a 3)
  42. C     NXX      : nombre de contraintes internes (0 a 3)
  43. C     NPSI     : existence (1) ou non (0) de variable de durcissement
  44. C                d'irradiation differente de la fluence neutronique
  45. C     tableaux des parametres du modele :
  46. C     PCOHI    : relatif aux coefficients de Hill
  47. C     PECOU    : relatif a l'ecoulement viscoplastique
  48. C     PEDIR    : relatif a la contrainte seuil
  49. C     PRVCE    : relatif aux lois d'evolution des deformations equivalentes
  50. C     PECRX    : relatif aux lois d'evolution des contraintes internes
  51. C     PDVDI    : relatif au durcissement d'irradiation
  52. C     PCROI    : relatif a la croissance d'irradiation
  53. C     NPINCR   : nombre de parametres de PINCR
  54. C     PINCR    : tableau des increments maximaux autorises
  55. C
  56. C     Sorties :
  57. C     ---------
  58. C     SIGF()   : contraintes a la fin du pas de temps *
  59. C     VARF()   : variables internes materiau a la fin du pas de temps **
  60. C     EPINF()  : deformation inelastique globale (somme des deformations
  61. C                plastiques de toutes natures) a la fin du pas de temps *
  62. C
  63. C     *  : dans la base generale de CAST3M
  64. C          (donc composantes 4, 5 et 6 des deformations = glissements)
  65. C     ** : dans la base generale de CAST3M ou dans base d'orthotropie MISTRAL
  66. C          pour les composantes de deformations et contraintes internes
  67. C          selon valeur de ICBASE extrait de XMAT
  68. C-----------------------------------------------------------------------
  69.  
  70.       IMPLICIT INTEGER (I-N)
  71.       IMPLICIT REAL*8 (A-H, O-Z)
  72.  
  73. C     Nombre total de variables internes :
  74. C        NIJ*(4+NDVPMA)+2*(1+NDVPMA)+NIJ*NXXMAX+8+2*NIJ+1+NDVPMA = 101
  75.       PARAMETER ( NIJ = 6 , NDVPMA = 4 , NXXMAX = 3 )
  76.       PARAMETER ( FEPSIJ = 2. )
  77.       PARAMETER ( NPCELM = 19 )
  78.  
  79.       DIMENSION SIG0(1:*),VAR0(1:*),DEPST(1:*),XMAT(1:*),
  80.      &          TXR(1:IDIM,1:*)
  81.       DIMENSION PDILT(1:*),PCOHI(1:*),PECOU(1:*),
  82.      &          PEDIR(1:*),PRVCE(1:*),PECRX(1:*),PDVDI(1:*),
  83.      &          PCROI(1:*),PINCR(1:*)
  84.       DIMENSION PCOEL(1:NPCELM)
  85.       DIMENSION PBLOC(1:3,1:3),VEC(1:3,1:3),
  86.      &          IP(1:3),SENS(1:3),PVEC(1:3,1:3),
  87.      &          PBPO(1:3,1:3),PINV(1:3,1:3),IC(6),JC(6),
  88.      &          VECT(1:3,1:3),VECTP(1:3,1:3),
  89.      &          PROMAT(1:3,1:3)
  90.       DIMENSION SIG(1:6),EPS(1:6),EPSTH(1:6),EPSEL(1:6),EPSCR(1:6),
  91.      &          EPSP(1:6,0:NDVPMA),EPSE(0:NDVPMA),EPSEQ(0:NDVPMA),
  92.      &          XX(1:6,1:NXXMAX)
  93.       DIMENSION CSIG(1:6,1:6),CEPS(1:6,1:6),CSIEP0(1:6)
  94.       DIMENSION SEQ(0:NDVPMA),RR(0:NDVPMA)
  95.       DIMENSION SIGP(1:6),VEPS(1:6),VEPSTH(1:6),VEPSEL(1:6),VEPSCR(1:6),
  96.      &          VEPSP(1:6,0:NDVPMA),VEPSE(0:NDVPMA),
  97.      &          EPSEQP(0:NDVPMA),XXP(1:6,1:NXXMAX)
  98.       DIMENSION DELSIG(1:6),
  99.      &          DELEPS(1:6),DELETH(1:6),DELEEL(1:6),DELECR(1:6),
  100.      &          DELEP(1:6,0:NDVPMA),DELEE(0:NDVPMA),DELEEQ(0:NDVPMA),
  101.      &          DELXX(1:6,1:NXXMAX)
  102.       DIMENSION VEPSIN(1:6)
  103.       REAL*8 NU120,NU230,NU310,MU120,MU130,MU230
  104.       CHARACTER FICH*1
  105.       DIMENSION SIGF(1:*),VARF(1:*),EPINF(1:*)
  106. C     ------------------------------------------------------------------
  107.       DATA IC / 1 , 2 , 3 , 1 , 1 , 2 /
  108.       DATA JC / 1 , 2 , 3 , 2 , 3 , 3 /
  109.       DATA FICH / ' ' /
  110. C     ------------------------------------------------------------------
  111. C     Variables diverses
  112. C
  113.       T = TEMP0
  114.       TT = T0
  115.       FI = FI0
  116.       IF (DT.EQ.0.D0) THEN
  117.         DTI = 1.D-20
  118.       ELSE
  119.         DTI = DT
  120.       END IF
  121.       TTP = (TF-T0)/DTI
  122.       FIP = (FIF-FI0)/DTI
  123. C     ------------------------------------------------------------------
  124. C     Determination de la matrice PBPO de passage de la base generale
  125. C     CAST3M a la base principale d'orthotropie pour MISTRAL
  126. C     et de son inverse PINV
  127. C
  128.       IF (IFOURB.EQ.2) THEN
  129.         PBLOC(1,1) = TXR(1,1)
  130.         PBLOC(1,2) = TXR(2,1)
  131.         PBLOC(1,3) = TXR(3,1)
  132.         PBLOC(2,1) = TXR(1,2)
  133.         PBLOC(2,2) = TXR(2,2)
  134.         PBLOC(2,3) = TXR(3,2)
  135.         PBLOC(3,1) = TXR(1,3)
  136.         PBLOC(3,2) = TXR(2,3)
  137.         PBLOC(3,3) = TXR(3,3)
  138.         VEC(1,1) = XMAT(10)
  139.         VEC(1,2) = XMAT(11)
  140.         VEC(1,3) = XMAT(12)
  141.         VEC(2,1) = XMAT(13)
  142.         VEC(2,2) = XMAT(14)
  143.         VEC(2,3) = XMAT(15)
  144.         VEC(3,1) = VEC(1,2)*VEC(2,3)-VEC(1,3)*VEC(2,2)
  145.         VEC(3,2) = VEC(1,3)*VEC(2,1)-VEC(1,1)*VEC(2,3)
  146.         VEC(3,3) = VEC(1,1)*VEC(2,2)-VEC(1,2)*VEC(2,1)
  147.         SENIP1 = XMAT(16)
  148.         SENIP2 = XMAT(17)
  149.       ELSE IF ((IFOURB.GE.-3).AND.(IFOURB.LE.0)) THEN
  150.         PBLOC(1,1) = TXR(1,1)
  151.         PBLOC(1,2) = TXR(2,1)
  152.         PBLOC(1,3) = 0.
  153.         PBLOC(2,1) = TXR(1,2)
  154.         PBLOC(2,2) = TXR(2,2)
  155.         PBLOC(2,3) = 0.
  156.         PBLOC(3,1) = 0.
  157.         PBLOC(3,2) = 0.
  158.         PBLOC(3,3) = 1.
  159.         IF (IFOURB.EQ.-2) THEN
  160.           VEC(1,1) = XMAT(5)
  161.           VEC(1,2) = XMAT(6)
  162.         ELSE
  163.           VEC(1,1) = XMAT(8)
  164.           VEC(1,2) = XMAT(9)
  165.         END IF
  166.         VEC(1,3) = 0.
  167.         VEC(3,1) = 0.
  168.         VEC(3,2) = 0.
  169.         VEC(3,3) = 1.
  170. C*      VEC(2,1) = VEC(3,2)*VEC(1,3)-VEC(3,3)*VEC(1,2)
  171. C*      VEC(2,2) = VEC(3,3)*VEC(1,1)-VEC(3,1)*VEC(1,3)
  172. C*      VEC(2,3) = VEC(3,1)*VEC(1,2)-VEC(3,2)*VEC(1,1)
  173.         VEC(2,1) = -VEC(1,2)
  174.         VEC(2,2) = VEC(1,1)
  175.         VEC(2,3) = 0.
  176.         SENIP1 = XMAT(10)
  177.         SENIP2 = XMAT(11)
  178.       ELSE IF ((IFOURB.GE.3).AND.(IFOURB.LE.14)) THEN
  179.         PBLOC(1,1) = 1.
  180.         PBLOC(1,2) = 0.
  181.         PBLOC(1,3) = 0.
  182.         PBLOC(2,1) = 0.
  183.         PBLOC(2,2) = 1.
  184.         PBLOC(2,3) = 0.
  185.         PBLOC(3,1) = 0.
  186.         PBLOC(3,2) = 0.
  187.         PBLOC(3,3) = 1.
  188.         VEC(1,1) = 1.
  189.         VEC(1,2) = 0.
  190.         VEC(1,3) = 0.
  191.         VEC(2,1) = 0.
  192.         VEC(2,2) = 1.
  193.         VEC(2,3) = 0.
  194.         VEC(3,1) = 0.
  195.         VEC(3,2) = 0.
  196.         VEC(3,3) = 1.
  197.         SENIP1 = XMAT(7)
  198.         SENIP2 = XMAT(8)
  199.       END IF
  200.       IP(1) = ABS(nint(SENIP1))
  201.       IP(2) = ABS(nint(SENIP2))
  202.       SENS(1) = SENIP1/ABS(SENIP1)
  203.       SENS(2) = SENIP2/ABS(SENIP2)
  204.       DO I = 1,2
  205.         DO J = 1,3
  206.           IF (J.EQ.IP(I)) THEN
  207.             PVEC(I,J) = SENS(I)
  208.           ELSE
  209.             PVEC(I,J) = 0.
  210.           END IF
  211.         END DO
  212.       END DO
  213.       PVEC(3,1) = PVEC(1,2)*PVEC(2,3)-PVEC(1,3)*PVEC(2,2)
  214.       PVEC(3,2) = PVEC(1,3)*PVEC(2,1)-PVEC(1,1)*PVEC(2,3)
  215.       PVEC(3,3) = PVEC(1,1)*PVEC(2,2)-PVEC(1,2)*PVEC(2,1)
  216.       CALL PRODMA(PBLOC,VEC,3,3,3, PROMAT)
  217.       CALL PRODMA(PROMAT,PVEC,3,3,3, PBPO)
  218.       DO I = 1,3
  219.         DO J = 1,3
  220.           PINV(I,J) = PBPO(J,I)
  221.         END DO
  222.       END DO
  223. C     ------------------------------------------------------------------
  224. C     Changement de base pour contraintes debut de pas de temps
  225. C
  226.       IF (IFOURB.EQ.2) THEN
  227.         VECT(1,1) = SIG0(1)
  228.         VECT(2,1) = SIG0(4)
  229.         VECT(3,1) = SIG0(5)
  230.         VECT(1,2) = SIG0(4)
  231.         VECT(2,2) = SIG0(2)
  232.         VECT(3,2) = SIG0(6)
  233.         VECT(1,3) = SIG0(5)
  234.         VECT(2,3) = SIG0(6)
  235.         VECT(3,3) = SIG0(3)
  236.       ELSE IF ((IFOURB.GE.-3).AND.(IFOURB.LE.0)) THEN
  237.         VECT(1,1) = SIG0(1)
  238.         VECT(2,1) = SIG0(4)
  239.         VECT(3,1) = 0.
  240.         VECT(1,2) = SIG0(4)
  241.         VECT(2,2) = SIG0(2)
  242.         VECT(3,2) = 0.
  243.         VECT(1,3) = 0.
  244.         VECT(2,3) = 0.
  245.         IF (IFOURB.EQ.-2) THEN
  246.           VECT(3,3) = 0.
  247.         ELSE
  248.           VECT(3,3) = SIG0(3)
  249.         END IF
  250.       ELSE IF ((IFOURB.GE.3).AND.(IFOURB.LE.14)) THEN
  251.         VECT(1,1) = SIG0(1)
  252.         VECT(2,1) = 0.
  253.         VECT(3,1) = 0.
  254.         IF ((IFOURB.EQ.5).OR.(IFOURB.EQ.6).OR.(IFOURB.EQ.10).OR.
  255.      &      (IFOURB.EQ.13)) THEN
  256.           VECT(2,2) = 0.
  257.         ELSE
  258.           VECT(2,2) = SIG0(2)
  259.         END IF
  260.         VECT(3,2) = 0.
  261.         VECT(1,3) = 0.
  262.         VECT(2,3) = 0.
  263.         IF ((IFOURB.EQ.4).OR.(IFOURB.EQ.6).OR.(IFOURB.EQ.8)) THEN
  264.           VECT(3,3) = 0.
  265.         ELSE
  266.           VECT(3,3) = SIG0(3)
  267.         END IF
  268.       END IF
  269.       CALL PRODMA(VECT,PBPO,3,3,3, PROMAT)
  270.       CALL PRODMA(PINV,PROMAT,3,3,3, VECTP)
  271.       SIG(1) = VECTP(1,1)
  272.       SIG(2) = VECTP(2,2)
  273.       SIG(3) = VECTP(3,3)
  274.       SIG(4) = VECTP(1,2)
  275.       SIG(5) = VECTP(1,3)
  276.       SIG(6) = VECTP(2,3)
  277. C     ------------------------------------------------------------------
  278. C     Coefficients des equations entre derivees des contraintes
  279. C     et vitesses de deformations
  280. C
  281.       IF (IFOURB.EQ.2) THEN
  282.         CSIEP0(1) = DEPST(1)/DTI
  283.         CSIEP0(2) = DEPST(2)/DTI
  284.         CSIEP0(3) = DEPST(3)/DTI
  285.         CSIEP0(4) = DEPST(4)/FEPSIJ/DTI
  286.         CSIEP0(5) = DEPST(5)/FEPSIJ/DTI
  287.         CSIEP0(6) = DEPST(6)/FEPSIJ/DTI
  288.         DO IJ = 1,6
  289.           DO IPJP = 1,3
  290.             CSIG(IPJP,IJ) = 0.
  291.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  292.           END DO
  293.           DO IPJP = 4,6
  294.             CSIG(IPJP,IJ) = 0.
  295.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  296.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  297.           END DO
  298.         END DO
  299.       ELSE IF ((IFOURB.GE.-3).AND.(IFOURB.LE.0)) THEN
  300.         IF (IFOURB.EQ.-2) THEN
  301.           CSIEP0(3) = 0.
  302.           IJ = 3
  303.           DO IPJP = 1,3
  304.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  305.             CEPS(IPJP,IJ) = 0.
  306.           END DO
  307.           DO IPJP = 4,6
  308.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  309.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  310.             CEPS(IPJP,IJ) = 0.
  311.           END DO
  312.         ELSE
  313.           CSIEP0(3) = DEPST(3)/DTI
  314.           IJ = 3
  315.           DO IPJP = 1,3
  316.             CSIG(IPJP,IJ) = 0.
  317.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  318.           END DO
  319.           DO IPJP = 4,6
  320.             CSIG(IPJP,IJ) = 0.
  321.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  322.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  323.           END DO
  324.         END IF
  325.         CSIEP0(1) = DEPST(1)/DTI
  326.         CSIEP0(2) = DEPST(2)/DTI
  327.         CSIEP0(4) = DEPST(4)/FEPSIJ/DTI
  328.         CSIEP0(5) = 0.
  329.         CSIEP0(6) = 0.
  330.         DO IJ = 1,2
  331.           DO IPJP = 1,3
  332.             CSIG(IPJP,IJ) = 0.
  333.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  334.           END DO
  335.           DO IPJP = 4,6
  336.             CSIG(IPJP,IJ) = 0.
  337.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  338.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  339.           END DO
  340.         END DO
  341.         IJ = 4
  342.         DO IPJP = 1,3
  343.           CSIG(IPJP,IJ) = 0.
  344.           CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  345.         END DO
  346.         DO IPJP = 4,6
  347.           CSIG(IPJP,IJ) = 0.
  348.           CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  349.      &                  + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  350.         END DO
  351.         DO IJ = 5,6
  352.           DO IPJP = 1,3
  353.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  354.             CEPS(IPJP,IJ) = 0.
  355.           END DO
  356.           DO IPJP = 4,6
  357.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  358.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  359.             CEPS(IPJP,IJ) = 0.
  360.           END DO
  361.         END DO
  362.       ELSE IF ((IFOURB.GE.3).AND.(IFOURB.LE.14)) THEN
  363.         IF ((IFOURB.EQ.5).OR.(IFOURB.EQ.10).OR.(IFOURB.EQ.13)) THEN
  364.           CSIEP0(2) = 0.
  365.           CSIEP0(3) = DEPST(3)/DTI
  366.           IJ = 2
  367.           DO IPJP = 1,3
  368.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  369.             CEPS(IPJP,IJ) = 0.
  370.           END DO
  371.           DO IPJP = 4,6
  372.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  373.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  374.             CEPS(IPJP,IJ) = 0.
  375.           END DO
  376.           IJ = 3
  377.           DO IPJP = 1,3
  378.             CSIG(IPJP,IJ) = 0.
  379.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  380.           END DO
  381.           DO IPJP = 4,6
  382.             CSIG(IPJP,IJ) = 0.
  383.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  384.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  385.           END DO
  386.         ELSE IF ((IFOURB.EQ.4).OR.(IFOURB.EQ.8)) THEN
  387.           CSIEP0(2) = DEPST(2)/DTI
  388.           CSIEP0(3) = 0.
  389.           IJ = 2
  390.           DO IPJP = 1,3
  391.             CSIG(IPJP,IJ) = 0.
  392.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  393.           END DO
  394.           DO IPJP = 4,6
  395.             CSIG(IPJP,IJ) = 0.
  396.             CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  397.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  398.           END DO
  399.           IJ = 3
  400.           DO IPJP = 1,3
  401.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  402.             CEPS(IPJP,IJ) = 0.
  403.           END DO
  404.           DO IPJP = 4,6
  405.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  406.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  407.             CEPS(IPJP,IJ) = 0.
  408.           END DO
  409.         ELSE IF (IFOURB.EQ.6) THEN
  410.           CSIEP0(2) = 0.
  411.           CSIEP0(3) = 0.
  412.           DO IJ = 2,3
  413.             DO IPJP = 1,3
  414.               CSIG(IPJP,IJ)=PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  415.               CEPS(IPJP,IJ)=0.
  416.             END DO
  417.             DO IPJP = 4,6
  418.               CSIG(IPJP,IJ)= PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  419.      &                     + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  420.               CEPS(IPJP,IJ) = 0.
  421.             END DO
  422.           END DO
  423.         ELSE
  424.           CSIEP0(2) = DEPST(2)/DTI
  425.           CSIEP0(3) = DEPST(3)/DTI
  426.           DO IJ = 2,3
  427.             DO IPJP = 1,3
  428.               CSIG(IPJP,IJ)= 0.
  429.               CEPS(IPJP,IJ)= PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  430.             END DO
  431.             DO IPJP = 4,6
  432.               CSIG(IPJP,IJ)= 0.
  433.               CEPS(IPJP,IJ)= PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  434.      &                     + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  435.             END DO
  436.           END DO
  437.         END IF
  438.         CSIEP0(1) = DEPST(1)/DTI
  439.         CSIEP0(4) = 0.
  440.         CSIEP0(5) = 0.
  441.         CSIEP0(6) = 0.
  442.         IJ = 1
  443.         DO IPJP = 1,3
  444.           CSIG(IPJP,IJ) = 0.
  445.           CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  446.         END DO
  447.         DO IPJP = 4,6
  448.           CSIG(IPJP,IJ) = 0.
  449.           CEPS(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  450.      &                  + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  451.         END DO
  452.         DO IJ = 4,6
  453.           DO IPJP = 1,3
  454.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  455.             CEPS(IPJP,IJ) = 0.
  456.           END DO
  457.           DO IPJP = 4,6
  458.             CSIG(IPJP,IJ) = PINV(IC(IJ),IC(IPJP))*PINV(JC(IJ),JC(IPJP))
  459.      &                    + PINV(IC(IJ),JC(IPJP))*PINV(JC(IJ),IC(IPJP))
  460.             CEPS(IPJP,IJ) = 0.
  461.           END DO
  462.         END DO
  463.       END IF
  464. C     ------------------------------------------------------------------
  465. C     Recuperation des variables internes au debut du pas de temps
  466. C
  467.       IVAR = NIJ*(3+NDVPMA)
  468.       DO IJ = 1,NIJ
  469.         EPSTH(IJ) = VAR0(IJ)
  470.         EPSEL(IJ) = VAR0(NIJ+IJ)
  471.         EPSCR(IJ) = VAR0(IVAR+IJ)
  472.       END DO
  473.       IVAR = NIJ*(4+NDVPMA)+1
  474.       JVAR = IVAR+NDVPMA+1
  475.       DO IDP = 1-NDPI,NDVP
  476.         KVAR = NIJ*(2+IDP)
  477.         DO IJ = 1,NIJ
  478.           EPSP(IJ,IDP) = VAR0(KVAR+IJ)
  479.         END DO
  480.         EPSE(IDP) = VAR0(IVAR+IDP)
  481.         EPSEQ(IDP) = VAR0(JVAR+IDP)
  482.       END DO
  483.       IVAR = NIJ*(3+NDVPMA)+2*(1+NDVPMA)
  484.       DO IXX = 1,NXX
  485.         JVAR = IVAR+NIJ*IXX
  486.         DO IJ = 1,NIJ
  487.           XX(IJ,IXX) = VAR0(JVAR+IJ)
  488.         END DO
  489.       END DO
  490.       IFVAR = NIJ*(4+NDVPMA)+2*(1+NDVPMA)+NIJ*NXXMAX
  491.       FIT = VAR0(IFVAR+1)
  492.       PSI = VAR0(IFVAR+2)
  493.       DTA = VAR0(IFVAR+3)
  494.       SEQA0 = VAR0(IFVAR+4)
  495.       RRA0 = VAR0(IFVAR+5)
  496.       VEEA0 = VAR0(IFVAR+6)
  497.       KPLAST = nint(VAR0(IFVAR+7))
  498. C     ------------------------------------------------------------------
  499. C     Changement eventuel de base pour deformations de toutes natures
  500. C     et contraintes internes
  501. C
  502.       IF (IFOURB.EQ.2) THEN
  503.         ICBASE = nint(XMAT(18))
  504.       ELSE IF ((IFOURB.GE.-3).AND.(IFOURB.LE.0)) THEN
  505.         ICBASE = nint(XMAT(12))
  506.       ELSE IF ((IFOURB.GE.3).AND.(IFOURB.LE.14)) THEN
  507.         ICBASE = nint(XMAT(9))
  508.       END IF
  509.       IF (ICBASE.GE.1) THEN
  510.         VECT(1,1) = EPSTH(1)
  511.         VECT(2,2) = EPSTH(2)
  512.         VECT(3,3) = EPSTH(3)
  513.         VECT(1,2) = EPSTH(4)/FEPSIJ
  514.         VECT(2,1) = VECT(1,2)
  515.         VECT(1,3) = EPSTH(5)/FEPSIJ
  516.         VECT(3,1) = VECT(1,3)
  517.         VECT(2,3) = EPSTH(6)/FEPSIJ
  518.         VECT(3,2) = VECT(2,3)
  519.         CALL PRODMA(VECT,PBPO,3,3,3, PROMAT)
  520.         CALL PRODMA(PINV,PROMAT,3,3,3, VECTP)
  521.         EPSTH(1) = VECTP(1,1)
  522.         EPSTH(2) = VECTP(2,2)
  523.         EPSTH(3) = VECTP(3,3)
  524.         EPSTH(4) = VECTP(1,2)
  525.         EPSTH(5) = VECTP(1,3)
  526.         EPSTH(6) = VECTP(2,3)
  527.         VECT(1,1) = EPSEL(1)
  528.         VECT(2,2) = EPSEL(2)
  529.         VECT(3,3) = EPSEL(3)
  530.         VECT(1,2) = EPSEL(4)/FEPSIJ
  531.         VECT(2,1) = VECT(1,2)
  532.         VECT(1,3) = EPSEL(5)/FEPSIJ
  533.         VECT(3,1) = VECT(1,3)
  534.         VECT(2,3) = EPSEL(6)/FEPSIJ
  535.         VECT(3,2) = VECT(2,3)
  536.         CALL PRODMA(VECT,PBPO,3,3,3, PROMAT)
  537.         CALL PRODMA(PINV,PROMAT,3,3,3, VECTP)
  538.         EPSEL(1) = VECTP(1,1)
  539.         EPSEL(2) = VECTP(2,2)
  540.         EPSEL(3) = VECTP(3,3)
  541.         EPSEL(4) = VECTP(1,2)
  542.         EPSEL(5) = VECTP(1,3)
  543.         EPSEL(6) = VECTP(2,3)
  544.         VECT(1,1) = EPSCR(1)
  545.         VECT(2,2) = EPSCR(2)
  546.         VECT(3,3) = EPSCR(3)
  547.         VECT(1,2) = EPSCR(4)/FEPSIJ
  548.         VECT(2,1) = VECT(1,2)
  549.         VECT(1,3) = EPSCR(5)/FEPSIJ
  550.         VECT(3,1) = VECT(1,3)
  551.         VECT(2,3) = EPSCR(6)/FEPSIJ
  552.         VECT(3,2) = VECT(2,3)
  553.         CALL PRODMA(VECT,PBPO,3,3,3, PROMAT)
  554.         CALL PRODMA(PINV,PROMAT,3,3,3, VECTP)
  555.         EPSCR(1) = VECTP(1,1)
  556.         EPSCR(2) = VECTP(2,2)
  557.         EPSCR(3) = VECTP(3,3)
  558.         EPSCR(4) = VECTP(1,2)
  559.         EPSCR(5) = VECTP(1,3)
  560.         EPSCR(6) = VECTP(2,3)
  561.         DO IDP = 1-NDPI,NDVP
  562.           VECT(1,1) = EPSP(1,IDP)
  563.           VECT(2,2) = EPSP(2,IDP)
  564.           VECT(3,3) = EPSP(3,IDP)
  565.           VECT(1,2) = EPSP(4,IDP)/FEPSIJ
  566.           VECT(2,1) = VECT(1,2)
  567.           VECT(1,3) = EPSP(5,IDP)/FEPSIJ
  568.           VECT(3,1) = VECT(1,3)
  569.           VECT(2,3) = EPSP(6,IDP)/FEPSIJ
  570.           VECT(3,2) = VECT(2,3)
  571.           CALL PRODMA(VECT,PBPO,3,3,3, PROMAT)
  572.           CALL PRODMA(PINV,PROMAT,3,3,3, VECTP)
  573.           EPSP(1,IDP) = VECTP(1,1)
  574.           EPSP(2,IDP) = VECTP(2,2)
  575.           EPSP(3,IDP) = VECTP(3,3)
  576.           EPSP(4,IDP) = VECTP(1,2)
  577.           EPSP(5,IDP) = VECTP(1,3)
  578.           EPSP(6,IDP) = VECTP(2,3)
  579.         END DO
  580.         IF (ICBASE.EQ.2) THEN
  581.           DO IXX = 1,NXX
  582.             VECT(1,1) = XX(1,IXX)
  583.             VECT(2,2) = XX(2,IXX)
  584.             VECT(3,3) = XX(3,IXX)
  585.             VECT(1,2) = XX(4,IXX)
  586.             VECT(2,1) = XX(4,IXX)
  587.             VECT(1,3) = XX(5,IXX)
  588.             VECT(3,1) = XX(5,IXX)
  589.             VECT(2,3) = XX(6,IXX)
  590.             VECT(3,2) = XX(6,IXX)
  591.             CALL PRODMA(VECT,PBPO,3,3,3, PROMAT)
  592.             CALL PRODMA(PINV,PROMAT,3,3,3, VECTP)
  593.             XX(1,IXX) = VECTP(1,1)
  594.             XX(2,IXX) = VECTP(2,2)
  595.             XX(3,IXX) = VECTP(3,3)
  596.             XX(4,IXX) = VECTP(1,2)
  597.             XX(5,IXX) = VECTP(1,3)
  598.             XX(6,IXX) = VECTP(2,3)
  599.           END DO
  600.         END IF
  601.       END IF
  602. C     ------------------------------------------------------------------
  603. C     Variables et parametres divers(es)
  604. C
  605.       DO IJ = 1,NIJ
  606.         EPS(IJ) = EPSTH(IJ)+EPSEL(IJ)+EPSCR(IJ)
  607.         DO IDP = 1-NDPI,NDVP
  608.           EPS(IJ) = EPS(IJ)+EPSP(IJ,IDP)
  609.         END DO
  610.       END DO
  611.       DTMI = 0.
  612.       DTMAX = DTI
  613.       TMAX = 1.D+4*DTI
  614.       NPAS = 0
  615.       NCHPL = 0
  616.       EPSMAX = 10.
  617.       IES = 0
  618.       IUM = 8
  619.       IUI = IUM-2
  620.       ISTOP = 0
  621. C     ------------------------------------------------------------------
  622. C     Passage de la base locale d'orthotropie CAST3M a celle de MISTRAL
  623. C     pour les coefficients d'elasticite
  624. C
  625.       PYG1 = XMAT(1)
  626.       PYG2 = XMAT(2)
  627.       IF (IFOURB.EQ.-2) THEN
  628.         PYG3 = XMAT(7)
  629.         PNU12 = XMAT(3)
  630.         PNU23 = XMAT(8)
  631.         PNU13 = XMAT(9)
  632.         PG12 = XMAT(4)
  633.         PG23 = PG12
  634.         PG13 = PG12
  635.       ELSE
  636.         PYG3 = XMAT(3)
  637.         PNU12 = XMAT(4)
  638.         PNU23 = XMAT(5)
  639.         PNU13 = XMAT(6)
  640.         IF (IFOURB.EQ.2) THEN
  641.           PG12 = XMAT(7)
  642.           PG23 = XMAT(8)
  643.           PG13 = XMAT(9)
  644.         ELSE IF ((IFOURB.GE.-3).AND.(IFOURB.LE.0)) THEN
  645.           PG12 = XMAT(7)
  646.           PG23 = PG12
  647.           PG13 = PG12
  648.         ELSE IF ((IFOURB.GE.3).AND.(IFOURB.LE.14)) THEN
  649.           PG12 = 0.
  650.           PG23 = 0.
  651.           PG13 = 0.
  652.         END IF
  653.       END IF
  654. C
  655.       IP(3) = 6-IP(1)-IP(2)
  656.       IF (IP(1).EQ.1) THEN
  657. C Cas 1.1 : IP(1)=1 , IP(2)=2 , IP(3)=3
  658.         IF (IP(2).EQ.2) THEN
  659.           E10 = PYG1
  660.           E20 = PYG2
  661.           E30 = PYG3
  662.           NU120 = PNU12
  663.           NU230 = PNU23
  664.           NU310 = PNU13*PYG3/PYG1
  665.           MU120 = PG12
  666.           MU130 = PG13
  667.           MU230 = PG23
  668. C Cas 1.2 : IP(1)=1 , IP(2)=3 , IP(3)=2
  669.         ELSE
  670.           E10 = PYG1
  671.           E20 = PYG3
  672.           E30 = PYG2
  673.           NU120 = PNU13
  674.           NU230 = PNU23*PYG3/PYG2
  675.           NU310 = PNU12*PYG2/PYG1
  676.           MU120 = PG13
  677.           MU130 = PG12
  678.           MU230 = PG23
  679.         END IF
  680.       ELSE IF (IP(1).EQ.2) THEN
  681. C Cas 2.1 : IP(1)=2 , IP(2)=1 , IP(3)=3
  682.         IF (IP(2).EQ.1) THEN
  683.           E10 = PYG2
  684.           E20 = PYG1
  685.           E30 = PYG3
  686.           NU120 = PNU12*PYG2/PYG1
  687.           NU230 = PNU13
  688.           NU310 = PNU23*PYG3/PYG2
  689.           MU120 = PG12
  690.           MU130 = PG23
  691.           MU230 = PG13
  692. C Cas 2.2 : IP(1)=2 , IP(2)=3 , IP(3)=1
  693.         ELSE
  694.           E10 = PYG2
  695.           E20 = PYG3
  696.           E30 = PYG1
  697.           NU120 = PNU23
  698.           NU230 = PNU13*PYG3/PYG1
  699.           NU310 = PNU12
  700.           MU120 = PG23
  701.           MU130 = PG12
  702.           MU230 = PG13
  703.         END IF
  704.       ELSE IF (IP(1).EQ.3) THEN
  705. C Cas 3.1 : IP(1)=3 , IP(2)=1 , IP(3)=2
  706.        IF (IP(2).EQ.1) THEN
  707.         E10 = PYG3
  708.         E20 = PYG1
  709.         E30 = PYG2
  710.         NU120 = PNU13*PYG3/PYG1
  711.         NU230 = PNU12
  712.         NU310 = PNU23
  713.         MU120 = PG13
  714.         MU130 = PG23
  715.         MU230 = PG12
  716. C Cas 3.2 : IP(1)=3 , IP(2)=2 , IP(3)=1
  717.        ELSE
  718.         E10 = PYG3
  719.         E20 = PYG2
  720.         E30 = PYG1
  721.         NU120 = PNU23*PYG3/PYG2
  722.         NU230 = PNU12*PYG2/PYG1
  723.         NU310 = PNU13
  724.         MU120 = PG23
  725.         MU130 = PG13
  726.         MU230 = PG12
  727.        END IF
  728.       END IF
  729.       PCOEL( 1) = 0
  730.       PCOEL( 2) = E10
  731.       PCOEL( 3) = 0.
  732.       PCOEL( 4) = E20
  733.       PCOEL( 5) = 0.
  734.       PCOEL( 6) = E30
  735.       PCOEL( 7) = 0.
  736.       PCOEL( 8) = NU120
  737.       PCOEL( 9) = 0.
  738.       PCOEL(10) = NU230
  739.       PCOEL(11) = 0.
  740.       PCOEL(12) = NU310
  741.       PCOEL(13) = 0.
  742.       PCOEL(14) = MU120
  743.       PCOEL(15) = 0.
  744.       PCOEL(16) = MU130
  745.       PCOEL(17) = 0.
  746.       PCOEL(18) = MU230
  747.       PCOEL(19) = 0.
  748. C     ------------------------------------------------------------------
  749. C     Appel du module MISTRAL pour calculer l'evolution durant DTI = DT
  750. C     (ou DTI = 1.E-20 si DT = 0.)
  751. C
  752.       CALL MISTR1 (T, TT,FI,FIT, SIG,
  753.      &             EPS,EPSTH,EPSEL,EPSP,EPSCR,EPSE,EPSEQ,XX,PSI,
  754.      &             DTA,SEQA0,RRA0,VEEA0, KPLAST,
  755.      &             DTI, TTP,FIP, CSIG,CEPS,CSIEP0,
  756.      &             PDILT,PCOEL, NDPI,NDVP,NXX,NPSI,
  757.      &             PCOHI,PECOU,PEDIR,PRVCE,PECRX,PDVDI, PCROI,
  758.      &             NPINCR,PINCR, DTMI,DTMAX,TMAX,
  759.      &             NPAS,NCHPL, EPSMAX, IES,FICH,IUM, ISTOP,
  760.      &             DELTT,DELFI,DELFIT, DELSIG,DELEPS,DELETH,
  761.      &             DELEEL,DELEP,DELECR,DELEE,DELEEQ,DELXX,DELPSI)
  762.       IF (ISTOP.NE.0) THEN
  763.        STOP
  764.       END IF
  765. C
  766. C     Vitesses de deformations a la fin du pas de temps
  767. C
  768.       CALL MISTR0 (TT,TTP,FI,FIT, SIG, XX,EPSE,EPSEQ,PSI,KPLAST,
  769.      &             CSIG,CEPS,CSIEP0, IUI,
  770.      &             PDILT,PCOEL, NDPI,NDVP,NXX,
  771.      &             PCOHI,PECOU,PEDIR,PRVCE,PECRX,PDVDI, PCROI,
  772.      &             SEQ,RR,
  773.      &             SIGP, VEPS,VEPSTH,VEPSEL,VEPSP,VEPSCR,
  774.      &             VEPSE,EPSEQP,XXP,PSIP,
  775.      &             ISTOP)
  776.       IF (ISTOP.NE.0) THEN
  777.        STOP
  778.       END IF
  779.       DO IJ = 1,NIJ
  780.         VEPSIN(IJ) = 0.
  781.         DO IDP = 1-NDPI,NDVP
  782.           VEPSIN(IJ) = VEPSIN(IJ)+VEPSP(IJ,IDP)
  783.         END DO
  784.       END DO
  785. C     ------------------------------------------------------------------
  786. C     Changement de base inverse pour contraintes a fin de pas de temps
  787. C
  788.       VECTP(1,1) = SIG(1)
  789.       VECTP(2,2) = SIG(2)
  790.       VECTP(3,3) = SIG(3)
  791.       VECTP(1,2) = SIG(4)
  792.       VECTP(2,1) = SIG(4)
  793.       VECTP(1,3) = SIG(5)
  794.       VECTP(3,1) = SIG(5)
  795.       VECTP(2,3) = SIG(6)
  796.       VECTP(3,2) = SIG(6)
  797.       CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  798.       CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  799.       SIGF(1) = VECT(1,1)
  800.       SIGF(2) = VECT(2,2)
  801.       SIGF(3) = VECT(3,3)
  802.       IF (IFOURB.EQ.2) THEN
  803.         SIGF(4) = VECT(1,2)
  804.         SIGF(5) = VECT(1,3)
  805.         SIGF(6) = VECT(2,3)
  806.       ELSE IF ((IFOURB.GE.-3).AND.(IFOURB.LE.0)) THEN
  807.         SIGF(4) = VECT(1,2)
  808.       END IF
  809. C     ------------------------------------------------------------------
  810. C     Changement eventuel de base pour deformations de toutes natures
  811. C     et contraintes internes a fin de pas de temps
  812. C
  813.       IF (ICBASE.GE.1) THEN
  814.         VECTP(1,1) = EPSTH(1)
  815.         VECTP(2,2) = EPSTH(2)
  816.         VECTP(3,3) = EPSTH(3)
  817.         VECTP(1,2) = EPSTH(4)
  818.         VECTP(2,1) = EPSTH(4)
  819.         VECTP(1,3) = EPSTH(5)
  820.         VECTP(3,1) = EPSTH(5)
  821.         VECTP(2,3) = EPSTH(6)
  822.         VECTP(3,2) = EPSTH(6)
  823.         CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  824.         CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  825.         EPSTH(1) = VECT(1,1)
  826.         EPSTH(2) = VECT(2,2)
  827.         EPSTH(3) = VECT(3,3)
  828.         EPSTH(4) = VECT(1,2)*FEPSIJ
  829.         EPSTH(5) = VECT(1,3)*FEPSIJ
  830.         EPSTH(6) = VECT(2,3)*FEPSIJ
  831.         VECTP(1,1) = EPSEL(1)
  832.         VECTP(2,2) = EPSEL(2)
  833.         VECTP(3,3) = EPSEL(3)
  834.         VECTP(1,2) = EPSEL(4)
  835.         VECTP(2,1) = EPSEL(4)
  836.         VECTP(1,3) = EPSEL(5)
  837.         VECTP(3,1) = EPSEL(5)
  838.         VECTP(2,3) = EPSEL(6)
  839.         VECTP(3,2) = EPSEL(6)
  840.         CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  841.         CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  842.         EPSEL(1) = VECT(1,1)
  843.         EPSEL(2) = VECT(2,2)
  844.         EPSEL(3) = VECT(3,3)
  845.         EPSEL(4) = VECT(1,2)*FEPSIJ
  846.         EPSEL(5) = VECT(1,3)*FEPSIJ
  847.         EPSEL(6) = VECT(2,3)*FEPSIJ
  848.         VECTP(1,1) = EPSCR(1)
  849.         VECTP(2,2) = EPSCR(2)
  850.         VECTP(3,3) = EPSCR(3)
  851.         VECTP(1,2) = EPSCR(4)
  852.         VECTP(2,1) = EPSCR(4)
  853.         VECTP(1,3) = EPSCR(5)
  854.         VECTP(3,1) = EPSCR(5)
  855.         VECTP(2,3) = EPSCR(6)
  856.         VECTP(3,2) = EPSCR(6)
  857.         CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  858.         CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  859.         EPSCR(1) = VECT(1,1)
  860.         EPSCR(2) = VECT(2,2)
  861.         EPSCR(3) = VECT(3,3)
  862.         EPSCR(4) = VECT(1,2)*FEPSIJ
  863.         EPSCR(5) = VECT(1,3)*FEPSIJ
  864.         EPSCR(6) = VECT(2,3)*FEPSIJ
  865.         DO IDP = 1-NDPI,NDVP
  866.           VECTP(1,1) = EPSP(1,IDP)
  867.           VECTP(2,2) = EPSP(2,IDP)
  868.           VECTP(3,3) = EPSP(3,IDP)
  869.           VECTP(1,2) = EPSP(4,IDP)
  870.           VECTP(2,1) = EPSP(4,IDP)
  871.           VECTP(1,3) = EPSP(5,IDP)
  872.           VECTP(3,1) = EPSP(5,IDP)
  873.           VECTP(2,3) = EPSP(6,IDP)
  874.           VECTP(3,2) = EPSP(6,IDP)
  875.           CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  876.           CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  877.           EPSP(1,IDP) = VECT(1,1)
  878.           EPSP(2,IDP) = VECT(2,2)
  879.           EPSP(3,IDP) = VECT(3,3)
  880.           EPSP(4,IDP) = VECT(1,2)*FEPSIJ
  881.           EPSP(5,IDP) = VECT(1,3)*FEPSIJ
  882.           EPSP(6,IDP) = VECT(2,3)*FEPSIJ
  883.         END DO
  884.         VECTP(1,1) = VEPS(1)
  885.         VECTP(2,2) = VEPS(2)
  886.         VECTP(3,3) = VEPS(3)
  887.         VECTP(1,2) = VEPS(4)
  888.         VECTP(2,1) = VEPS(4)
  889.         VECTP(1,3) = VEPS(5)
  890.         VECTP(3,1) = VEPS(5)
  891.         VECTP(2,3) = VEPS(6)
  892.         VECTP(3,2) = VEPS(6)
  893.         CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  894.         CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  895.         VEPS(1) = VECT(1,1)
  896.         VEPS(2) = VECT(2,2)
  897.         VEPS(3) = VECT(3,3)
  898.         VEPS(4) = VECT(1,2)*FEPSIJ
  899.         VEPS(5) = VECT(1,3)*FEPSIJ
  900.         VEPS(6) = VECT(2,3)*FEPSIJ
  901.         IF ((NDPI+NDVP).GT.0) THEN
  902.           VECTP(1,1) = VEPSIN(1)
  903.           VECTP(2,2) = VEPSIN(2)
  904.           VECTP(3,3) = VEPSIN(3)
  905.           VECTP(1,2) = VEPSIN(4)
  906.           VECTP(2,1) = VEPSIN(4)
  907.           VECTP(1,3) = VEPSIN(5)
  908.           VECTP(3,1) = VEPSIN(5)
  909.           VECTP(2,3) = VEPSIN(6)
  910.           VECTP(3,2) = VEPSIN(6)
  911.           CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  912.           CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  913.           VEPSIN(1) = VECT(1,1)
  914.           VEPSIN(2) = VECT(2,2)
  915.           VEPSIN(3) = VECT(3,3)
  916.           VEPSIN(4) = VECT(1,2)*FEPSIJ
  917.           VEPSIN(5) = VECT(1,3)*FEPSIJ
  918.           VEPSIN(6) = VECT(2,3)*FEPSIJ
  919.         END IF
  920.         IF (ICBASE.EQ.2) THEN
  921.           DO IXX = 1,NXX
  922.             VECTP(1,1) = XX(1,IXX)
  923.             VECTP(2,2) = XX(2,IXX)
  924.             VECTP(3,3) = XX(3,IXX)
  925.             VECTP(1,2) = XX(4,IXX)
  926.             VECTP(2,1) = XX(4,IXX)
  927.             VECTP(1,3) = XX(5,IXX)
  928.             VECTP(3,1) = XX(5,IXX)
  929.             VECTP(2,3) = XX(6,IXX)
  930.             VECTP(3,2) = XX(6,IXX)
  931.             CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  932.             CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  933.             XX(1,IXX) = VECT(1,1)
  934.             XX(2,IXX) = VECT(2,2)
  935.             XX(3,IXX) = VECT(3,3)
  936.             XX(4,IXX) = VECT(1,2)
  937.             XX(5,IXX) = VECT(1,3)
  938.             XX(6,IXX) = VECT(2,3)
  939.           END DO
  940.         END IF
  941.       END IF
  942. C     ------------------------------------------------------------------
  943. C     Stockage des variables internes a la fin du pas de temps
  944. C
  945.       IVAR = NIJ*(3+NDVPMA)
  946.       DO IJ = 1,NIJ
  947.         VARF(IJ) = EPSTH(IJ)
  948.         VARF(NIJ+IJ) = EPSEL(IJ)
  949.         VARF(IVAR+IJ) = EPSCR(IJ)
  950.       END DO
  951.       IVAR = NIJ*(4+NDVPMA)+1
  952.       JVAR = IVAR+NDVPMA+1
  953.       LVAR = IFVAR+8
  954.       VARF(LVAR) = 0.
  955.       DO IDP = 1-NDPI,NDVP
  956.         KVAR = NIJ*(2+IDP)
  957.         DO IJ = 1,NIJ
  958.           VARF(KVAR+IJ) = EPSP(IJ,IDP)
  959.         END DO
  960.         VARF(IVAR+IDP) = EPSE(IDP)
  961.         VARF(JVAR+IDP) = EPSEQ(IDP)
  962.         VARF(LVAR) = VARF(LVAR) + EPSE(IDP)
  963.       END DO
  964.       IVAR = NIJ*(3+NDVPMA)+2*(1+NDVPMA)
  965.       DO IXX = 1,NXX
  966.         JVAR = IVAR+NIJ*IXX
  967.         DO IJ = 1,NIJ
  968.           VARF(JVAR+IJ) = XX(IJ,IXX)
  969.         END DO
  970.       END DO
  971.       VARF(IFVAR+1) = FIT
  972.       VARF(IFVAR+2) = PSI
  973.         IF (DT.NE.0.D0) THEN
  974.       VARF(IFVAR+3) = DTA
  975.       VARF(IFVAR+4) = SEQA0
  976.       VARF(IFVAR+5) = RRA0
  977.       VARF(IFVAR+6) = VEEA0
  978.       VARF(IFVAR+7) = KPLAST
  979.         ENDIF
  980.       IVAR = LVAR+NIJ
  981.       DO IJ = 1,NIJ
  982.         VARF(LVAR+IJ) = VEPS(IJ)
  983.         VARF(IVAR+IJ) = VEPSIN(IJ)
  984.       END DO
  985.       IVAR = IVAR+NIJ+1
  986.       DO IDP = 1-NDPI,NDVP
  987.         VARF(IVAR+IDP) = VEPSE(IDP)
  988.       END DO
  989. C     ------------------------------------------------------------------
  990. C     Deformation inelastique = somme des deformations plastiques
  991. C     a la fin du pas de temps
  992. C
  993.       IF (ICBASE.EQ.0) THEN
  994.         VECTP(1,1) = 0.
  995.         VECTP(2,2) = 0.
  996.         VECTP(3,3) = 0.
  997.         VECTP(1,2) = 0.
  998.         VECTP(1,3) = 0.
  999.         VECTP(2,3) = 0.
  1000.         DO IDP = 1-NDPI,NDVP
  1001.           VECTP(1,1) = VECTP(1,1) + EPSP(1,IDP)
  1002.           VECTP(2,2) = VECTP(2,2) + EPSP(2,IDP)
  1003.           VECTP(3,3) = VECTP(3,3) + EPSP(3,IDP)
  1004.           VECTP(1,2) = VECTP(1,2) + EPSP(4,IDP)
  1005.           VECTP(1,3) = VECTP(1,3) + EPSP(5,IDP)
  1006.           VECTP(2,3) = VECTP(2,3) + EPSP(6,IDP)
  1007.         END DO
  1008.         VECTP(2,1) = VECTP(1,2)
  1009.         VECTP(3,1) = VECTP(1,3)
  1010.         VECTP(3,2) = VECTP(2,3)
  1011.         CALL PRODMA(VECTP,PINV,3,3,3, PROMAT)
  1012.         CALL PRODMA(PBPO,PROMAT,3,3,3, VECT)
  1013.         EPINF(1) = VECT(1,1)
  1014.         EPINF(2) = VECT(2,2)
  1015.         EPINF(3) = VECT(3,3)
  1016.         IF (IFOURB.EQ.2) THEN
  1017.           EPINF(4) = VECT(1,2)*FEPSIJ
  1018.           EPINF(5) = VECT(1,3)*FEPSIJ
  1019.           EPINF(6) = VECT(2,3)*FEPSIJ
  1020.         ELSE IF ((IFOURB.GE.-3).AND.(IFOURB.LE.0)) THEN
  1021.           EPINF(4) = VECT(1,2)*FEPSIJ
  1022.         END IF
  1023.       ELSE
  1024.         DO IJ =1,NSTRS
  1025.           EPINF(IJ) = 0.
  1026.           DO IDP = 1-NDPI,NDVP
  1027.             EPINF(IJ) = EPINF(IJ)+EPSP(IJ,IDP)
  1028.           END DO
  1029.         END DO
  1030.       END IF
  1031. C     ------------------------------------------------------------------
  1032.       RETURN
  1033.       END
  1034.  
  1035.  
  1036.  
  1037.  
  1038.  
  1039.  

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