Sources Esope | Cast3M

copain
C COPAIN    SOURCE    FD144363  26/03/18    21:15:03     12501          
      SUBROUTINE COPAIN(IPT1,IPT2,IPT3,IPT4,
     &  MPOVA1,MPOVA2,MPOVA3,MPOVA4,MPOVA5,MPOVA6,
     &  MPOVA7,MPOVA8,MPOVA9,
     &  MPOVAB,MPOVAC,MPOVAD,MPOVAE,
     &  GRAV,HMAX,TTAINF,TTASUP,MPOVACP,MPOVAF,MPOVAI,MPOVAJ,
     &  MPOUT1,MPOUT2,MPOUT3,MPOUT4)
C--------------------------------------------------------------------
C Evaluation du coef. d'echange thermique en paroi et du debit de
C condensation lorsque la correlation retenue est celle de COPAIN
C--------------------------------------------------------------------
C
C  E/   IPT1   : Maillage des points CELL
C  E/   IPT2   : Maillage des points LAYERW
C  E/   IPT3   : Maillage des points WALL
C  E/   IPT4   : Maillage de SEG2 CELL/LAYERW
C  E/   MPOVA1 : Pression totale
C  E/   MPOVA2 : Temperature du melange gazeux
C  E/   MPOVA3 : Densitee du melange
C  E/   MPOVA4 : Densitee partielle de vapeur
C  E/   MPOVA5 : Pression partielle de vapeur
C  E/   MPOVA6 : Temperature de saturation Tsat(Pv)
C  E/   MPOVA7 : Coefficient d'echange par conduction
C  E/   MPOVA8 : Temperature des murs
C  E/   MPOVA9 : Surface des parois
C  E/   MPOVAB : Lambda du melange
C  E/   MPOVAC : Nu du melange
C  E/   MPOVAD : Constante des gaz parfait des incondensables
C  E/   MPOVAE : Coeff de diffusion de la vapeur
C  E/   MPOVAF : Constante des gaz parfait du melange
C  E/   MPOVAI : Masse molaire des incondensables
C  E/   MPOVAJ : Masse molaire de la vapeur d'eau
C  E/   GRAV   : Gravite (acceleration de la pesanteur)  
C  E/   HMAX   : Borne sup pour le coefficient de transfert thermique
C  E/   TTAINF : Borne inf pour le coefficient theta
C  E/   TTASUP : Borne sup pour le coefficient theta
C  E/   MPOVACP: Coefficient calorifique de l'eau
C   /S  MPOUT1 : Coefficient d'echange global en paroi
C   /S  MPOUT2 : Debit de condensation / flux thermique
C   /S  MPOUT3 : Enthalpie liquide
C   /S  MPOUT4 : Coefficient d'échange convectif en paroi
C
C--------------------------------------------------------------------
C Les segments MPOUT... ne sont pas desactives
C--------------------------------------------------------------------
C
C 13/12/16 : Création
C 2026/03  : FD : On borne THETA (+ nettoyage)
C
C--------------------------------------------------------------------
C
C Langage : ESOPE + FORTRAN77
C Auteurs : 
C
C-----------------------------------------------------------------------
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)                             
      POINTEUR MPOVA7.MPOVAL,MPOVA8.MPOVAL,MPOVA9.MPOVAL
      POINTEUR MPOVAB.MPOVAL,MPOVAC.MPOVAL,MPOVAD.MPOVAL,MPOVAE.MPOVAL
      POINTEUR MPOVAF.MPOVAL,MPOVAI.MPOVAL,MPOVAJ.MPOVAL,MPOVACP.MPOVAL
      POINTEUR MPOUT1.MPOVAL,MPOUT2.MPOVAL,MPOUT3.MPOVAL,MPOUT4.MPOVAL
C
 
-INC PPARAM
-INC CCOPTIO
-INC SMCHPOI
-INC SMCOORD
-INC SMELEME
-INC SMTABLE
      SEGMENT REDIR
        INTEGER IND1(NNGOT),IND2(NNGOT),IND3(NNGOT)
      ENDSEGMENT
C
C- Tableau d'adressage indirect
C
      SEGACT IPT1,IPT2,IPT3,IPT4
      SEGACT MPOVA1,MPOVA2,MPOVA3,MPOVA4,MPOVA5,MPOVA6,MPOVA7,MPOVA8,
     &       MPOVA9,MPOVAB,MPOVAC,MPOVAD,MPOVAE,MPOVACP,MPOVAF,
     &       MPOVAI,MPOVAJ
      NNGOT  = nbpts
      SEGINI REDIR
      NBLAY  = IPT2.NUM(/2)
      DO 10 I=1,NBLAY
         K = IPT2.NUM(1,I)
         IND2(K) = I
 10   CONTINUE
      NBWAL  = IPT3.NUM(/2)
      DO 20 I=1,NBWAL
         K = IPT3.NUM(1,I)
         IND3(K) = I
 20   CONTINUE
      NBCEL  = IPT1.NUM(/2)
      DO 30 I=1,NBCEL
         K = IPT1.NUM(1,I)
         IND1(K) = I
 30   CONTINUE
      UST    = 1.D0 / 3.D0
      DST    = 2.D0 / 3.D0
C
C- Pour chaque couple CELL/LAYERW on calcule les grandeurs liees a la
C  condensation
C
      SEGINI, MPOUT1=MPOVA7
      N  = MPOVA7.VPOCHA(/1)
      NC = 1
      SEGINI MPOUT2,MPOUT3,MPOUT4
      N1 = IPT4.NUM(/2)
      DO 70 I=1,N1
       IPOS = IND2(IPT4.NUM(2,I))
       IF (IPOS.NE.0) THEN
          II1     = IND1(IPT4.NUM(1,I))
          JPOS    = IND3(IPT4.NUM(2,I))
C
          T_G     = MPOVA2.VPOCHA(II1,1)
          Tsat_G  = MPOVA6.VPOCHA(II1,1)
          P_G     = MPOVA1.VPOCHA(II1,1)
          Pv_G    = MPOVA5.VPOCHA(II1,1)
          RHO_G   = MPOVA3.VPOCHA(II1,1)
          RHOv_G  = MPOVA4.VPOCHA(II1,1)     
C
          Rnc     = MPOVAD.VPOCHA(II1,1)
          T_W     = MPOVA8.VPOCHA(JPOS,1)
          Pv_W    = PSATT0(T_W)
          Pnc_W   = P_G - Pv_W
          RHOnc_W = Pnc_W / (Rnc * T_W)
          RHOv_W  = ROVAP0(Pv_W,T_W)
          RHO_W   = RHOv_W + RHOnc_W
C          
          HL_W    = HLS0(P_G,T_W)
          XLAMBDA = MPOVAB.VPOCHA(II1,1)
          XNU     = MPOVAC.VPOCHA(II1,1)
          CPmel   = MPOVACP.VPOCHA(II1,1)
          ALFA    = XLAMBDA / (RHO_G * CPmel) 
          WNC     = MPOVAI.VPOCHA(II1,1)
          WVAP    = MPOVAJ.VPOCHA(II1,1)
C
          Ync_G = (RHO_G - RHOv_G) / RHO_G
          Ync_W = RHOnc_W / RHO_W
          X1    = (Ync_W - Ync_G) / 0.01D0
          X1    = MAX(X1,0.D0)
          X1    = MIN(X1,1.D0)
          FX1   = (X1**2.D0) * (3.D0 - 2.D0*X1)
          FACT1 = (T_G - T_W) / T_G
          FACT2 = FX1 * (Ync_W - Ync_G) / ((WNC/(WNC-WVAP)) - Ync_G)
          GBDRSR = ABS(FACT1 + FACT2) * GRAV
C
C En absence de condensation HT n'est pas multiplié par theta
          HT     = 0.13D0 * XLAMBDA * ((GBDRSR/(XNU*ALFA)) ** UST)
          IF ((T_G - T_W) .GE. 0.02D0) THEN
             S_W   = MPOVA9.VPOCHA(IPOS,1)
             Dvmel = MPOVAE.VPOCHA(II1,1)
             HV_G  = HVS0(Pv_G,T_G)
C
             IF (Ync_w .LT. 0.01D0) THEN
                DYSY  = -1.D0
             ELSE
                DYSY  = (Ync_W - Ync_G) / Ync_W
             ENDIF 
             IF (DYSY .LT. 0.D0) THEN
                DC  = 0.D0
                HDC = 0.D0
             ELSE
                THETA = 0.8254D0 + 0.616D0 * (Pv_G-Pv_W) / (P_G-Pv_W)
                THETA = MAX(THETA,TTAINF)
                THETA = MIN(THETA,TTASUP)
                HT    = MIN(THETA*HT,HMAX)
                XKC   = (Dvmel ** DST) * (ALFA ** UST) * HT / XLAMBDA
C
C Lissage de Dc en deça d'un écart de 0.2°C entre Tsat_G et Tsat_W
                X2    = (Tsat_G - T_W) / 0.2D0
                X2    = MAX(X2,0.D0)
                X2    = MIN(X2,1.D0)
                FX2   = (X2**2.D0) * (3.D0 - 2.D0*X2)
                DC    = RHO_G * XKC * S_W * DYSY * FX2
                DC    = MAX(DC,0.D0)
                HDC   = (DC / S_W) * (HV_G - HL_W) / (T_G - T_W)
                HDC   = MAX(HDC,0.D0)
             ENDIF 
             IF ((HT + HDC) . LT . MPOVA7.VPOCHA(IPOS,1)) THEN
                HT = MPOVA7.VPOCHA(IPOS,1) - HDC
             ENDIF
             QW2 = (HT + HDC) * (T_G - T_W)
             DCSPHI = (DC / QW2) * (-1.D0)
          ELSE
             HT = MIN(HT,HMAX)
             DC  = 0.D0
             HDC = 0.D0
             IF ((HT + HDC) . LT . MPOVA7.VPOCHA(IPOS,1)) THEN
                HT = MPOVA7.VPOCHA(IPOS,1) - HDC
             ENDIF
             DCSPHI = 0.D0
          ENDIF
          HT = MAX(HT,0.D0)
          MPOUT1 . VPOCHA(IPOS,1) = HT + HDC
          MPOUT2 . VPOCHA(IPOS,1) = DCSPHI
          MPOUT3 . VPOCHA(IPOS,1) = HL_W
          MPOUT4 . VPOCHA(IPOS,1) = HT
       ENDIF
  70  CONTINUE
C
C- Menage
C
      SEGSUP REDIR
      SEGDES IPT1,IPT2,IPT3,IPT4,MPOVA1,MPOVA2,MPOVA3,MPOVA4,MPOVA5,
     &       MPOVA6,MPOVA7,MPOVA8,MPOVA9,MPOVAB,MPOVAC,MPOVAD,
     &       MPOVAE,MPOVAF,MPOVAI,MPOVAJ,MPOVACP
C
      RETURN
      END