Sources Esope | Cast3M

cldeta
C CLDETA    SOURCE    CHAT      05/01/12    22:04:36     5004
      SUBROUTINE CLDETA (YUNG,XNU,RT,XLTR,XLTT,EPTT,EPTR,EPRS,OUVER,
     .                   CONTR,DDEFT,DCONPR)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
 
-INC PPARAM
-INC CCOPTIO
C
      DIMENSION RT(*),XLTR(*),XLTT(*),EPTT(*),EPTR(*),EPRS(*),OUVER(*)
      DIMENSION CONTR(*),DDEFT(*),DCONPR(*)
      DIMENSION DDR(6,6),ET(3),EPPT(3),EPSPR(3)
C
C***********************************************************************
C************* CALCUL DE L ETAT DE CHARGE OU DE DECHARGE ***************
C***********************************************************************
C
C     INITIALISATION
C
      CALL ZDANUL(DDR,36)
      RFE=1.D-8
      RFS=YUNG*RFE
C
      DO 10 I=1,3
        ET(I)=YUNG
        IF(XLTR(I).LE.0.D0) ET(I)=0.D0
        EPPT(I)=EPTT(I)-XLTT(I)/YUNG
        EPSPR(I)=EPRS(I)
        IF(XLTR(I).GT.0.D0.AND.RT(I).LE.RFS) EPSPR(I)=EPTR(I)
C
        IF(XLTR(I).GT.0.D0.AND.RT(I).GT.RFS) THEN
C
          IF(XLTT(I).GT.0.D0) THEN
            IF(RT(I).GE.XLTT(I)) THEN
              EPSPR(I)=((XLTR(I)-RT(I))/(XLTR(I)-XLTT(I)))*EPPT(I)
            ELSE
              EPSPR(I)=EPPT(I)+((1.D0-RT(I)/XLTT(I))*(EPTR(I)-EPPT(I)))
            ENDIF
          ELSE
            EPSPR(I)=(1.D0-RT(I)/XLTR(I))*EPTR(I)
          ENDIF
C
          EPPMAX=(1.D0-RT(I)/XLTR(I))*EPSPR(I)
          EPRM=(1.D0-RT(I)/XLTR(I))*EPRS(I)
          EPRM=MIN(EPPMAX,EPRM)
        ELSE
            EPPMAX=EPSPR(I)
            EPRM=EPRS(I)
        ENDIF
C
        IF(XLTR(I).GT.0.D0)
     .     ET(I)=(YUNG*XLTR(I))/(XLTR(I)+(YUNG*EPSPR(I)))
 
C
        IF(ABS(CONTR(I)).LE.RFS) THEN
          IF(RT(I).LE.RFS) THEN
            IF((ABS(OUVER(I)-EPRM)).LE.RFE.OR.OUVER(I).LE.EPRM) THEN
              IF(DDEFT(I).LT.0.D0) THEN
                ET(I)=YUNG
              ELSE
                ET(I)=0.D0
              ENDIF
            ELSE
              ET(I)=0.D0
            ENDIF
          ELSE
            IF((ABS(OUVER(I)-EPPMAX)).LE.RFE.OR.OUVER(I).GE.EPPMAX) THEN
              IF(DDEFT(I).LT.0.D0) ET(I)=0.D0
            ELSE
              IF((ABS(OUVER(I)-EPRM)).LE.RFE.OR.OUVER(I).LE.EPRM) THEN
                IF(DDEFT(I).LT.0.D0) ET(I)=YUNG
              ELSE
                IF(DDEFT(I).LT.0.D0) ET(I)=0.D0
              ENDIF
            ENDIF
          ENDIF
        ENDIF
C
        IF((ABS(XLTR(I)-RT(I))).LE.RFS.AND.XLTR(I).GT.0.D0) ET(I)=YUNG
C
   10 CONTINUE
C
      ET1=ET(1)
      ET2=ET(2)
      ET3=ET(3)
      UPUN=1.D0+XNU
      UMUN=1.D0-XNU
      UMDN=1.D0-2.D0*XNU
      EMET1=YUNG-ET1
      EMET2=YUNG-ET2
      EMET3=YUNG-ET3
C
      DENO=(EMET1*EMET2*EMET3)+
     .   (ET1*EMET2*EMET3)+(EMET1*ET2*EMET3)+(EMET1*EMET2*ET3)+
     .   (UPUN*UMUN*((EMET1*ET2*ET3)+(ET1*EMET2*ET3)+(ET1*ET2*EMET3)))+
     .   (UPUN*UPUN*UMDN*(ET1*ET2*ET3))
      USDENO=1.D0/DENO
C
      DDR(1,1)=((EMET2*EMET3)+(ET2*EMET3)+(EMET2*ET3)+
     .      (ET2*ET3*UPUN*UMUN))*YUNG*ET1*USDENO
      DDR(1,2)=XNU*YUNG*ET1*ET2*(EMET3+(ET3*UPUN))*USDENO
      DDR(1,3)=XNU*YUNG*ET1*ET3*(EMET2+(ET2*UPUN))*USDENO
      DDR(2,1)=DDR(1,2)
      DDR(2,2)=((EMET1*EMET3)+(ET1*EMET3)+(EMET1*ET3)+
     .      (ET1*ET3*UPUN*UMUN))*YUNG*ET2*USDENO
      DDR(2,3)=XNU*YUNG*ET2*ET3*(EMET1+(ET1*UPUN))*USDENO
      DDR(3,1)=DDR(1,3)
      DDR(3,2)=DDR(2,3)
      DDR(3,3)=((EMET1*EMET2)+(ET1*EMET2)+(EMET1*ET2)+
     .      (ET1*ET2*UPUN*UMUN))*YUNG*ET3*USDENO
      DDR(4,4)=YUNG/(2.D0*UPUN)
      DDR(5,5)=DDR(4,4)
      DDR(6,6)=DDR(4,4)
C
      CALL MULMAT(DCONPR,DDR,DDEFT,6,1,6)
C
      IF(IIMPI.EQ.9) THEN
        WRITE(IOIMP,*) 'DDR'
        WRITE(IOIMP,1000) ((DDR(I,J),J=1,6),I=1,6)
        WRITE(IOIMP,1001) (DDEFT(I),I=1,6)
        WRITE(IOIMP,1002) (DCONPR(I),I=1,6)
 1000   FORMAT(6(1X,1PE12.5))
 1001   FORMAT(1X,'DDEFT =',6(1X,1PE12.5))
 1002   FORMAT(1X,'DCONPR=',6(1X,1PE12.5))
      ENDIF
C
      RETURN
C
      END