Esope source files

   1 : C BSUR      SOURCE    DURA      06/07/11    21:17:29     5481           
   2 :       SUBROUTINE BSUR (X,DX,XL,RUG,XW,XN,TN,EN,BN,KIMP,PSLIM,REL,RINDEX,
   3 :      &         P1,T1,QAE,QEE1,PHI1,P2,T2,U2,QEE2,PHI2,QW2,RE,H,PSQ,
   4 :      &         NPP,ITP,PF,PP,DPF,DPP,RECU,XKUL,XKUT1,XKUT2,XKUT3,XKUT4)
   5 : 
   6 :       IMPLICIT INTEGER(I-N)
   7 :       IMPLICIT REAL*8 (A-H,O-Z)
   8 :       DIMENSION XN(NPP),TN(NPP),EN(NPP),BN(NPP)
   9 : 
  10 : C     operateur FUITE
  11 : C>>>  superheated vapour
  12 : C     version with wall condensation
  13 : C
  14 : C     QAE:  QA.E  constant air flowrate
  15 : C     QEEi: QEi.E water flowrate
  16 : C     indice 1,2 : entree et sortie de troncon
  17 : C     RECU,XKUL,XKUT1,XKUT2,XKUT3,XKUT4 : coef lois de frot utilisateur
  18 : C
  19 : C  H correspond au flux d'energie total cede a la paroi
  20 : C  HM (convection) est utilise pour le calcul de T2
  21 : C
  22 :       HW=0.D0
  23 :       RINDEX = 1.0
  24 :       DQEE = 0.D0
  25 : 
  25 :       CALL BPHYS(T0,P0,RA,RV,CA,CV,CL,XLAT0,ROL,XKL,XKT,XREL)
  26 : 
  27 : C     interpolation de la temperature de paroi
  27 :       CALL BTPAR(XN,TN,X,TP,NPP,ITP)
  28 :       TP=TP+T0
  29 : C     epaisseur a l'entree
  29 :       CALL BTPAR(XN,EN,X,E1,NPP,ITP)
  30 : C     etendue a l'entree
  30 :       CALL BTPAR(XN,BN,X,B1,NPP,ITP)
  31 : 
  32 :       QA1=QAE/E1/B1
  33 :       QE1=QEE1/E1/B1
  34 :       Q = QA1 + QE1
  35 : C     calcul de R et CP
  35 :       CALL BCAR(RA,RV,PHI1,CV,CA,Q,RQ,R,QA0,QV,CP)
  36 :       IF(KIMP.GE.1) THEN
  37 :         write(6,*) 'entree bsur : E1,B1,PHI1= ',E1,B1,PHI1
  38 :       ENDIF
  39 : 
  40 :       RO=P1/R/T1
  41 :       U1=Q/RO
  42 : 
  43 :       U=U1
  44 :       T=T1
  45 :       P=P1
  46 :       QE=QE1
  47 :       QA=QA1
  48 :       QX=QA+QE
  49 :       PHI=PHI1
  50 : 
  51 :       IF ((X+DX).GT.1.) DX=1-X
  52 :       X=X+DX
  53 : 
  54 :       DU=0
  55 :       T2=0
  56 :       P2=0
  57 : c     ouverture en sortie de maille
  57 :       CALL BTPAR(XN,EN,X,E2,NPP,ITP)
  58 : c     etendue en sortie de maille
  58 :       CALL BTPAR(XN,BN,X,B2,NPP,ITP)
  59 :       QA2=QAE/E2/B2
  60 :       E=E1
  61 :       B=B1
  62 : 
  63 :       NITER=0
  64 : 
  65 :   10  CONTINUE
  66 : 
  67 :       NITER=NITER+1
  68 : 
  69 : c      write(6,*) ' bsur X,NITER,T1,P1,RE',X,NITER,T1,P1,RE
  70 :       TT2=T2
  71 :       PP2=P2
  72 : 
  73 :       XMU=BMUG(T,PHI,T0)
  74 :       RE=QX*2*E/XMU
  75 :       PR=XMU*CP/BLA(T,T0)
  76 :       DPP = XL*QX*QX/4/E
  77 :       BK = BKFRO(RE,REL,XKL,XKT,2*E,RUG,RECU,XKUL,XKUT1,XKUT2,XKUT3,
  78 :      *XKUT4)
  79 :       DP=DX*BK*DPP/RO
  80 : 
  81 :       P2=P1-DP
  82 :          IF (P2.LT.PSLIM) THEN
  83 :             PSQ = -1.
  84 :             RETURN
  85 :          ENDIF
  86 : 
  87 :       S=0.D0
  88 : C*
  89 :       XLA=BLA(T,T0)
  90 : 
  91 :       HM = BHECH(RE,REL,PR,XLA,E)
  92 : 
  93 : C if vapour at the inlet
  94 : 
  95 :       IF(PHI1.GT.1.D-2) THEN
  96 : C
  97 : C wall condensation
  98 :          D=2.5D-5
  99 :          A=2.3D-5
 100 :          PP=BPSAT(TP)
 101 :          PV=PHI1*P1
 102 :          RAPP=(P1-PP)/(P1-PV)
 103 : 
 104 :          IF (RAPP.GT.(1.0001)) THEN
 105 :          XJ=XW*(HM/CV)*(D/A)**(2./3.)*DLOG(RAPP)
 106 :          ELSE
 107 :          XJ=0.D0
 108 :          ENDIF
 109 : 
 110 :          DQEE=-2.*XJ*DX*XL*B
 111 :          XLAT=XLAT0+(CV-CL)*T
 112 :          SW=(XLAT/QX/CP/E)*DQEE
 113 : C**
 114 :          DELTAT = T-TP
 115 :          IF (ABS(DELTAT).GT.1D-3) THEN
 116 :          HW=XJ*XLAT/(T-TP)
 117 : C        write(6,*) ' T TP HW ',T,TP,HW
 118 :          ENDIF
 119 : 
 120 :          QEE2=QEE1+DQEE
 121 :          QE2=QEE2/E2/B2
 122 : 
 123 :          IF ((QEE2/QEE1).LT.(-1.D-5)) THEN
 124 :            write(6,*) 'bsur : relative water flowrate ',QE2/QE1
 125 :            RINDEX = 0.5
 126 :            RETURN
 127 :          ELSE
 128 :            A=RA*QA2/RV/QE2
 129 :            PHI2=1./(1.+A)
 130 :          ENDIF
 131 : * no vapour at the inlet
 132 :       ELSE
 133 :          XLAT=XLAT0+(CV-CL)*T
 134 :          HW=0.D0
 135 :          PHI2=0.D0
 136 :          QEE2=0.D0
 137 :          QE2=0.D0
 138 :       ENDIF
 139 :       IF(KIMP.GE.1) THEN
 140 :          write(6,*) 'bsur: NITER,PHI2,QA2,QE2= ',NITER,PHI2,QA2,QE2
 141 :       ENDIF
 142 :       QEEX=(QEE1+QEE2)/2.
 143 : 
 144 :       T2 = BT2(T1,TP,S,HM,PR,RE,XLA,XL,DX)
 145 :       IF(KIMP.GE.1) WRITE(*,*) 'bsur : T2= ',T2
 146 : 
 147 :       U2 = (QA2+QE2)*R*T2/P2
 148 :       RO=(P1+P2)/R/(T1+T2)
 149 :       T=(T1+T2)/2
 150 :       P=(P1+P2)/2
 151 :       PHI=(PHI1+PHI2)/2
 152 : 
 153 :       E=(E1+E2)/2
 154 :       B=(B1+B2)/2
 155 :       QA=QAE/E/B
 156 :       QE=QEEX/E/B
 157 :       QX=QA+QE
 158 : 
 159 :       QW2=(QEE2-QEE1)/DX/XL/B
 160 : 
 161 :       ERT=DABS((T2-TT2)/T2)
 162 :       ERP=DABS((P2-PP2)/P2)
 163 : 
 164 :       IF(((ERT.GT.1E-4).OR.(ERP.GT.1E-4))
 165 :      $  .AND.(NITER.LE.10))    GOTO 10
 166 : 
 167 :          IF(KIMP.GE.2) THEN
 168 :          write(6,*) ' HM HW HW/HM ',HM,HW,HW/HM
 169 :          ENDIF
 170 : 
 171 : C     calcul des puissances echangees
 172 : 
 173 :       RAPH=(CV*T+XLAT0)/XLAT
 174 :       H=HM+HW*RAPH
 175 : 
 176 : C     WRITE(6,*) ' bsur ',H
 177 : 
 178 :       HH1=(QAE*CA*T1)+(QEE1*(CV*T1+XLAT0))
 179 :       HH2=(QAE*CA*T2)+(QEE2*(CV*T2+XLAT0))
 180 : C !   PF = HH2-HH1
 181 :       PF = (HH2-HH1) - DQEE*(CV*T1+XLAT0)
 182 :       PP = H*(T-TP)*2*DX*XL*B
 183 :       IF(KIMP.EQ.-3) THEN
 184 :         write(6,*) ' DQEE CV T1 XLAT0 ',DQEE,CV,T1,XLAT0
 185 :         write(6,1000)  X,RAPH,PF,PP,HH1,HH2
 186 : 1000    FORMAT(1X,' bsur X RH DH EC HH1 HH2 ',6E12.5)
 187 :       ENDIF
 188 :       IF(KIMP.NE.0.AND.X.EQ.1.) THEN
 189 :         write(6,2110) X,RE,BK
 190 : 2110    format(1X,'bsur X RE BKRO ',3E12.5)
 191 :       ENDIF
 192 :       DPF=PF/DX/XL/B
 193 :       DPP=PP/DX/XL/B
 194 : 
 195 :       RETURN
 196 :       END
 197 :  
 198 :  
 199 :  
 200 :  
 201 :  
 202 :  
 203 :