Télécharger la plaquette Cast3M


fr en

Télécharger incre7.eso

Retour à la liste

Numérotation des lignes : 

incre7

  1. C INCRE7    SOURCE    CHAT      05/01/13    00:35:30     5004
  2.        SUBROUTINE INCRE7(SIG,VAR,DSPT,EPSVPT,VARPT,XMAT,NSTRS0,
  3.      &     MFR,NVARI,NCOMAT,IFOUR)
  4.       IMPLICIT INTEGER(I-N)
  5.        IMPLICIT REAL*8(A-H,O-Z)
  6.        DIMENSION SIG(*),VAR(*),DSPT(*)
  7.        DIMENSION EPSVPT(*),VARPT(*),XMAT(*)
  8.        DIMENSION XX(6),YY(6),EVPT(6),EPPT(6),SIGPT(6)
  9.        DIMENSION SIG0(6),EPS0(6)
  10.        PARAMETER (AMAX = 1.0D20 , AMIN = 1.D-10)
  11.        DETIER = 2.0D0/3.0D0
  12.        ROOT = SQRT(3.0D0/2.0D0)
  13.        CALL ZDANUL(VARPT,NVARI)
  14.        IF (MFR.EQ.5) THEN
  15.          NSTRS=6
  16.          SIG0(1)=SIG(1)
  17.          SIG0(2)=SIG(2)
  18.          SIG0(3)=0.D0
  19.          SIG0(4)=SIG(3)
  20.          SIG0(5)=SIG(4)
  21.          SIG0(6)=SIG(5)
  22.        ELSE
  23.          NSTRS=NSTRS0
  24.          DO 10 I=1,NSTRS
  25.           SIG0(I)=SIG(I)
  26.   10     CONTINUE
  27.        ENDIF
  28.        TRACE =(SIG0(1)+SIG0(2)+SIG0(3))/3.0D0
  29. C-------------------------------------------------------------------|
  30. C                      COEFFICIENTS MATERIAU                        |
  31. C-------------------------------------------------------------------|
  32.        YOU =XMAT( 1)
  33.        XNU =XMAT( 2)
  34.        RP0 =XMAT( 5)
  35.        QP  =XMAT( 6)
  36.        BP  =XMAT( 7)
  37.        CP1 =XMAT( 8)
  38.        CP2 =XMAT( 9)
  39.        DP1 =XMAT(10)
  40.        DP2 =XMAT(11)
  41.        XK  =XMAT(12)
  42.        XN  =XMAT(13)
  43.        RV0 =XMAT(14)
  44.        QV  =XMAT(15)
  45.        BV  =XMAT(16)
  46.        CV1 =XMAT(17)
  47.        CV2 =XMAT(18)
  48.        DV1 =XMAT(19)
  49.        DV2 =XMAT(20)
  50.        CVP1=XMAT(21)
  51.        CVP2=XMAT(22)
  52.  
  53.        IF (CP1.EQ.0.D0) THEN
  54.          DCP1 = 0.D0
  55.        ELSE
  56.          DCP1 = DP1 / CP1
  57.        ENDIF
  58.  
  59.        IF (CP2.EQ.0.D0) THEN
  60.          DCP2 = 0.D0
  61.        ELSE
  62.          DCP2 = DP2 / CP2
  63.        ENDIF
  64.  
  65.        IF (CV1.EQ.0.D0) THEN
  66.          DCV1 = 0.D0
  67.        ELSE
  68.          DCV1 = DV1 / CV1
  69.        ENDIF
  70.  
  71.        IF (CV2.EQ.0.D0) THEN
  72.          DCV2 = 0.D0
  73.        ELSE
  74.          DCV2 = DV2 / CV2
  75.        ENDIF
  76.  
  77. C-------------------------------------------------------------------|
  78. C                        PARTIE VISCOPLASTIQUE                      |
  79. C-------------------------------------------------------------------I
  80. C
  81. C *** CALCUL DE XX=DEVIATEUR(SIGMA-XP)
  82. C
  83.        TRACE =(SIG0(1)+SIG0(2)+SIG0(3))/3.0D0
  84.        DO 71 I=1,NSTRS
  85.        A = 0.0D0
  86.        IF (I.LE.3) A=1.0D0
  87.        XX(I) = SIG0(I)-A*TRACE - VAR(5*NSTRS+I)-VAR(6*NSTRS+I)
  88.   71   CONTINUE
  89. C
  90. C *** CRITERE VISCOPLASTIQUE
  91. C
  92.        AJ2 = PROCON(XX,XX,NSTRS)
  93.        AJ2 = SQRT(1.5D0*AJ2)
  94.        CRIT= AJ2-VAR(10*NSTRS+4)-RV0
  95. C
  96. C *** ECOULEMENT
  97. C
  98.        IF(CRIT.GT.0.D0) THEN
  99.          VPT=(CRIT/XK)**XN
  100.        ELSE
  101.          VPT=0.D0
  102.        ENDIF
  103.        VARPT(10*NSTRS+3)=VPT
  104.        DO 80 I=1,NSTRS
  105.        IF(VPT.EQ.0.D0) THEN
  106.          EVPT(I)=0.D0
  107.          VARPT(7*NSTRS+I)=0.D0
  108.          VARPT(8*NSTRS+I)=0.D0
  109.          VARPT(9*NSTRS+I)=0.D0
  110.        ELSE
  111.          EVPT(I)=1.5*XX(I)*VPT/AJ2
  112.          VARPT(9*NSTRS+I)=EVPT(I)
  113.          VARPT(7*NSTRS+I)=EVPT(I)-(1.5*DCV1*VAR(5*NSTRS+I)*VPT)
  114.          VARPT(8*NSTRS+I)=EVPT(I)-(1.5*DCV2*VAR(6*NSTRS+I)*VPT)
  115.        ENDIF
  116.   80   CONTINUE
  117. C-------------------------------------------------------------------|
  118. C                         PARTIE PLASTIQUE                          |
  119. C-------------------------------------------------------------------I
  120. C
  121. C *** CALCUL DE XX=DEVIATEUR(SIGMA-XP)/J2(SIGMA-X)
  122. C
  123.        DO 70 I=1,NSTRS
  124.        A = 0.0D0
  125.        IF (I.LE.3) A=1.0D0
  126.        XX(I) = SIG0(I)-A*TRACE - VAR(I)-VAR(NSTRS+I)
  127.   70   CONTINUE
  128.        AJ2 = PROCON (XX,XX,NSTRS)
  129.        AJ2 = SQRT(1.5D0*AJ2)
  130.        IF(AJ2.NE.0.D0) THEN
  131.          DO 81 I=1,NSTRS
  132.          XX(I)=1.5*XX(I)/AJ2
  133.   81     CONTINUE
  134.        ENDIF
  135. C
  136. C *** CRITERE PLASTIQUE
  137. C
  138.        CRIT=AJ2-VAR(10*NSTRS+2)-RP0
  139. C
  140. C *** CALCUL DE PPT1=XX:(E:(EPSPT-EVPT)-2/3*SIGVPPT)
  141. C
  142.        IF(CRIT.LE.0.D0) THEN
  143.          PPT=0.D0
  144.          VARPT(10*NSTRS+1)=PPT
  145.          DO 85 I=1,NSTRS
  146.          VARPT(2*NSTRS+I)=0.D0
  147.          VARPT(3*NSTRS+I)=0.D0
  148.          VARPT(4*NSTRS+I)=0.D0
  149.          EPPT(I)=0.D0
  150.    85    CONTINUE
  151.        ELSE
  152.        IF(IFOUR.EQ.-2) THEN
  153.        X2MU=YOU/(1.0+XNU)
  154.        XCO = X2MU/(1.0-XNU)
  155.        SIGPT(1)=DSPT(1)-XCO*(EVPT(1)+XNU*EVPT(2))
  156.        SIGPT(2)=DSPT(2)-XCO*(EVPT(2)+XNU*EVPT(1))
  157.        SIGPT(4)=DSPT(4)-X2MU*EVPT(4)
  158.        SIGPT(3)=0.0D0
  159.        ELSEIF (IFOUR.EQ.2.OR.IFOUR.EQ.0.OR.IFOUR.EQ.1.OR.
  160.      &         IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
  161.        X2MU = YOU/(1+XNU)
  162.        ALMB = X2MU*XNU/(1-2*XNU)
  163.        IF (MFR.EQ.5) THEN
  164.        X2MU=YOU/(1.0+XNU)
  165.        XCO = X2MU/(1.0-XNU)
  166.        SIGPT(1)=DSPT(1)-XCO*(EVPT(1)+XNU*EVPT(2))
  167.        SIGPT(2)=DSPT(2)-XCO*(EVPT(2)+XNU*EVPT(1))
  168.        SIGPT(3)=DSPT(3)-X2MU*EVPT(3)
  169.        SIGPT(4)=DSPT(4)-X2MU*EVPT(4)
  170.        SIGPT(5)=DSPT(5)-X2MU*EVPT(5)
  171.        ELSE
  172.         TRACE=EVPT(1)+EVPT(2)+EVPT(3)
  173.         DO 20 I=1,NSTRS,1
  174.         A = 0.0
  175.         IF (I.LE.3) A=1.0
  176.         SIGPT(I) = DSPT(I)-ALMB*A*TRACE-X2MU*EVPT(I)
  177.    20   CONTINUE
  178.        ENDIF
  179. C= Modes de calcul 1D contraintes planes suivant z (DYCZ et GYCZ)
  180.        ELSE IF (IFOUR.EQ.4.OR.IFOUR.EQ.8) THEN
  181.          XCO = YOU/(1.-XNU*XNU)
  182.          SIGPT(1) = DSPT(1)-XCO*(EVPT(1)+XNU*EVPT(2))
  183.          SIGPT(2) = DSPT(2)-XCO*(EVPT(2)+XNU*EVPT(1))
  184.          SIGPT(3) = 0.D0
  185. C= Modes de calcul 1D contraintes planes suivant y (CYDZ et CYGZ) et
  186. C=                 1D axisymetrique contraintes axiales nulles (AXCZ)
  187.        ELSE IF (IFOUR.EQ.5.OR.IFOUR.EQ.10.OR.IFOUR.EQ.13) THEN
  188.          XCO = YOU/(1.-XNU*XNU)
  189.          SIGPT(1) = DSPT(1)-XCO*(EVPT(1)+XNU*EVPT(3))
  190.          SIGPT(2) = 0.D0
  191.          SIGPT(3) = DSPT(3)-XCO*(EVPT(3)+XNU*EVPT(1))
  192. C= Modes de calcul 1D contraintes planes suivant y et z (CYCZ)
  193.        ELSE IF (IFOUR.EQ.6) THEN
  194.          SIGPT(1) = DSPT(1)-YOU*EVPT(1)
  195.          SIGPT(2) = 0.D0
  196.          SIGPT(3) = 0.D0
  197. C= Autres modes de calcul 1D deformations planes (DYDZ GYDZ DYGZ GYGZ)
  198. C=                        1D axisymetrique (AXDZ AXGZ) et 1D spherique
  199.        ELSE IF (IFOUR.EQ. 3.OR.IFOUR.EQ. 7.OR.IFOUR.EQ. 9.OR.IFOUR.EQ.11
  200.      .      .OR.IFOUR.EQ.12.OR.IFOUR.EQ.14.OR.IFOUR.EQ.15) THEN
  201.          X2MU = YOU/(1.+XNU)
  202.          TRACE = EVPT(1)+EVPT(2)+EVPT(3)
  203.          XCO = X2MU*XNU*TRACE/(1.-XNU)
  204.          SIGPT(1) = DSPT(1)-X2MU*EVPT(1)-XCO
  205.          SIGPT(2) = DSPT(2)-X2MU*EVPT(2)-XCO
  206.          SIGPT(3) = DSPT(3)-X2MU*EVPT(3)-XCO
  207.        ENDIF
  208.        DO 30 I=1,NSTRS
  209.        SIGPT(I)=SIGPT(I)-DETIER*(
  210.      &         CVP1*VARPT(7*NSTRS+I)-CVP2*VARPT(8*NSTRS+I))
  211.    30  CONTINUE
  212.        PPT1=PROCON(SIGPT,XX,NSTRS)
  213. C
  214. C *** CALCUL DE
  215. C     H1=2/3*XX:(C1P(NP-3/2*(D1P/C1P)*X1P)+C2P(NP-3/2*(D2P/C2P)*X2P))
  216. C
  217.        DO 40 I=1,NSTRS
  218.        YY(I)=CP1*(XX(I)-1.5*DCP1*VAR(I))+
  219.      &       CP2*(XX(I)-1.5*DCP2*VAR(NSTRS+I))
  220.  40    CONTINUE
  221.        H1=DETIER*PROCON(XX,YY,NSTRS)
  222. C
  223. C *** CALCUL DE H2=XX:(E:XX)
  224. C
  225.        IF(IFOUR.EQ.-2) THEN
  226.        X2MU=YOU/(1.0+XNU)
  227.        XCO = X2MU/(1.0-XNU)
  228.        YY(1)=XCO*(XX(1)+XNU*XX(2))
  229.        YY(2)=XCO*(XX(2)+XNU*XX(1))
  230.        YY(4)=X2MU*XX(4)
  231.        YY(3)=0.0D0
  232.        ELSEIF (IFOUR.EQ.2.OR.IFOUR.EQ.0.OR.IFOUR.EQ.1.OR.
  233.      &         IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
  234.        X2MU = YOU/(1+XNU)
  235.        ALMB = X2MU*XNU/(1-2*XNU)
  236.        IF (MFR.EQ.5) THEN
  237.        X2MU=YOU/(1.0+XNU)
  238.        XCO = X2MU/(1.0-XNU)
  239.        YY(1)=XCO*(XX(1)+XNU*XX(2))
  240.        YY(2)=XCO*(XX(2)+XNU*XX(1))
  241.        YY(3)=X2MU*XX(3)
  242.        YY(4)=X2MU*XX(4)
  243.        YY(5)=X2MU*XX(5)
  244.        ELSE
  245.         TRACE=XX(1)+XX(2)+XX(3)
  246.         DO 50 I=1,NSTRS,1
  247.         A = 0.0
  248.         IF (I.LE.3) A=1.0
  249.         YY(I) = ALMB*A*TRACE+X2MU*XX(I)
  250.    50   CONTINUE
  251.        ENDIF
  252. C= Modes de calcul 1D contraintes planes suivant z (DYCZ et GYCZ)
  253.        ELSE IF (IFOUR.EQ.4.OR.IFOUR.EQ.8) THEN
  254.          XCO = YOU/(1.-XNU*XNU)
  255.          YY(1) = XCO*(XX(1)+XNU*XX(2))
  256.          YY(2) = XCO*(XX(2)+XNU*XX(1))
  257.          YY(3) = 0.D0
  258. C= Modes de calcul 1D contraintes planes suivant y (CYDZ et CYGZ) et
  259. C=                 1D axisymetrique contraintes axiales nulles (AXCZ)
  260.        ELSE IF (IFOUR.EQ.5.OR.IFOUR.EQ.10.OR.IFOUR.EQ.13) THEN
  261.          XCO = YOU/(1.-XNU*XNU)
  262.          YY(1) = XCO*(XX(1)+XNU*XX(3))
  263.          YY(2) = 0.D0
  264.          YY(3) = XCO*(XX(3)+XNU*XX(1))
  265. C= Modes de calcul 1D contraintes planes suivant y et z (CYCZ)
  266.        ELSE IF (IFOUR.EQ.6) THEN
  267.          YY(1) = YOU*XX(1)
  268.          YY(2) = 0.D0
  269.          YY(3)=0.D0
  270. C= Autres modes de calcul 1D deformations planes (DYDZ GYDZ DYGZ GYGZ)
  271. C=                        1D axisymetrique (AXDZ AXGZ) et 1D spherique
  272.        ELSE IF (IFOUR.EQ. 3.OR.IFOUR.EQ. 7.OR.IFOUR.EQ. 9.OR.IFOUR.EQ.11
  273.      .      .OR.IFOUR.EQ.12.OR.IFOUR.EQ.14.OR.IFOUR.EQ.15) THEN
  274.          X2MU = YOU/(1.+XNU)
  275.          TRACE = XX(1)+XX(2)+XX(3)
  276.          XCO = X2MU*XNU*TRACE/(1.-XNU)
  277.          YY(1) = X2MU*XX(1)+XCO
  278.          YY(2) = X2MU*XX(2)+XCO
  279.          YY(3) = X2MU*XX(3)+XCO
  280.        ENDIF
  281.        H2=PROCON(XX,YY,NSTRS)
  282. C
  283. C *** CALCUL DE H3=BP*(RP-RP0-QP)
  284. C
  285.       H3=BP*(VAR(10*NSTRS+2)-QP)
  286. C
  287. C *** MODULE PLASTIQUE
  288. C
  289.       H=H1+H2-H3
  290. C
  291. C *** CALCUL DE PPT,EPPT,ALPHAPT
  292. C
  293.       PPT=PPT1/H
  294.       IF(PPT.LT.0.D0) PPT=0.D0
  295.       VARPT(10*NSTRS+1)=PPT
  296.       DO 60 I=1,NSTRS
  297.       EPPT(I)=XX(I)*PPT
  298.       VARPT(2*NSTRS+I)=EPPT(I)-(1.5*DCP1*VAR(I)*PPT)
  299.       VARPT(3*NSTRS+I)=EPPT(I)-(1.5*DCP2*VAR(NSTRS+I)*PPT)
  300.       VARPT(4*NSTRS+I)=EPPT(I)
  301.    60 CONTINUE
  302.       ENDIF
  303. C
  304.        IF (MFR.EQ.5) THEN
  305.          EPSVPT(1)=EVPT(1)+EPPT(1)
  306.          EPSVPT(2)=EVPT(2)+EPPT(2)
  307.          EPSVPT(3)=EVPT(4)+EPPT(4)
  308.          EPSVPT(4)=EVPT(5)+EPPT(5)
  309.          EPSVPT(5)=EVPT(6)+EPPT(6)
  310.        ELSE
  311.          DO 11 I=1,NSTRS
  312.           EPSVPT(I)=EVPT(I)+EPPT(I)
  313.   11     CONTINUE
  314.        ENDIF
  315. C
  316.        VARPT(10*NSTRS+5)=VPT+PPT
  317.        RETURN
  318.        END
  319.  
  320.  
  321.  
  322.  
  323.  
  324.  
  325.  
  326.  
  327.  
  328.  
  329.  
  330.  
  331.  

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