Télécharger la plaquette Cast3M


fr en

Télécharger acti1.eso

Retour à la liste

Numérotation des lignes : 

acti1

  1. C ACTI1     SOURCE    PV        22/04/19    16:17:58     11344          
  2. C ACTI1    SOURCE
  3.       SUBROUTINE ACTI1(SIG0,SIGF,D,NSTRSS,BETINSA)
  4. C
  5. C ===============================================================
  6. C           Un seul critere de traction: DRUCKER PRAGER 3D
  7. C ================================================================
  8. C CE SOUS-PROGRAMME EST APPELE DANS "BONE3D".
  9. C
  10.       IMPLICIT INTEGER(I-N)
  11.       IMPLICIT REAL*8(A-H,O-Z)
  12.       DIMENSION SIG0(NSTRSS),SIGF(NSTRSS)
  13.       DIMENSION DFSIG(6),VEC1(6),DEPSI(6),SIGE(6)
  14.       DIMENSION D(6,6),AC(6),A(6,6),AI(6,6),AIM(4,4)
  15. C
  16.       SEGMENT BETINSA
  17.        REAL*8 RT,RC,YOUN,XNU,GFT,GFC,CAR
  18.        REAL*8 DKT,DKC,SEQT,SEQC,ENDT,ENDC
  19.        INTEGER IFIS,IPLA,IBB,IGAU
  20.       ENDSEGMENT
  21. C
  22.       CALL ZERO(SIGE,6,1)
  23.       DO 2 I=1,NSTRSS
  24.        SIGE(I)=SIGF(I)
  25.     2 CONTINUE
  26.       IAPEX=0
  27.       PRB=1.D-5
  28.       PRB2=1.D-2
  29.       ITER=1
  30.       ITR = 1500
  31.       ITANG=0
  32.       IBROY=0
  33.       CRIMAX=0.D0
  34.       SEQ = 0.D0
  35.       SEQ1 = 0.D0
  36. C
  37.     8 CONTINUE
  38. C
  39.       DO 4 I=1,NSTRSS
  40.        SIGF(I)=SIGE(I)
  41.     4 CONTINUE
  42.       DK=DKT
  43.       DLAM0=0.D0
  44.       CALL DRUTRA(SIGF,SEQTE,BETINSA)
  45.       FCRI0 = SEQTE - SEQT
  46.       CRIMAX=ABS(100.D0*FCRI0)
  47.     9 CONTINUE
  48. C
  49. C ************ CALCUL DU JACOBIEN INITIAL  ********************
  50. C
  51.       CALL DRUTR1(SIGF,DFSIG,BETINSA)
  52. C
  53.       DO    I=1,NSTRSS
  54.       AC(I)=0.D0
  55.       DO    J=1,NSTRSS
  56.       AC(I)=AC(I)+D(I,J)*DFSIG(J)
  57.       enddo
  58.       enddo
  59.       FDF=0.D0
  60.       DO 30 J=1,NSTRSS
  61.       FDF=FDF+DFSIG(J)*AC(J)
  62.    30 CONTINUE
  63. C
  64.       CALL ENDAME(1,BETINSA)
  65.       CALL FORECR(DK,PAEC,1,SEQ,BETINSA)
  66. C
  67.       DJAC0=-(PAEC+FDF)
  68. C
  69. C ************ DEBUT ITERATION INTERNES ************************
  70. C
  71.    40 CONTINUE
  72. C
  73. C *************** Determination de DK  *************************
  74. C
  75.       DLAM1=-FCRI0/DJAC0+DLAM0
  76.       DK=DKT+DLAM1
  77. C
  78. C *********** Estimation contrainte quivalente *****************
  79. C
  80.       CALL ENDAME(1,BETINSA)
  81.       CALL FORECR(DK,PAEC,1,SEQ,BETINSA)
  82. C
  83. C ************** Determination de DPHI *************************
  84. C
  85.       CALL DRUTR2(SIGE,SEQ,DPHI,DLAM1,VEC1,BETINSA)
  86. C
  87. C **************** Cas de l'apex *******************************
  88. C
  89.       IF (ABS(DPHI).LE.10E-10) THEN
  90.       IAPEX=1
  91.       DO    I=1,NSTRSS
  92.       DO    J=1,NSTRSS
  93.       AI(I,J)=0.D0
  94.       enddo
  95.       enddo
  96.         AI(1,1)=1./3.
  97.         AI(1,2)=AI(1,1)
  98.         AI(1,2)=AI(1,1)
  99.         AI(1,3)=AI(1,1)
  100.         AI(2,1)=AI(1,1)
  101.         AI(2,2)=AI(1,1)
  102.         AI(2,3)=AI(1,1)
  103.         AI(3,1)=AI(1,1)
  104.         AI(3,2)=AI(1,1)
  105.         AI(3,3)=AI(1,1)
  106.         GOTO 75
  107.       ENDIF
  108. C
  109. C ************** Mise a jour des contraintes ***************
  110. C
  111. C ---------------- calcul de la matrice A ------------------
  112. C
  113.       DO    I=1,NSTRSS
  114.       DO    J=1,NSTRSS
  115.       A(I,J)=0.D0
  116.       enddo
  117.       enddo
  118. C
  119.       DG=YOUN/(1.D0+XNU)
  120. C
  121.       A(1,1)=1.D0+2.*(DLAM1*DG)/2.D0/DPHI
  122.       A(2,2)=A(1,1)
  123.       A(3,3)=A(1,1)
  124.       A(1,2)=-(DLAM1*DG)/2.D0/DPHI
  125.       A(1,3)=A(1,2)
  126.       A(2,1)=A(1,2)
  127.       A(2,3)=A(1,2)
  128.       A(3,1)=A(1,2)
  129.       A(3,2)=A(1,2)
  130.       A(4,4)=1.D0+3.*(DLAM1*DG)/2.D0/DPHI
  131.       A(5,5)=A(4,4)
  132.       A(6,6)=A(4,4)
  133. C
  134. C -------------- invertion de la matrice A -----------------
  135. C
  136.       CALL ZERO(AI,6,6)
  137.       CALL ZERO(AIM,4,4)
  138. C
  139.       DO    I=1,3
  140.       DO    J=1,3
  141.       AIM(I,J)=A(I,J)
  142.       enddo
  143.       enddo
  144.       CALL INVMA2(AIM,3,ISING)
  145.       IF (ISING.EQ.1) THEN
  146.        WRITE(*,*)'MATRICE AIM singuliere ds ACTI1'
  147.       ENDIF
  148.       DO    I=1,3
  149.       DO    J=1,3
  150.       AI(I,J)=AIM(I,J)
  151.       enddo
  152.       enddo
  153.       AI(4,4) = 1./A(4,4)
  154.       AI(5,5) = 1./A(5,5)
  155.       AI(6,6) = 1./A(6,6)
  156. C
  157. C -------------- mise a jour des contraintes  ------------
  158. C
  159.    75 CONTINUE
  160. C
  161.       DO 80 I=1,NSTRSS
  162.       DEPSI(I)=SIGE(I)-DLAM1*VEC1(I)
  163.    80 CONTINUE
  164. C
  165.       DO    I=1,NSTRSS
  166.       SIGF(I)=0.
  167.       DO    J=1,NSTRSS
  168.       SIGF(I)=SIGF(I)+AI(I,J)*DEPSI(J)
  169.       enddo
  170.       enddo
  171. C
  172. C ******** VERIFICATION DU CRITERE ****************
  173. C
  174.       CALL DRUTRA(SIGF,SEQI,BETINSA)
  175.       FCRI1 = SEQI - SEQ
  176. C
  177.       IF (IBROY.EQ.0.AND.ABS(FCRI1).GE.CRIMAX) THEN
  178. C      WRITE(*,*)'****************************************'
  179. C      WRITE(*,*)'LE RESIDU DIVERGE AVEC BROYDEN'
  180. C      WRITE(*,*)'on passe donc a la secante'
  181. C      WRITE(*,*)'Dans l element',IBB
  182. C      WRITE(*,*)'et au point d intégration',IGAU
  183. C      WRITE(*,*)'CRIMAX=',CRIMAX
  184. C      WRITE(*,*)'****************************************'
  185.       ITER=ITR
  186.       ENDIF
  187. C
  188. C *******  Compteur sur la methode de resolution ****
  189. C
  190.       IF (IBROY.EQ.0.AND.ITER.EQ.ITR) THEN
  191.         IBROY=1
  192.         ITANG=1
  193.         ITER=1
  194.         IAPEX=0
  195.         GOTO 8
  196.       ENDIF
  197. C
  198. C ******* non convergence **************************
  199. C
  200.       IF (ABS(FCRI1).GT.PRB.AND.ITER.LT.ITR) THEN
  201.        IF (IBROY.EQ.0) THEN
  202.        DJAC1=(FCRI0-FCRI1)/(DLAM0-DLAM1)
  203.        DLAM0=DLAM1
  204.        DJAC0=DJAC1
  205.        FCRI0=FCRI1
  206.        ITER=ITER+1
  207.        GOTO 40
  208.        ENDIF
  209.        IF (IBROY.EQ.1.AND.ITANG.EQ.1) THEN
  210.        DLAM0=DLAM1
  211.        FCRI0=FCRI1
  212.        ITER=ITER+1
  213.         IF (ITER.GE.(ITR-5)) THEN
  214.         WRITE(*,*)'ITER=',ITER
  215.         WRITE(*,*)'FCRI=',FCRI1
  216.         ENDIF
  217.        GOTO 9
  218.        ENDIF
  219.       ENDIF
  220.       IF (ITER.GE.ITR.AND.ABS(FCRI1).GT.PRB2) THEN
  221.         WRITE(*,*)'NON CONVERGENCE INTERNE dans ACTI1'
  222.         WRITE(*,*)'FCRI=',FCRI1
  223.         WRITE(*,*)'Dans l element',IBB
  224.         WRITE(*,*)'et au point d intégration',IGAU
  225.         STOP
  226.       ENDIF
  227. C
  228. C *******  convergence **************************
  229. C
  230.       SEQT = SEQ
  231.       DKT = DK
  232. C
  233. C ***********************************************
  234. C
  235.   100 CONTINUE
  236. C
  237.       RETURN
  238.       END
  239.  
  240.  
  241.  
  242.  
  243.  
  244.  
  245.  
  246.  
  247.  
  248.  

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