Télécharger la plaquette Cast3M


fr en

Télécharger cmicro.eso

Retour à la liste

Numérotation des lignes : 

cmicro

  1. C CMICRO    SOURCE    OF166741  25/11/04    21:15:29     12349          
  2.  
  3.       SUBROUTINE CMICRO (WRK52,WRK53,WRK54,NVARI,iecou)
  4. *
  5. *   modele d'endommagement microplan couple a la plasticite
  6. *   C. La Borderie + S. Fichant Oct. 95
  7. *   routines utilisees:
  8. *   micro1: plasticite nadai
  9. *            IDECAL=3 DANS LE CAS ISO IDECAL=8 DANS LE CAS ANISO
  10. *   jacob3: diagonalisation:
  11. *            attention jacob3 modifie la matrice a diagonaliser!!
  12. *   prodt et prodt2
  13. *   attention prodt2 ne fonctionne qu'avec la matrice des V. P. !!
  14. *
  15.       IMPLICIT INTEGER(I-N)
  16.       IMPLICIT REAL*8(A-H,O-Z)
  17.  
  18.       INTEGER NVARI
  19.       INTEGER ISTRS,I,J
  20.       REAL*8 YOUNG,XNU,EPSD0,BT,DEUXMU,ALFA
  21.       REAL*8 DOM(6),SIGAN(6),TRSIG,DEF33(3,3),EPSIPP(3),VECP(3,3)
  22.       REAL*8 DEFPT(6),DEFTOT(6),EPSITOP(3),EPSE1,EPSE
  23.       REAL*8 VECPT(3,3),DOM33(3,3),DEFRPDE(3,3),DOMRPDE(3,3)
  24.       REAL*8 D1,D2,D3,DOM3(3)
  25.       REAL*8 SIGPP(3),SIGPM(3),SIG33(3,3),SIG33P(3,3),SIG33M(3,3)
  26.       REAL*8 S33PRD(3,3),S33MRD(3,3)
  27.       REAL*8 LAMBDAP(6), LAMBDAM(6)
  28.       LOGICAL COMP
  29.       INTEGER IDECAL
  30.  
  31. -INC PPARAM
  32. -INC CCOPTIO
  33. -INC DECHE
  34. -INC TECOU
  35.  
  36.       SEGMENT WRKK1
  37.         REAL*8 DEFELA(NSTRS1)
  38.         REAL*8 DDEFEL(NSTRS1)
  39.       ENDSEGMENT
  40.  
  41. *      print*,'dans microp'
  42. *
  43. *   on recupere les variable materielles
  44. *
  45.       YOUNG=XMAT(1)
  46.       XNU=XMAT(2)
  47.       EPSD0=XMAT(5)
  48.       BT=XMAT(6)
  49.       ALFA=XMAT(9)
  50.       DEUXMU=YOUNG/(1.D0+XNU)
  51.  
  52.       NSTRS1 = iecou.NSTRSS
  53. *
  54. *     recuperation des variables internes d'endommagement
  55. *
  56.       DO ISTRS=1,NSTRS1
  57.         DOM(ISTRS)=VAR0(2+ISTRS)
  58.       END DO
  59. *
  60. *     on ecoule plastiquement sur la contrainte effective
  61. *
  62.       IDECAL=8
  63. *      print*,'appel a micro1'
  64.       CALL CMICR1(wrk52,wrk53,wrk54,nstrs1,NVARI,IDECAL,
  65.      &      .false.,DEFPT,EPSE1,EPSE)
  66.       IF (KERRE .NE. 0) THEN
  67.            print*,'on n''a pas converge dans micro1'
  68.  
  69.         CALL CMICR1(wrk52,wrk53,wrk54,nstrs1,NVARI,IDECAL,
  70.      &              .true.,DEFPT,EPSE1,EPSE)
  71.          RETURN
  72.       ENDIF
  73. *
  74. *     on ecoule en endommagement sur les deformations elastiques
  75. *
  76. *      print*,'apres micro1'
  77.  
  78.       SEGINI WRKK1
  79. *
  80. *     calcul de l'increment deformations elastiques DDEFEL
  81. *     1) on calcule l'increment de deformations totales avec
  82. *        l'increment de contraintes elastique DSIGT
  83. *     2) on retranche l'increment de deformations plastiques DEFP
  84. *
  85.         TRSIG= DSIGT(1)+DSIGT(2)+DSIGT(3)
  86.         DO ISTRS=1,3
  87.           DDEFEL(ISTRS)=( (1.D0+XNU)*DSIGT(ISTRS)-XNU*TRSIG)/YOUNG
  88.      1                 - DEFP(ISTRS)
  89.         END DO
  90.         DO ISTRS=4,NSTRS1
  91.           DDEFEL(ISTRS)= (1.D0+XNU)*DSIGT(ISTRS)/YOUNG
  92.      1                 - 0.5d0*DEFP(ISTRS)
  93.         END DO
  94. *
  95. *      on diagonalise l'increment de deformations elastiques:
  96. *      1) on met sous forme 3x3 avec endoca
  97. *      2) on diagonalise avec jacob3
  98. *
  99. *       print*,'incr de def el dans rpg'
  100. *       print*,ddefel
  101. *       print*,'avant endoca ddefl'
  102.          CALL ENDOCA (DDEFEL,DEF33,1)
  103. *       print*,'apres endoca ddefl'
  104.          CALL JACOB3 (DEF33,IDIM,EPSIPP,VECP)
  105. *       print*,'apres jacob3 def33'
  106. *       print*,'incr deformations principales'
  107. *       print*,epsipp
  108. *
  109. *     calcul des deformations elastiques
  110. *
  111. *      print*,'-----sigf------'
  112. *      print*,sigf(1),sigf(2),sigf(3)
  113.          TRSIG=SIGF(1)+SIGF(2)+SIGF(3)
  114.          DO ISTRS=1,3
  115.            DEFELA(ISTRS)=( (1.D0+XNU)*SIGF(ISTRS)-XNU*TRSIG)/YOUNG
  116.          END DO
  117. *      print*,'-----defela------'
  118. *      print*,defela(1),defela(2),defela(3)
  119.          DO ISTRS=4,NSTRS1
  120.            DEFELA(ISTRS)= (1.D0+XNU)*SIGF(ISTRS)/YOUNG
  121.          END DO
  122. *         print*,'deformations elastiques dans rpg'
  123. *         print*,defela
  124. *
  125. *     on met les deformations sous forme de matrice 3x3
  126. *     puis on ecrit la matrice dans le repere de depsilon:DEFRPDE
  127. *
  128.          CALL ENDOCA (DEFELA,DEF33,1)
  129. *         print*,'deformations elastiques dans rpg(3x3)'
  130. *         print*,def33
  131.          CALL PRODT (DEFRPDE,DEF33,VECP,3,3)
  132. *          print*,'def elast dans rpddeps'
  133. *         print*,defrpde
  134. *
  135. *     on calcule l'endommagement resultant
  136. *
  137.          IF (DEFRPDE(1,1) .GT. EPSD0) THEN
  138.            D1=1.D0-EPSD0/DEFRPDE(1,1)*EXP(BT*( EPSD0 - DEFRPDE(1,1)))
  139.          ELSE
  140.            D1=0.D0
  141.          END IF
  142.          IF (DEFRPDE(2,2) .GT. EPSD0) THEN
  143.            D2=1.D0-EPSD0/DEFRPDE(2,2)*EXP(BT*( EPSD0 - DEFRPDE(2,2)))
  144.          ELSE
  145.            D2=0.D0
  146.          END IF
  147.          IF (DEFRPDE(3,3) .GT. EPSD0) THEN
  148.            D3=1.D0-EPSD0/DEFRPDE(3,3)*EXP(BT*( EPSD0 - DEFRPDE(3,3)))
  149.          ELSE
  150.            D3=0.D0
  151.          END IF
  152. *      print*,defrpde(1,1),epsd0,D1
  153. *      print*,defrpde(2,2),epsd0,D2
  154. *      print*,defrpde(3,3),epsd0,D3
  155. *
  156. *     on met l'endommagement initial dans le meme repere
  157. *
  158.          CALL ENDOCA(DOM,DOM33,1)
  159. *      print*,'prodt DOMRPDE'
  160.          CALL PRODT (DOMRPDE,DOM33,VECP,3,3)
  161. *       print*,'endommagement initial dans rpddeps'
  162. *       print*,domrpde
  163. *
  164. *     et on en deduit l'increment d'endommagement dans RPDE
  165. *
  166.          IF(d1.gt.domrpde(1,1))then
  167.             domrpde(1,1)=d1
  168.          endif
  169.          IF(d2.gt.domrpde(2,2))then
  170.             domrpde(2,2)=d2
  171.          endif
  172.          IF(d3.gt.domrpde(3,3))then
  173.             domrpde(3,3)=d3
  174.          endif
  175. *       print*,'endommagement final dans rpddeps'
  176. *       print*,domrpde
  177. *
  178. *     on remet D dans RPG
  179. *
  180.       DO J=1,3
  181.         DO I=1,3
  182.            VECPT(I,J)=VECP(J,I)
  183.         END DO
  184.       END DO
  185.       call prodt(DOM33,DOMRPDE,VECPT,3,3)
  186. *       print*,'endommagement final dans rpg'
  187. *       print*,dom33
  188. *
  189. *     on separe les contraintes effectives en + et - dans rpsigma
  190. *
  191.          CALL ENDOCA (SIGF,SIG33,1)
  192.          CALL JACOB3 (SIG33,3,SIGPP,VECP)
  193. *         print*,'contraintes ppales'
  194. *         print*,sigpp
  195.          DO I=1,3
  196.            IF (SIGPP(I) .LT. 0.D0)THEN
  197.               SIGPM(I)=SIGPP(I)
  198.               SIGPP(I)=0.D0
  199.            ELSE
  200.               SIGPM(I)=0.D0
  201.            END IF
  202.          END DO
  203.          CALL PRODT2(SIG33P,SIGPP,VECP,3)
  204.          CALL PRODT2(SIG33M,SIGPM,VECP,3)
  205. *         print*,'contraintes dans rpg'
  206. *         print*,sig33p
  207. *         print*,sig33m
  208. *
  209. *      on met le tout dans le repere ppal d'endo
  210. *      attention jacob3 modifie la matrice fournie
  211. *       --> on passe une copie
  212.          DO J=1,3
  213.            DO I=1,3
  214.              SIG33(I,J)=DOM33(I,J)
  215.            END DO
  216.          END DO
  217.          CALL JACOB3(SIG33,3,DOM3,VECP)
  218. *         print*,'endom dans rpd'
  219. *         print*,dom3
  220. **************
  221. *  A REVOIR  EN ATENDANT MIEUX ON BORNE LES VALEURS PROPRES DE D ENTRE 0. ET 1.*
  222. **************
  223.          DO I=1,3
  224.            DOM3(I)=MAX(DOM3(I),0.D0)
  225.            DOM3(I)=MIN(DOM3(I),1.D0-1.d-6)
  226.          END DO
  227.          CALL PRODT2 (DOM33,DOM3,VECP,3)
  228. **************
  229. *  FIN A REVOIR
  230. **************
  231. *      print*,'dom3=',dom3
  232. *      print*,'dom33=',dom33
  233.          CALL PRODT (S33PRD,SIG33P,VECP,3,3)
  234.          CALL PRODT (S33MRD,SIG33M,VECP,3,3)
  235.  
  236. *      print*,'s33prd',s33prd(3,3)
  237. *      print*,'s33mrd',s33mrd(3,3)
  238. *      print*,'callambdap'
  239.          COMP=.FALSE.
  240.          CALL CLMBDA(DOM3(1),DOM3(2),DOM3(3),LAMBDAP,YOUNG,XNU,ALFA,
  241.      1                                                            COMP)
  242.          COMP=.TRUE.
  243.          CALL CLMBDA(DOM3(1),DOM3(2),DOM3(3),LAMBDAM,YOUNG,XNU,ALFA,
  244.      1                                                            COMP)
  245. *       print*,'apres callambdap'
  246.        SIG33P(1,1)=LAMBDAP(1)*S33PRD(1,1)+LAMBDAM(1)*S33MRD(1,1)
  247.        SIG33P(1,2)=LAMBDAP(4)*S33PRD(1,2)+LAMBDAM(4)*S33MRD(1,2)
  248.        SIG33P(1,3)=LAMBDAP(5)*S33PRD(1,3)+LAMBDAM(5)*S33MRD(1,3)
  249.        SIG33P(2,1)=SIG33P(1,2)
  250.        SIG33P(2,2)=LAMBDAP(2)*S33PRD(2,2)+LAMBDAM(2)*S33MRD(2,2)
  251.        SIG33P(2,3)=LAMBDAP(6)*S33PRD(2,3)+LAMBDAM(6)*S33MRD(2,3)
  252.        SIG33P(3,1)=SIG33P(1,3)
  253.        SIG33P(3,2)=SIG33P(2,3)
  254.        SIG33P(3,3)=LAMBDAP(3)*S33PRD(3,3)+LAMBDAM(3)*S33MRD(3,3)
  255. *       print*,'contraintes dans RPD'
  256. *       print*,sig33p
  257. *       print*,'sig33p(1,2)',sig33p(1,2)
  258. *      print*,'dsigt',dsigt(1),dsigt(2),dsigt(3)
  259. *      print*,'dsigt',dsigt(4),dsigt(5),dsigt(6)
  260. *      print*,'sigf',sigf(1),sigf(2),sigf(3)
  261. *      print*,'sigf',sigf(4),sigf(5),sigf(6)
  262. *
  263. *     on remet le tout dans le repere global
  264. *
  265. *      print*,'matrice des vect proprs'
  266. *      print*,vecp
  267.       DO J=1,3
  268.         DO I=1,3
  269.            VECPT(I,J)=VECP(J,I)
  270.         END DO
  271.       END DO
  272.       CALL PRODT (SIG33,SIG33P,VECPT,3,3)
  273. *      print*,'contraintes dans rpg'
  274. *      print*,sig33
  275. *
  276. *     on rend les contraintes et les variables internes finales
  277. *
  278.       SIGAN(1)=SIGF(1)-SIG33(1,1)
  279.       SIGF(1)=SIG33(1,1)
  280.       VARF(3)=MAX(DOM33(1,1),0.d0)
  281.       SIGAN(2)=SIGF(2)-SIG33(2,2)
  282.       SIGF(2)=SIG33(2,2)
  283.       VARF(4)=MAX(DOM33(2,2),0.d0)
  284.       SIGAN(3)=SIGF(3)-SIG33(3,3)
  285.       SIGF(3)=SIG33(3,3)
  286.       VARF(5)=MAX(DOM33(3,3),0.d0)
  287.       SIGAN(4)=SIGF(4)-SIG33(1,2)
  288.       SIGF(4)=SIG33(1,2)
  289.       VARF(6)=MAX(DOM33(1,2),0.d0)
  290.       VARF(7)=MAX(DOM33(1,3),0.d0)
  291.       VARF(8)=MAX(DOM33(2,3),0.d0)
  292.       IF(IFOUR.GE.1.OR.IFOUR.LE.-3) THEN
  293.         SIGAN(5)=SIGF(5)-SIG33(1,3)
  294.         SIGF(5)=SIG33(1,3)
  295.         SIGAN(6)=SIGF(6)-SIG33(2,3)
  296.         SIGF(6)=SIG33(2,3)
  297.       ELSE
  298.         SIGAN(5)=0.D0
  299.         SIGAN(6)=0.D0
  300.       END IF
  301.       DO ISTRS=1,6
  302.          VARF(ISTRS+8)=SIGAN(ISTRS)
  303.       END DO
  304. *      print*,'sigf',sigf(1),sigf(2),sigf(3)
  305. *      print*,'sigf',sigf(4),sigf(5),sigf(6)
  306. *      print*,'sigan',sigan(1),sigan(2),sigan(3)
  307. *      print*,'sigan',sigan(4),sigan(5),sigan(6)
  308.       SEGSUP WRKK1
  309.  
  310.       RETURN
  311.       END
  312.  
  313.  
  314.  

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