Télécharger la plaquette Cast3M


fr en

Télécharger dp_sol.eso

Retour à la liste

Numérotation des lignes : 

dp_sol

  1. C DP_SOL    SOURCE    CB215821  26/06/25    21:15:07     12581          
  2. C DP_SOIL   SOURCE    AF221230  12/10/04    21:15:18     7520
  3.       SUBROUTINE DP_SOL (XMAT,VAR0,VARF,SIG0,SIGF,DEPST,XCAR)
  4. C
  5. C====&===1=========2=========3=========4=========5=========6=========7==
  6. C Commentaires :   Loi Drucker Prager avec ecrouissage non lineaire et
  7. C                  loi non associée
  8. C Traits       : -
  9. C
  10. C Auteur       : A. Frau (Dr - Ing.) - CEA/DEN/DANS/DM2S/SEMT/EMSI
  11. C====&===1=========2=========3=========4=========5=========6=========7==
  12. C
  13. C----DECLARATION GENERALES----------------------------------------------
  14. C
  15.       IMPLICIT INTEGER(I-N)
  16.       IMPLICIT REAL*8(A-H,O-Z)
  17.  
  18. -INC PPARAM
  19. -INC CCOPTIO
  20. C
  21. C----DECLARATION SEGMENTS-----------------------------------------------
  22. C
  23.       REAL*8 XMAT(*),SIG0(*),SIGF(*),VAR0(*),VARF(*),DEPST(*),XCAR(*)
  24.       REAL*8 XNU,XYOUN,XSY
  25.       REAL*8 XMU,XLAM,XALPHA,XSIGI,XBET1,XGAM1,XDEL1
  26.       REAL*8 XDEEPP(6),XEPP0(6)
  27.       REAL*8 XSIGT(6),XCHIT,XQT,XFHIT,XDSIGT(6),XNDSIGT
  28.       REAL*8 XNUT(6),XPIT(6)
  29.       REAL*8 XTOL1,XVPL1
  30.       REAL*8 XGAMN,XSIGN(6),XEPPN(6),XCHIN,XQN,XFHIN,XDSIGN(6),XNDSIGN
  31.       REAL*8 XGAMN1,XSIGN1(6),XEPPN1(6),XCHIN1,XQN1,XFHIN1
  32.       REAL*8 XDFDG,XDEGAM,XDSIGN1(6),XNDSIGN1,XTSIGN1
  33.       REAL*8 XGAMF,XCHIF,XQF,XEPPF(6),XVALT1,XVALT2
  34. C
  35. C----CONTROLE OPTION CALCUL --------------------------------------------
  36.       IF ((IFOUR.NE.2).AND.(IFOUR.NE.0).AND.(IFOUR.NE.-1)) THEN
  37.         WRITE(IOIMP,10000)
  38. 10000   FORMAT('Option calcul mauvais')
  39.         WRITE(IOIMP,10001)
  40. 10001   FORMAT('  seulement:')
  41.         WRITE(IOIMP,10002)
  42. 10002   FORMAT('    3D MODE PLAN')
  43.         WRITE(IOIMP,10003)
  44. 10003   FORMAT('    2D MODE AXIS')
  45.         WRITE(IOIMP,10004)
  46. 10004   FORMAT('    2D MODE PLAN')
  47.         STOP
  48.       ENDIF
  49. C
  50. C----MISE EN DONNEES----------------------------------------------------
  51. C
  52. C     Module de young
  53.       XYOUN = XMAT(1)
  54. C     Coef Poisson
  55.       XNU  = XMAT(2)
  56. C     alpha
  57.       XALPHA  =  XMAT(5)
  58. C     Sigma_y
  59.       XSY  =  XMAT(6)
  60. C     Sig_inf
  61.       XSIGI  =  XMAT(7)
  62. C     Beta
  63.       XBET1  =  XMAT(8)
  64. C     Gamma
  65.       XGAM1  =  XMAT(9)
  66. C     Delta
  67.       XDEL1  =  XMAT(10)
  68. C
  69. C----VARIABLES INTERNES-------------------------------------------------
  70.  
  71. c      Coef mu elastique
  72.        XMU = (0.5D0)*((XYOUN)/(1.0D0 + XNU))
  73.  
  74. c      Coef lambda elastique
  75.        XLAM = ((XNU)/(1.0D0-((2.0D0)*(XNU))))
  76. C
  77. C----ACTUALISATION DEFORMATION------------------------------------------
  78. C
  79. C
  80. C----PREDICTION ELASTIQUE-----------------------------------------------
  81. C
  82. c      gamma_test = 0 donc le delta_plastique = 0
  83. c      Initialissation Eps_pla_initiales
  84.        DO 11 I1=1,6
  85.          XDEEPP(I1) = 0.0D0
  86.          XEPP0(I1) = VAR0(I1)
  87. 11     CONTINUE
  88. c
  89. c      Calcul de la variable q_test = q_ini -> gamma_test=0
  90.        XCHIT = VAR0(7)
  91. c
  92. c      Calcul de la variable q_test = q_ini -> gamma_test=0
  93.        XQT = ((-1.0D0)*(XSIGI - XSY))*(1.0D0
  94.      &     - EXP(((-1.0D0)*(XBET1))*(XCHIT)))
  95. c
  96. c      Calcul Sigma_test a partir du step initiale
  97.        DO 12 I1=1,6
  98.          IF (I1.LE.3) THEN
  99.            XSIGT(I1) = SIG0(I1) +
  100.      &     ((2.0D0)*(XMU))*(DEPST(I1) + ((XLAM)*(TRACE(DEPST))))
  101.          ELSE
  102.            IF (IFOUR.EQ.2) THEN
  103.              XSIGT(I1) = SIG0(I1) + ((1.0D0)*(XMU))*(DEPST(I1))
  104.            ELSE
  105.              IF (I1.EQ.4) THEN
  106.                XSIGT(I1) = SIG0(I1) + ((1.0D0)*(XMU))*(DEPST(I1))
  107.              ELSE
  108.                XSIGT(I1) = 0.0D0
  109.              ENDIF
  110.            ENDIF
  111.          ENDIF
  112. 12     CONTINUE
  113. C
  114. C      Deviatoire de Sigma_test
  115.        DO 13 I1=1,6
  116.          IF (I1.LE.3) THEN
  117.            XDSIGT(I1) = XSIGT(I1) - ((TRACE(XSIGT)/(3.0D0)))
  118.          ELSE
  119.            IF (IFOUR.EQ.2) THEN
  120.              XDSIGT(I1) = XSIGT(I1)
  121.            ELSE
  122.              IF (I1.EQ.4) THEN
  123.                XDSIGT(I1) = XSIGT(I1)
  124.              ELSE
  125.                XDSIGT(I1) = 0.0D0
  126.              ENDIF
  127.            ENDIF
  128.          ENDIF
  129. 13     CONTINUE
  130. c
  131. c      Calcul de la norme partie deviatoire
  132.        XNDSIGT = 0.0D0
  133.        DO 14 I1=1,6
  134.          IF (I1.LE.3) THEN
  135.            XNDSIGT = XNDSIGT + (XDSIGT(I1)*XDSIGT(I1))
  136.          ELSE
  137.            XNDSIGT = XNDSIGT + ((2.0D0)*(XDSIGT(I1)*XDSIGT(I1)))
  138.          ENDIF
  139. 14     CONTINUE
  140.        XNDSIGT = ((XNDSIGT)**(0.5D0))
  141. C
  142. c      Calcul du tensuer pi
  143.        Xtmp    = 3.0D0**(-0.5D0)
  144.        XPIT(1) = Xtmp
  145.        XPIT(2) = Xtmp
  146.        XPIT(3) = Xtmp
  147.        XPIT(4) = 0.0D0
  148.        XPIT(5) = 0.0D0
  149.        XPIT(6) = 0.0D0
  150.  
  151. C
  152. c      Calcul ((3**0.5)*(trace(sigma_test)))
  153.        XTSIGT = 0.0D0
  154.        DO 23 I1=1,6
  155.          XTSIGT = XTSIGT + (XPIT(I1)*(XSIGT(I1)))
  156. 23     CONTINUE
  157. c
  158. c      Calcul du vecteur nu_test
  159.        IF ((ABS(XNDSIGT/XTSIGT)).LT.(1.0E-8)) THEN
  160.          XNUT(1) = 1.0D0/3.0D0
  161.          XNUT(2) = 1.0D0/3.0D0
  162.          XNUT(3) =-2.0D0/3.0D0
  163.          XNUT(4) = 0.0D0
  164.          XNUT(5) = 0.0D0
  165.          XNUT(6) = 0.0D0
  166.          XVPL1   = 0.0D0
  167.        ELSE
  168.          XNUT(1) = XDSIGT(1)/XNDSIGT
  169.          XNUT(2) = XDSIGT(2)/XNDSIGT
  170.          XNUT(3) = XDSIGT(3)/XNDSIGT
  171.          XNUT(4) = XDSIGT(4)/XNDSIGT
  172.          XNUT(5) = XDSIGT(5)/XNDSIGT
  173.          XNUT(6) = XDSIGT(6)/XNDSIGT
  174.          XVPL1   = 1.0D0
  175.        ENDIF
  176. c      Calcul du critere fhi_test
  177.        XFHIT = XNDSIGT + (((2.0D0)**(0.5))*(XTSIGT))*(XALPHA)
  178.      &  - (((2.0D0/3.0D0)**(0.5D0))*(XSY - XQT))
  179. C
  180. C-----TEST---------------------------------------------
  181. C
  182.        IF (XFHIT.LT.0) THEN
  183.          GOTO 101
  184. c        fhi_test < 0.0D0  pas de plasticité
  185.        ELSE IF (XFHIT.GE.0) THEN
  186.          GOTO 102
  187. c        fhi_test >= 0.0D0  on est en plasticité
  188.        END IF
  189. C
  190. C-----CALCUL ELASTIQUE---------------------------------------------
  191. C
  192.  
  193. 101   CONTINUE
  194. c
  195. c     si fhi_test < 0  gamma=0
  196.        DO 15 I1=1,6
  197.          SIGF(I1) = XSIGT(I1)
  198. 15     CONTINUE
  199. C
  200. c     pas des variations des variables internes
  201.        DO 16 I1=1,8
  202.          IF (I1.EQ.8) THEN
  203.            VARF(I1) = 0.0D0
  204.          ELSE
  205.            VARF(I1) = VAR0(I1)
  206.          ENDIF
  207. 16     CONTINUE
  208.       RETURN
  209. C
  210. C-----CALCUL PLASTIQUE---------------------------------------------
  211. C
  212. 102   CONTINUE
  213. c
  214. c        si fhi_test >= 0  determination de gamma
  215. c        Nombre max interaction
  216.          NINT1 = 10000
  217. C        Tollerance
  218.          XTOL1 = 1.E-8
  219.       XFHIN = XFHIT
  220.       XQN = XQT
  221.       XGAMN = 0.0D0
  222.       XCHIN = XCHIT
  223. c
  224.       DO 17 I1=1,6
  225.         XEPPN(I1) = XEPP0(I1)
  226.         XSIGN(I1) = XSIGT(I1)
  227. 17    CONTINUE
  228. c
  229.       DO 110 II1=1,NINT1
  230.         IF (II1.EQ.NINT1) THEN
  231. c          PRINT*,'!!!!!   ATTENTION   !!!!'
  232. c          PRINT*,'Nombre des iteractions'
  233. c          PRINT*,'       maximales depasse'
  234. c          STOP
  235.         ENDIF
  236.         XDFDG = (((2.0D0)*(XMU*XVPL1))*(XDEL1))
  237.      &  + ((2.0D0)/(3.0D0))*(((XSIGI - XSY)*(
  238.      &  EXP(((-1.0D0)*(XBET1))*(XCHIN))))*(XBET1)) +
  239.      &  ((4.0D0)*((XALPHA)*(XGAM1)))*
  240.      &  (1.0D0 + ((3.0D0)*(XLAM)))*(XMU)
  241.         XDEGAM = XFHIN/XDFDG
  242.         XGAMN1 = XGAMN + XDEGAM
  243.         XCHIN1 = XCHIN + (((2.0D0/3.0D0)**(0.5D0))*(XDEGAM))
  244.         XQN1 = ((-1.0D0)*(XSIGI - XSY))*(1.0D0
  245.      &       - EXP(((-1.0D0)*(XBET1))*(XCHIN1)))
  246.         IF ((XGAMN1).LT.(0.0D0)) THEN
  247.            WRITE(IOIMP,10005)
  248. 10005      FORMAT('Error....multiplicateur negative')
  249.            STOP
  250.         ENDIF
  251.  
  252.         DO 18 I1=1,6
  253.          IF (I1.LE.3) THEN
  254.             XDEEPP(I1) = ((((XDEL1*XVPL1))*(XNUT(I1)))+
  255.      &       (((XGAM1))*((2.0D0)**(0.5D0)))*(XPIT(I1)))*XDEGAM
  256.          ELSE
  257.            IF (IFOUR.EQ.2) THEN
  258.               XDEEPP(I1) = ((((XDEL1*XVPL1))*(XNUT(I1)))+
  259.      &         (((XGAM1))*((2.0D0)**(0.5D0)))*(XPIT(I1)))*XDEGAM
  260.            ELSE
  261.              IF (I1.EQ.4) THEN
  262.                XDEEPP(I1) = ((((XDEL1*XVPL1))*(XNUT(I1)))+
  263.      &          (((XGAM1))*((2.0D0)**(0.5D0)))*(XPIT(I1)))*XDEGAM
  264.              ELSE
  265.                XDEEPP(I1) = 0.0D0
  266.              ENDIF
  267.            ENDIF
  268.          ENDIF
  269.  
  270.          IF (I1.LE.3) THEN
  271.             XEPPN1(I1) = XEPPN(I1) + XDEEPP(I1)
  272.          ELSE
  273.            IF (IFOUR.EQ.2) THEN
  274.               XEPPN1(I1) = XEPPN(I1) + XDEEPP(I1)
  275.            ELSE
  276.              IF (I1.EQ.4) THEN
  277.                 XEPPN1(I1) = XEPPN(I1) + XDEEPP(I1)
  278.              ELSE
  279.                 XEPPN1(I1) = 0.0D0
  280.              ENDIF
  281.            ENDIF
  282.          ENDIF
  283.  
  284. 18      CONTINUE
  285.  
  286.         DO 21 I1=1,6
  287.         IF (I1.LE.3) THEN
  288.             XSIGN1(I1) = XSIGN(I1) -
  289.      &    ((2.0D0)*(XMU))*(XDEEPP(I1)
  290.      &    + ((XLAM)*(TRACE(XDEEPP))))
  291.           ELSE
  292.            IF (IFOUR.EQ.2) THEN
  293.             XSIGN1(I1) = XSIGN(I1) -
  294.      &    ((1.0D0)*(XMU))*(XDEEPP(I1))
  295.            ELSE
  296.              IF (I1.EQ.4) THEN
  297.                XSIGN1(I1) = XSIGN(I1) -
  298.      &       ((1.0D0)*(XMU))*(XDEEPP(I1))
  299.              ELSE
  300.                XSIGN1(I1) = 0.0D0
  301.              ENDIF
  302.            ENDIF
  303.           ENDIF
  304. 21      CONTINUE
  305.  
  306.         DO 19 I1=1,6
  307.           IF (I1.LE.3) THEN
  308.             XDSIGN1(I1) = XSIGN1(I1)
  309.      &       - ((TRACE(XSIGN1)/(3.0D0)))
  310.           ELSE
  311.            IF (IFOUR.EQ.2) THEN
  312.             XDSIGN1(I1) = XSIGN1(I1)
  313.            ELSE
  314.              IF (I1.EQ.4) THEN
  315.                XDSIGN1(I1) = XSIGN1(I1)
  316.              ELSE
  317.                XDSIGN1(I1) = 0.0D0
  318.              ENDIF
  319.            ENDIF
  320.           ENDIF
  321. 19      CONTINUE
  322.  
  323.         XNDSIGN1 = 0.0D0
  324.         DO 20 I1=1,6
  325.  
  326.          IF (I1.LE.3) THEN
  327.            XNDSIGN1 = XNDSIGN1
  328.      &        + (XDSIGN1(I1)*XDSIGN1(I1))
  329.          ELSE
  330.            XNDSIGN1 = XNDSIGN1
  331.      &        + ((2.0D0)*(XDSIGN1(I1)*XDSIGN1(I1)))
  332.          ENDIF
  333. 20      CONTINUE
  334.         XNDSIGN1 = ((XNDSIGN1)**(0.5D0))
  335.  
  336.         XTSIGN1 = XSIGN1(1) + XSIGN1(2) + XSIGN1(3)
  337.         XTSIGN1 = ((XTSIGN1)/((3.0D0)**(0.5D0)))
  338.         XFHIN1 = XNDSIGN1 +
  339.      &     (((2.0D0)**(0.5))*(XTSIGN1))*(XALPHA)
  340.      &   - (((2.0D0/3.0D0)**(0.5D0))*(XSY - XQN1))
  341. c
  342.         IF (ABS((XFHIN1)/(XFHIT)).LT.(XTOL1)) THEN
  343.           GOTO 103
  344.         ELSE
  345.           XFHIN = XFHIN1
  346.           XGAMN = XGAMN1
  347.           XCHIN = XCHIN1
  348.           XQN = XQN1
  349.           XEPPN(1) = XEPPN1(1)
  350.           XEPPN(2) = XEPPN1(2)
  351.           XEPPN(3) = XEPPN1(3)
  352.           XEPPN(4) = XEPPN1(4)
  353.           XEPPN(5) = XEPPN1(5)
  354.           XEPPN(6) = XEPPN1(6)
  355.           XSIGN(1) = XSIGN1(1)
  356.           XSIGN(2) = XSIGN1(2)
  357.           XSIGN(3) = XSIGN1(3)
  358.           XSIGN(4) = XSIGN1(4)
  359.           XSIGN(5) = XSIGN1(5)
  360.           XSIGN(6) = XSIGN1(6)
  361.         ENDIF
  362. 110   CONTINUE
  363.  
  364.  
  365. 103   CONTINUE
  366.       XCHIF = XCHIN1
  367.       XGAMF = XGAMN1
  368.       SIGF(1) = XSIGN1(1)
  369.       SIGF(2) = XSIGN1(2)
  370.       SIGF(3) = XSIGN1(3)
  371.       SIGF(4) = XSIGN1(4)
  372.       SIGF(5) = XSIGN1(5)
  373.       SIGF(6) = XSIGN1(6)
  374.  
  375.       VARF(1) = XEPPN1(1)
  376.       VARF(2) = XEPPN1(2)
  377.       VARF(3) = XEPPN1(3)
  378.       VARF(4) = XEPPN1(4)
  379.       VARF(5) = XEPPN1(5)
  380.       VARF(6) = XEPPN1(6)
  381.       VARF(7) = XCHIF
  382.       VARF(8) = XGAMF
  383.  
  384.       RETURN
  385.  
  386.       END
  387.  
  388.  
  389.  

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