acti2
C ACTI2 SOURCE PV 22/04/19 16:17:59 11344 C ACTI2 SOURCE C C =============================================================== C Un seul critere de compression: DRUCKER PRAGER 3D C ================================================================ C CE SOUS-PROGRAMME EST APPELE DANS "BONE3D". C IMPLICIT INTEGER(I-N) IMPLICIT REAL*8(A-H,O-Z) DIMENSION SIG0(NSTRSS),SIGF(NSTRSS) DIMENSION D(6,6),AC(6),A(6,6),AI(6,6),AIM(4,4) C SEGMENT BETINSA REAL*8 RT,RC,YOUN,XNU,GFT,GFC,CAR REAL*8 DKT,DKC,SEQT,SEQC,ENDT,ENDC INTEGER IFIS,IPLA,IBB,IGAU ENDSEGMENT C DO 2 I=1,NSTRSS SIGE(I)=SIGF(I) 2 CONTINUE IAPEX=0 PRB=1.D-5 PRB2=1.D-2 ITER=1 ITR = 1500 ITANG=0 IBROY=0 CRIMAX=0.D0 SEQ = 0.D0 SEQ1 = 0.D0 C 8 CONTINUE C DO 4 I=1,NSTRSS SIGF(I)=SIGE(I) 4 CONTINUE DK=DKC DLAM0=0.D0 FCRI0 = SEQCE - SEQC CRIMAX=ABS(100.D0*FCRI0) 9 CONTINUE C C ************ CALCUL DU JACOBIEN INITIAL ******************** C C DO I=1,NSTRSS AC(I)=0.D0 DO J=1,NSTRSS enddo enddo FDF=0.D0 DO 30 J=1,NSTRSS 30 CONTINUE C C DJAC0=-(PAEC+FDF) C C ************ DEBUT ITERATION INTERNES ************************ C 40 CONTINUE C C *************** Determination de DK ************************* C DLAM1=-FCRI0/DJAC0+DLAM0 DK=DKC+DLAM1 C C *********** Estimation contrainte quivalente ***************** C C C ************** Determination de DPHI ************************* C C C **************** Cas de l'apex ******************************* C IF (ABS(DPHI).LE.10E-10) THEN IAPEX=1 DO I=1,NSTRSS DO J=1,NSTRSS AI(I,J)=0.D0 enddo enddo AI(1,1)=1./3. AI(1,2)=AI(1,1) AI(1,2)=AI(1,1) AI(1,3)=AI(1,1) AI(2,1)=AI(1,1) AI(2,2)=AI(1,1) AI(2,3)=AI(1,1) AI(3,1)=AI(1,1) AI(3,2)=AI(1,1) AI(3,3)=AI(1,1) GOTO 75 ENDIF C C ************** Mise a jour des contraintes *************** C C ---------------- calcul de la matrice A ------------------ C DO I=1,NSTRSS DO J=1,NSTRSS A(I,J)=0.D0 enddo enddo C DG=YOUN/(1.D0+XNU) C A(1,1)=1.D0+2.*(DLAM1*DG)/2.D0/DPHI A(2,2)=A(1,1) A(3,3)=A(1,1) A(1,2)=-(DLAM1*DG)/2.D0/DPHI A(1,3)=A(1,2) A(2,1)=A(1,2) A(2,3)=A(1,2) A(3,1)=A(1,2) A(3,2)=A(1,2) A(4,4)=(1.D0)+(3.*(DLAM1*DG)/2.D0/DPHI) A(5,5)=A(4,4) A(6,6)=A(4,4) C C -------------- invertion de la matrice A ----------------- C C DO I=1,3 DO J=1,3 AIM(I,J)=A(I,J) enddo enddo IF (ISING.EQ.1) THEN WRITE(*,*)'MATRICE AIM singuliere ds ACTI2' ENDIF DO I=1,3 DO J=1,3 AI(I,J)=AIM(I,J) enddo enddo AI(4,4) = 1. / A(4,4) AI(5,5) = 1. / A(5,5) AI(6,6) = 1. / A(6,6) C C -------------- mise a jour des contraintes ------------ C 75 CONTINUE C DO 80 I=1,NSTRSS DEPS(I)= SIGE(I) - (DLAM1 * VEC1(I)) 80 CONTINUE C DO I=1,NSTRSS SIGF(I)=0. DO J=1,NSTRSS SIGF(I)=SIGF(I)+AI(I,J)*DEPS(J) enddo enddo C C ******** VERIFICATION DU CRITERE **************** C FCRI1 = SEQI - SEQ C IF (IBROY.EQ.0.AND.ABS(FCRI1).GE.CRIMAX) THEN C WRITE(*,*)'****************************************' C WRITE(*,*)'LE RESIDU DIVERGE AVEC BROYDEN' C WRITE(*,*)'on passe donc a la secante' C WRITE(*,*)'Dans l element',IBB C WRITE(*,*)'et au point d intégration',IGAU C WRITE(*,*)'CRIMAX=',CRIMAX C WRITE(*,*)'****************************************' ITER=ITR ENDIF C C ******* Compteur sur la methode de resolution **** C IF (IBROY.EQ.0.AND.ITER.EQ.ITR) THEN IBROY=1 ITANG=1 ITER=1 IAPEX=0 GOTO 8 ENDIF C C ******* non convergence ************************** C IF (ABS(FCRI1).GT.PRB.AND.ITER.LT.ITR) THEN IF (IBROY.EQ.0) THEN DJAC1=(FCRI0-FCRI1)/(DLAM0-DLAM1) DLAM0=DLAM1 DJAC0=DJAC1 FCRI0=FCRI1 ITER=ITER+1 GOTO 40 ENDIF IF (IBROY.EQ.1.AND.ITANG.EQ.1) THEN DLAM0=DLAM1 FCRI0=FCRI1 ITER=ITER+1 IF (ITER.GE.(ITR-5)) THEN WRITE(*,*)'ITER=',ITER WRITE(*,*)'FCRI=',FCRI1 ENDIF GOTO 9 ENDIF ENDIF IF (ITER.GE.ITR.AND.ABS(FCRI1).GT.PRB2) THEN WRITE(*,*)'NON CONVERGENCE INTERNE dans ACTI2' WRITE(*,*)'FCRI=',FCRI1 WRITE(*,*)'Dans l element',IBB WRITE(*,*)'et au point d intégration',IGAU STOP ENDIF C C ******* convergence ************************** C SEQC = SEQ DKC = DK C C *********************************************** C 100 CONTINUE C RETURN END
© Cast3M 2003 - Tous droits réservés.
Mentions légales