Sources Esope | Cast3M

triple
C TRIPLE    SOURCE    CHAT      05/01/13    03:47:39     5004
      SUBROUTINE TRIPLE(ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,YUNG,XNU,SIGEL,DSIGP,DLAMBD,
     .XLAMBD,ICOUP,ICAS,ICRIT1,ICRIT2,ICRIME,ICRIMT,SIGMA,DSIGMA,XX,
     .PREC,RFSG,RFEP,RFPR,KERRE)
C
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
-INC PPARAM
-INC CCOPTIO
C
      DIMENSION SIGEL(*),DSIGP(*),DLAMBD(*),XLAMBD(*)
      DIMENSION SIGMA(*),DSIGMA(*)
      DIMENSION DF1(6),DF2(6),DF3(6),SIGMAT(6),SIG(6)
      DIMENSION A(3,3),C(3),DDLAM(3)
C
C  IC1  INDICE DU CRITERE 1
C  IC2  INDICE DU CRITERE 2
C  IC3  INDICE DU CRITERE 3
C  DLAMBD(IC1)  INCREMENT DLMBDA POUR LE CRITERE 1
C  DLAMBD(IC2)  INCREMENT DLMBDA POUR LE CRITERE 2
C  DLAMBD(IC3)  INCREMENT DLMBDA POUR LE CRITERE 3
C  DDLAM1  CORRECTION A DLAMBD(IC1) AU COURS DES ITERATIONS INTERNES
C  DDLAM2  CORRECTION A DLAMBD(IC2) AU COURS DES ITERATIONS INTERNES
C  DDLAM3  CORRECTION A DLAMBD(IC3) AU COURS DES ITERATIONS INTERNES
C  SIG  ETAT DE CONTRAINTES FINAL PROJETTE
C  SIGEL  ETAT DE CONTRAINTES INITIAL
C  DSIGP  INCREMENT DE CONTRAINTES A ECOULER
C  MAXITE  NOMBRE MAXIMAL D'ITERATIONS INTERNES
C  PREC  PRECISION POUR LA CONVERGENCE DES ITERATIONS INTERNES
C
      DO 1 I=1,6
      SIGMAT(I)=SIGEL(I)+DSIGP(I)
      SIGMA(I)=SIGEL(I)
    1 CONTINUE
C
      IF(IIMPI.EQ.9) THEN
        WRITE(IOIMP,3001) (SIGEL(I),I=1,6)
        WRITE(IOIMP,3002) (DSIGP(I),I=1,6)
        WRITE(IOIMP,3003) (SIGMAT(I),I=1,6)
      ENDIF
C
      DATA MAXITE/30/
      IC1=1
      IC2=2
      IC3=3
      DO 2 I=1,3
      DDLAM(I)=0.D0
    2 CONTINUE
      TRSIGE=SIGEL(1)+SIGEL(2)+SIGEL(3)
      SIGEQ2=AVM(SIGEL,SIGEL)
      SIGEQ=SQRT(SIGEQ2)
C
      SU=1.5*YUNG/(1.D0+XNU)
      DU=YUNG/(1.D0-2.D0*XNU)
C
C CALCUL DU: A11, A22, A33, A12, A23, A31
C
      A(1,1)=PENTE3+DU/3.D0
      A(2,2)=SU+ALFADV*ALFADV*DU*3.D0
      A(3,3)=PENTE1+SU
      A(1,2)=-ALFADV*DU
      A(2,1)=A(1,2)
      A(2,3)=SU
      A(3,2)=SU
      A(3,1)=0.D0
      A(1,3)=0.D0
C
C CALCUL DU: C(1), C(2), C(3)
C
      CALL KRITER(IC1,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,DDLAM(1),SIGEL,FSG,C(1),KERRE)
      CALL KRITER(IC2,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,DDLAM(2),SIGEL,FSG,C(2),KERRE)
      CALL KRITER(IC3,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,DDLAM(3),SIGEL,FSG,C(3),KERRE)
C
      IF(IIMPI.EQ.9) WRITE(IOIMP,3004) (C(I),I=1,3)
C
      TRDSIG=DSIGP(1)+DSIGP(2)+DSIGP(3)
      FF=SIGEL(1)*DSIGP(1)+SIGEL(2)*DSIGP(2)+SIGEL(3)*DSIGP(3)
      ZZ=2.D0*(SIGEL(4)*DSIGP(4)+SIGEL(5)*DSIGP(5)+SIGEL(6)*DSIGP(6))
      TRSIDS=FF+ZZ-TRSIGE*TRDSIG/3.D0
      C(1)=C(1)-TRDSIG/3.D0
      C(2)=C(2)+ALFADV*TRDSIG+1.5D0*TRSIDS/SIGEQ
      C(3)=C(3)+1.5D0*TRSIDS/SIGEQ
C
C CALCUL DU: DDLAM(1),DDLAM(2),DDLAM(3)
C
      CALL RESNEQ(A,C,DDLAM,3,KERRE)
C
      IF(IIMPI.EQ.9) THEN
        WRITE(IOIMP,3004) (C(I),I=1,3)
        WRITE(IOIMP,3005) (DDLAM(I),I=1,3)
      ENDIF
C
C  INITIALISATIONS
C
      ZR=0.D0
      ITER=0
      ICOUP=0
      ICRIT1=0
      ICRIT2=0
      ICRIME=0
      ICRIMT=0
      IBOU=6
      VMELA0=VMELAS
      DPEL20=DPELA2
      POREL0=PORELA
      DO 3 I=1,IBOU
    3 SIG(I)=SIGEL(I)
      DLAM01=DLAMBD(IC1)
      DLAM02=DLAMBD(IC2)
      DLAM03=DLAMBD(IC3)
      DLAMBD(IC1)=DDLAM(1)
      DLAMBD(IC2)=DDLAM(2)
      DLAMBD(IC3)=DDLAM(3)
C
      IF(ABS(DLAMBD(IC1)).LE.RFEP) DLAMBD(IC1)=0.D0
      IF(ABS(DLAMBD(IC2)).LE.RFEP) DLAMBD(IC2)=0.D0
      IF(ABS(DLAMBD(IC3)).LE.RFEP) DLAMBD(IC3)=0.D0
C
      IF(DLAMBD(IC1).GT.0.D0.AND.DLAMBD(IC2).GT.0.D0.
     .   AND.DLAMBD(IC3).GT.0.D0) GO TO 9
      IF(DLAMBD(IC1).GT.0.D0.OR.DLAMBD(IC2).GT.0.D0.
     .   OR.DLAMBD(IC3).GT.0.D0) GO TO 309
      IF(DLAMBD(IC1).EQ.0.D0.AND.DLAMBD(IC2).EQ.0.D0.
     .   AND.DLAMBD(IC3).EQ.0.D0) GO TO 313
      WRITE(IOIMP,901) IC1,DLAMBD(IC1),IC2,DLAMBD(IC2),IC3,DLAMBD(IC3)
      KERRE=640
      RETURN
C
    9 CONTINUE
C
C LE CALCUL DE DLAMBD PREMIERE ESTIMATION
C
      ITER=1
      CALL KRITER(IC1,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,ZR,SIGMAT,FSG,C(1),KERRE)
      CALL KRITER(IC2,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,ZR,SIGMAT,FSG,C(2),KERRE)
      CALL KRITER(IC3,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,ZR,SIGMAT,FSG,C(3),KERRE)
C
C CALCUL DU: A11, A22, A33, A12, A23, A31
C
      A(1,1)=PENTE3+DU/3.D0
      A(2,2)=SU+ALFADV*ALFADV*DU*3.D0
      A(3,3)=PENTE1+SU
      A(1,2)=-ALFADV*DU
      A(2,1)=A(1,2)
      A(2,3)=SU
      A(3,2)=SU
      A(3,1)=0.D0
      A(1,3)=0.D0
C
      CALL RESNEQ(A,C,DDLAM,3,KERRE)
C
      IF(IIMPI.EQ.9) THEN
        WRITE(IOIMP,3004) (C(I),I=1,3)
        WRITE(IOIMP,3005) (DDLAM(I),I=1,3)
      ENDIF
C
      DLAMBD(IC1)=DDLAM(1)
      DLAMBD(IC2)=DDLAM(2)
      DLAMBD(IC3)=DDLAM(3)
C
      CALL DDFDS(IC1,ALFADV,ALFAD1,ALFAD2,YUNG,XNU,SIGMAT,SIG,
     .DF1,KERRE)
      CALL DDFDS(IC2,ALFADV,ALFAD1,ALFAD2,YUNG,XNU,SIGMAT,SIG,
     .DF2,KERRE)
      CALL DDFDS(IC3,ALFADV,ALFAD1,ALFAD2,YUNG,XNU,SIGMAT,SIG,
     .DF3,KERRE)
      DO 4 I=1,IBOU
      SIG(I)=SIG(I)-DF1(I)*DDLAM(1)-DF2(I)*DDLAM(2)-DF3(I)*DDLAM(3)
      SIG(I)=SIG(I)+DSIGP(I)
    4 CONTINUE
C
      IF(IIMPI.EQ.9) THEN
        WRITE(IOIMP,3006) (DF1(I),I=1,6)
        WRITE(IOIMP,3007) (DF2(I),I=1,6)
        WRITE(IOIMP,3008) (DF3(I),I=1,6)
        WRITE(IOIMP,3009) (SIG(I),I=1,6)
      ENDIF
C
C  ITERATIONS INTERNES
C
  333 ITER=ITER+1
      CALL KRITER(IC1,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,DDLAM(1),SIG,FSG,C(1),KERRE)
      CALL KRITER(IC2,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,DDLAM(2),SIG,FSG,C(2),KERRE)
      CALL KRITER(IC3,ALFADV,ALFAD1,ALFAD2,VMELAS,DPELAS,DPELA1,
     .DPELA2,PORELA,PENTE1,PENTE2,PENTE3,DDLAM(3),SIG,FSG,C(3),KERRE)
C
C CALCUL DU: A11, A22, A33, A12, A23, A31
C
      A(1,1)=PENTE3+DU/3.D0
      A(2,2)=SU+ALFADV*ALFADV*DU*3.D0
      A(3,3)=PENTE1+SU
      A(1,2)=-ALFADV*DU
      A(2,1)=A(1,2)
      A(2,3)=SU
      A(3,2)=SU
      A(3,1)=0.D0
      A(1,3)=0.D0
C
      CALL RESNEQ(A,C,DDLAM,3,KERRE)
C
      IF(IIMPI.EQ.9) THEN
        WRITE(IOIMP,3004) (C(I),I=1,3)
        WRITE(IOIMP,3005) (DDLAM(I),I=1,3)
      ENDIF
C
      CALL DDFDS(IC1,ALFADV,ALFAD1,ALFAD2,YUNG,XNU,SIG,SIG,
     .DF1,KERRE)
      CALL DDFDS(IC2,ALFADV,ALFAD1,ALFAD2,YUNG,XNU,SIG,SIG,
     .DF2,KERRE)
      CALL DDFDS(IC3,ALFADV,ALFAD1,ALFAD2,YUNG,XNU,SIG,SIG,
     .DF3,KERRE)
      DO 5 I=1,IBOU
      SIG(I)=SIG(I)-DF1(I)*DDLAM(1)-DF2(I)*DDLAM(2)-DF3(I)*DDLAM(3)
    5 CONTINUE
C
      IF(IIMPI.EQ.9) THEN
        WRITE(IOIMP,3006) (DF1(I),I=1,6)
        WRITE(IOIMP,3007) (DF2(I),I=1,6)
        WRITE(IOIMP,3008) (DF3(I),I=1,6)
        WRITE(IOIMP,3009) (SIG(I),I=1,6)
      ENDIF
C
      DLAMBD(IC1)=DLAMBD(IC1)+DDLAM(1)
      DLAMBD(IC2)=DLAMBD(IC2)+DDLAM(2)
      DLAMBD(IC3)=DLAMBD(IC3)+DDLAM(3)
C
      IF(ABS(DLAMBD(IC1)).LE.RFEP) DLAMBD(IC1)=0.D0
      IF(ABS(DLAMBD(IC2)).LE.RFEP) DLAMBD(IC2)=0.D0
      IF(ABS(DLAMBD(IC3)).LE.RFEP) DLAMBD(IC3)=0.D0
C
      IF(DLAMBD(IC1).GT.0.D0.AND.DLAMBD(IC2).GT.0.D0.
     .   AND.DLAMBD(IC3).GT.0.D0) GO TO 310
      IF(DLAMBD(IC1).GT.0.D0.OR.DLAMBD(IC2).GT.0.D0.
     .   OR.DLAMBD(IC3).GT.0.D0) GO TO 309
      IF(DLAMBD(IC1).EQ.0.D0.AND.DLAMBD(IC2).EQ.0.D0.
     .   AND.DLAMBD(IC3).EQ.0.D0) GO TO 313
      WRITE(IOIMP,901) IC1,DLAMBD(IC1),IC2,DLAMBD(IC2),IC3,DLAMBD(IC3)
      KERRE=640
      RETURN
C
C  TESTS
C
  310 TEST1=ABS(DDLAM(1))/DLAMBD(IC1)
      TEST2=ABS(DDLAM(2))/DLAMBD(IC2)
      TEST3=ABS(DDLAM(3))/DLAMBD(IC3)
      IF(TEST1.LT.PREC.AND.TEST2.LT.PREC.AND.TEST3.LT.PREC) GO TO 311
      IF(ITER.LE.MAXITE) GO TO 333
      WRITE(IOIMP,900) IC1,TEST1,IC2,TEST2,IC3,TEST3
      KERRE=640
      RETURN
C
C IL Y A TRIPLAGE
C
  311 ICOUP=3
      ICRIT1=IC1
      ICRIT2=IC2
      ICRIT3=IC3
C
      CALL KRITER(IC1,ALFADV,ALFAD1,ALFAD2,VMELA0,DPELAS,DPELA1,
     .DPEL20,POREL0,PENTE1,PENTE2,PENTE3,DLAMBD(IC1),SIG,FSG,C(1),KERRE)
      CALL KRITER(IC2,ALFADV,ALFAD1,ALFAD2,VMELA0,DPELAS,DPELA1,
     .DPEL20,POREL0,PENTE1,PENTE2,PENTE3,DLAMBD(IC2),SIG,FSG,C(2),KERRE)
      CALL KRITER(IC3,ALFADV,ALFAD1,ALFAD2,VMELA0,DPELAS,DPELA1,
     .DPEL20,POREL0,PENTE1,PENTE2,PENTE3,DLAMBD(IC3),SIG,FSG,C(3),KERRE)
C
      VMELAS=VMELA0
      DPELA2=DPEL20
      PORELA=POREL0
C
      DLAMBD(IC1)=DLAMBD(IC1)+DLAM01
      DLAMBD(IC2)=DLAMBD(IC2)+DLAM02
      DLAMBD(IC3)=DLAMBD(IC3)+DLAM03
      DO 7 I=1,5
      XLAMBD(I)=XLAMBD(I)+DLAMBD(I)
      DLAMBD(I)=0.D0
    7 CONTINUE
      DO 8 I=1,IBOU
      SIGEL(I)=SIG(I)
      DSIGP(I)=0.D0
      DSIGMA(I)=SIGEL(I)-SIGMA(I)
    8 CONTINUE
      RETURN
C
  309 ICAS=0
      IF(DLAMBD(IC1).GT.0.D0.AND.DLAMBD(IC2).GT.0.D0) ICAS=IC1+IC2-2
      IF(DLAMBD(IC1).GT.0.D0.AND.DLAMBD(IC3).GT.0.D0) ICAS=IC1+IC3-2
      IF(DLAMBD(IC2).GT.0.D0.AND.DLAMBD(IC3).GT.0.D0) ICAS=IC2+IC3-2
      IF(ICAS.NE.0) GO TO 312
C
C IL N Y A PAS NI DE TRIPLAGE NI DE COUPLAGE
C
      IF(DLAMBD(IC1).GT.0.D0) ICRIT1=IC1
      IF(DLAMBD(IC2).GT.0.D0) ICRIT1=IC2
      IF(DLAMBD(IC3).GT.0.D0) ICRIT1=IC3
C
      ICRIME=IC1
      ICRIMT=IC2
      IF(ICRIT1.EQ.IC1) ICRIME=IC3
      IF(ICRIT1.EQ.IC2) ICRIMT=IC3
C
      ICOUP=1
      XX=0.D0
      DLAMBD(IC1)=DLAM01
      DLAMBD(IC2)=DLAM02
      DLAMBD(IC3)=DLAM03
      VMELAS=VMELA0
      DPELA2=DPEL20
      PORELA=POREL0
      RETURN
C
C IL N Y A PAS DE TRIPLAGE MAIS IL Y A COUPLAGE
C
  312 ICOUP=2
      IF(DLAMBD(IC1).EQ.0.D0) ICRIMT=IC1
      IF(DLAMBD(IC2).EQ.0.D0) ICRIMT=IC2
      IF(DLAMBD(IC3).EQ.0.D0) ICRIMT=IC3
C
      XX=0.D0
      DLAMBD(IC1)=DLAM01
      DLAMBD(IC2)=DLAM02
      DLAMBD(IC3)=DLAM03
      VMELAS=VMELA0
      DPELA2=DPEL20
      PORELA=POREL0
      RETURN
C
C IL N Y A PAS D ENDOMMAGEMENT
C
  313 ICOUP=0
      XX=0.D0
      DLAMBD(IC1)=DLAM01
      DLAMBD(IC2)=DLAM02
      DLAMBD(IC3)=DLAM03
      VMELAS=VMELA0
      DPELA2=DPEL20
      PORELA=POREL0
C
  900 FORMAT(1X,'ERREUR DANS TRIPLE - NON CONVERGENCE '/
     .       1X,'ICRIT1=',I4,1X,'TEST1 =',1PE12.5,/,
     .       1X,'ICRIT2=',I4,1X,'TEST2 =',1PE12.5,/,
     .       1X,'ICRIT3=',I4,1X,'TEST3 =',1PE12.5)
  901 FORMAT(1X,'ERREUR DANS TRIPLE: DLAMBD(',I1,')=',1PD12.5,1X,
     .'DLAMBD(',I1,')=',1PD12.5,1X,'DLAMBD(',I1,')=',1PD12.5)
 3001 FORMAT(1X,'SIGEL =',6(1X,1PD12.5))
 3002 FORMAT(1X,'DSIGP =',6(1X,1PD12.5))
 3003 FORMAT(1X,'SIGMAT=',6(1X,1PD12.5))
 3004 FORMAT(1X,'C     =',3(1X,1PD12.5))
 3005 FORMAT(1X,'DDLAM =',3(1X,1PD12.5))
 3006 FORMAT(1X,'DF1   =',6(1X,1PD12.5))
 3007 FORMAT(1X,'DF2   =',6(1X,1PD12.5))
 3008 FORMAT(1X,'DF3   =',6(1X,1PD12.5))
 3009 FORMAT(1X,'SIG   =',6(1X,1PD12.5))
C
      RETURN
      END