Télécharger bcond2.eso

Retour à la liste

Numérotation des lignes :

  1. C BCOND2 SOURCE CB215821 16/04/21 21:15:17 8920
  2. SUBROUTINE BCOND2 (X,DX,XL,RUG,XW,XN,TN,EN,BN,KIMP,PSLIM,REL,
  3. & P1,PS1,T1,Y1,QAE,QEE1,PHI1,
  4. & P2,T2,Y2,U2,QEE2,PHI2,QW2,RE,H,PSQ,RINDEX,
  5. & NPP,ITP,PF,PP,DPF,DPP,RECU,XKUL,XKUT1,XKUT2,XKUT3,XKUT4)
  6.  
  7. IMPLICIT INTEGER(I-N)
  8. IMPLICIT REAL*8 (A-H,O-Z)
  9. C operateur FUITE option gaz reel
  10. C>>> air + vapeur + liquide
  11. C QAE : QA*E*B invariant air flowrate
  12. C QEEi : QEEi*E*B invariant water flowrate
  13. C XW=1.(0.) with (without) wall condensation
  14. C indice 1,2 : entree et sortie de troncon
  15. C RECU,XKUL,XKUT1,XKUT2,XKUT3,XKUT4 : coef lois de frot utilisateur
  16. C H correspond a l'energie totale cedee a la paroi
  17. C H=HM+RAPH*HW transmis
  18. C H=HM+HW pour le calcul de T2
  19. C
  20. DIMENSION XN(NPP),TN(NPP),EN(NPP),BN(NPP)
  21.  
  22. IF(KIMP.GE.1) THEN
  23. write(6,*) 'entree bcond2 X T1 Y1 PHI1= ',X,T1,Y1,PHI1
  24. ENDIF
  25. RINDEX=1.0
  26.  
  27. ERPS1 = ABS (1.-(PS1/PSATT0(T1)))
  28. IF (ERPS1.GT.1D-4) THEN
  29. WRITE(6,*) ' bcond :erreur sur PS1 ',ERPS1
  30. ENDIF
  31.  
  32. C IMPLEMENTATION DES PROPRIETES PHYSIQUES
  33. CALL BPHYS(T0,P0,RA,RV,CA,XCV,XCL,XLAT0,XROL,XKL,XKT,XREL)
  34.  
  35. C TEMPERATURE PAROI A ENTREE (KELVIN)
  36. CALL BTPAR(XN,TN,X,TP,NPP,ITP)
  37. TP=TP+T0
  38.  
  39. C EPAISSEUR DE FISSURE EN ENTREE
  40. CALL BTPAR(XN,EN,X,E1,NPP,ITP)
  41. C etendue a l'entree
  42. CALL BTPAR(XN,BN,X,B1,NPP,ITP)
  43.  
  44. C-------------------------------------
  45. C DEFINITION DES DIFFERENTS DEBITS (homogenes a ro.u)
  46. C-------------------------------------
  47. C debit air
  48. QA=QAE/E1/B1
  49.  
  50. C debit liquide
  51. QE1=QEE1/E1/B1
  52.  
  53. C debit vapeur
  54. QV1=QE1*Y1
  55.  
  56. C debit liquide
  57. QL1= (1-Y1)*QE1
  58.  
  59. C debit total
  60. Q = QA + QE1
  61.  
  62. C-------------------------------------
  63. C DEFINITION DES MASSES VOLUMIQUES
  64. C CONSTITUANT:
  65. C : air
  66. C : vapeur
  67. C : eau
  68. C : eau + vapeur
  69. C : air + eau + vapeur
  70. C--------------------------------------
  71. C PS1=PHI1*P1
  72.  
  73. C masse volumique de air
  74. ROA1=(P1-PS1)/RA/T1
  75. C masse volumique vapeur
  76. ROV1=ROVAP0(PS1,T1)
  77.  
  78. C masse volumique eau liquide
  79. ROL1=ROLIQ0(P1,T1)
  80.  
  81. C calcul la fraction volumique de gaz
  82. AL =BALPHA(Y1,ROL1,ROV1)
  83.  
  84. C calcul la masse volumique air + eau + vapeur
  85. RO1=AL *(ROA1+ROV1)+(1.-AL )*ROL1
  86. C-------------------------------------------------
  87. C CALCUL DE LA VITESSE DE ECOULEMENT EN ENTREE
  88. C-------------------------------------------------
  89.  
  90. U1=QA/AL /ROA1
  91.  
  92. C-------------------------------------------------
  93. C CALCUL DU COEFFICIENT DE COMPRESSIBILITE Z
  94. C CONSTITUANT:
  95. C : vapeur
  96. C-------------------------------------------------
  97. ZV = ZVAP0(ROV1,T1)
  98.  
  99. C-------------------------------------------------
  100. C CALCUL DES VISCOSITES DYNAMIQUES
  101. C CONSTITUANT:
  102. C : air + vapeur
  103. C : air + eau + vapeur
  104. C-------------------------------------------------
  105.  
  106. C viscosite dynamique air+vapeur
  107. XMU=BMUG(T1,PS1/P1,T0)
  108.  
  109. C calcul de la viscosite dynamique air+vapeur+eau_liquide
  110. XMU=AL*XMU+BMUL(T1,T0)*(1.-AL)
  111. C WRITE(6,*) 'XMU',XMU
  112.  
  113. C--------------------------------------------------
  114. C CALCUL DES ENTHALPIES SPECIFIQUES
  115. C CONSTITUANT:
  116. C : eau liquide
  117. C : vapeur
  118. C--------------------------------------------------------
  119.  
  120. C calcul enthalpie specifique eau liquide
  121. HLS01 = HLS0(P1,T1)
  122. C WRITE(6,*) 'HL',HLS01
  123. C calcul enthalpie specifique de la vapeur
  124. HVS01 = HVS0(PS1,T1)
  125. C WRITE(6,*) 'HV',HVS01
  126.  
  127. C--------------------------------------------------------
  128. C CALCUL CHALEUR LATENTE DU MELANGE
  129. C CONSTITUANT:
  130. C : eau+vapeur
  131. C--------------------------------------------------------
  132.  
  133. C calcul la chaleur latente du melange
  134. XLAT = (HVS01-HLS01)
  135. C WRITE(6,*) 'XLAT',XLAT
  136. C--------------------------------------------------------
  137. C CALCUL DES CONDUCTIVITES THERMIQUES
  138. C CONSTITUANT:
  139. C : air
  140. C--------------------------------------------------------
  141.  
  142. XLA = BLA(T1,T0)
  143.  
  144. C--------------------------------------------------------
  145. C CALCUL DES CHALEURS SPECIFIQUES
  146. C CONSTITUANT:
  147. C : vapeur
  148. C : eau liquide
  149. C : air + vapeur + eau
  150. C--------------------------------------------------------
  151.  
  152. C calcul chaleur specifique de la vapeur
  153. CV1 = DHVDT0(PS1,T1)
  154.  
  155. C calcul chaleur specifique eau liquide
  156. CL1 = CPLIQ0(P1,HLS01)
  157.  
  158.  
  159. C calcul chaleur specifique du melange air vapeur liquide
  160. CP1 = (QA*CA+QE1*(Y1*CV1+(1-Y1)*CL1))/(QA+QE1)
  161.  
  162.  
  163. C WRITE(6,*) 'cv1,cl1,cp1',CV1,CL1,CP1
  164. C--------------------------------------------------------
  165. C ACTUALISATION DES VARIABLES
  166. C CONSTITUANT:
  167. C : vitesse U
  168. C : titre Y
  169. C : densite du melange RO
  170. C : pression vapeur PS
  171. C : pression totale P
  172. C : temperature T
  173. C : debit eau (vapeu+liquide) QE
  174. C : debit total QX
  175. C------------------------------------------------------
  176.  
  177. E=E1
  178. B=B1
  179. U=U1
  180. Y=Y1
  181. RO=RO1
  182. PS=PS1
  183. P=P1
  184. T=T1
  185. QE=QE1
  186. QX=QA+QE
  187.  
  188.  
  189.  
  190. IF ((X+DX).GT.1.D0) DX=1.D0-X
  191.  
  192. X=X+DX
  193. DU=0
  194. T2=0
  195. P2=0
  196.  
  197. CALL BTPAR(XN,EN,X,E2,NPP,ITP)
  198. CALL BTPAR(XN,BN,X,B2,NPP,ITP)
  199.  
  200. NITER = 0
  201. 10 CONTINUE
  202.  
  203. NITER = NITER + 1
  204. C write(6,*)
  205. C write(6,*) ' niter P2 T2 ',NITER,P2,T2
  206. TT2=T2
  207. PP2=P2
  208. YY2=Y2
  209.  
  210.  
  211. C--------------------------------------------------------
  212. C CALCUL DU NOMBRE DE REYNOLDS DE ECOULEMENT
  213. C CONSTITUANT:
  214. C : air + eau + vapeur
  215. C--------------------------------------------------------
  216.  
  217. C calcul du nombre de REYNOLDS
  218. RE=QX*2*E/XMU
  219. C--------------------------------------------------------
  220. C CALCUL DE LA PERTE CHARGE DU AU FORCE DE FROTTEMENT
  221. C CONSTITUANT:
  222. C : air + eau + vapeur
  223. C--------------------------------------------------------
  224.  
  225. DPP = XL*QX*QX /4/E
  226.  
  227. C calcul de la perte de charge DP
  228. BK=BKFRO(RE,REL,XKL,XKT,2*E,RUG,RECU,XKUL,XKUT1,XKUT2,XKUT3,
  229. * XKUT4)
  230. DP=DX*BK*DPP/RO
  231. IF(KIMP.GE.1) THEN
  232. WRITE(6,*) 'bcond2 ap BKFRO : RE,QX,DPP= ',RE,QX,DPP
  233. WRITE(6,*) 'bcond2 ap BKFRO : DX,BK,RO= ',DX,BK,RO
  234. ENDIF
  235. P2=P1-DP
  236. IF (P2.LT.PSLIM) THEN
  237. PSQ = -1.
  238. RETURN
  239. ENDIF
  240.  
  241.  
  242. C--------------------------------------------------------
  243. C CALCUL DU NOMBRE DE PRANDT DE ECOULEMENT
  244. C CONSTITUANT:
  245. C : air + eau + vapeur
  246. C--------------------------------------------------------
  247.  
  248. C calcul du nombre de Prandt estime en amont
  249. PR = XMU*CP1/BLA(T1,T0)
  250.  
  251.  
  252.  
  253. C--------------------------------------------------------
  254. C CALCUL DU COEFFICIENT ECHANGE THERMIQUE
  255. C CONSTITUANT:
  256. C : air + eau + vapeur
  257. C--------------------------------------------------------
  258.  
  259.  
  260. C calcul du coefficient echange thermique
  261. HM = BHECH(RE,REL,PR,XLA,E)
  262. H = HM
  263. HW=0.D0
  264.  
  265. C si condensation en film
  266. IF(T1.GE.TP) THEN
  267. HW=XW*BWALL(XLAT,CV1,PS/P,H,(DPSAT0(T1)/PSATT0(T1)))
  268. H = HM+HW
  269. ENDIF
  270.  
  271.  
  272. C calcul de la fraction molaire de vapeur
  273. PHI=PS/P
  274.  
  275. C--------------------------------------------------------
  276. C CALCUL DES COEFFICIENTS ALPHAi
  277. C CONSTITUANT:
  278. C : on se sert des fonctions de VARI
  279. C--------------------------------------------------------
  280.  
  281.  
  282. C definition des alphai pour resolution systeme
  283.  
  284. alpha1=DHVDT0(PS ,T )
  285. alpha2=DHVDP0(PS ,T )
  286. alpha3=DHLDP0(P ,T )
  287. HLIQ=HLS0(P ,T )
  288. alpha4=CPLIQ0(P ,HLIQ )
  289. alpha5=DPSAT0(T )
  290. alpha6=DZVDP0(PS ,T )
  291. alpha7=DZVDT0(PS ,T )
  292.  
  293.  
  294. C--------------------------------------------------------
  295. C INTEGRATION DES PARAMETRES DE EQUATION THERMIQUE
  296. C Qy dT/dx= (H)(Tw-T)+S resolution BT3
  297. C--------------------------------------------------------
  298.  
  299.  
  300. C definition de Qc
  301. Qc=QAE*CA+Y *QEE1*(alpha1+(alpha2*alpha5))+(1-Y )*alpha4*QEE1
  302. C WRITE(6,*) 'QA QE1 Qc ',QA,QE1,Qc
  303.  
  304.  
  305. C definition de a
  306. a=(RA/RV)*QAE
  307. C write(6,*) 'a ',a
  308.  
  309.  
  310. C DEFINITION DE ALPHA9
  311. x1=(PHI *alpha5)/(ZV*PSATT0(T )*(1-PHI )*(1-PHI ))
  312. x2=(PHI *(alpha5*alpha6+alpha7))/(1-PHI )
  313. x2=x2/ZV/ZV
  314. alpha9=x1-x2
  315.  
  316. C WRITE(6,*) 'ALPHA9 x1 x2 ',alpha9 ,x1,x2
  317.  
  318. C DEFINITION DE ALPHA10
  319. x3=(PHI *XLAT*a)/(ZV*P *(1-PHI )*(1-PHI ))
  320. x4=alpha3*(1-Y )*QEE1
  321. alpha10=x3-x4
  322. C WRITE(6,*) 'alpha10 x3 x4 ',alpha10,x3,x4
  323.  
  324.  
  325. C calcul du terme source de la thermique S
  326. C rappel DP=P1-P2 et non l'inverse
  327. S=-(alpha10)*(DP/(XL*DX))
  328.  
  329. Qy=Qc+a*XLAT*alpha9
  330.  
  331.  
  332. C----------------------------------------------------
  333. C Calcul la temperature T2 de sortie du troncon
  334. C----------------------------------------------------
  335.  
  336. HEFF = 2 * H
  337. C bt3 traite l echange de chaleur par unite d etendue
  338. Qy = Qy/B
  339. S = S/B
  340. T2 = BT3(T1,TP,S,HEFF,Qy,DX,XL)
  341.  
  342. IF(KIMP.GE.1) THEN
  343. WRITE(*,*) 'bcond2 ap BT3: niter Qy B HEFF S',niter,Qy,B,HEFF,S
  344. WRITE(6,*) 'bcond2 ap BT3: T1 T2 TP= ',T1,T2,TP
  345. ENDIF
  346.  
  347. C calcul de increment de temperature DT
  348. DT = T2-T1
  349.  
  350. C----------------------------------------------------
  351. C FIN DU TRONCON CALCULS DES NOUVELLES VALEURS
  352. C----------------------------------------------------
  353.  
  354. C calcul de nouvelle pression vapeur saturante
  355. PS2=PSATT0(T2)
  356.  
  357. C--------------------------------------------------------
  358. C DEFINITION DES NOUVELLES MASSES VOLUMIQUES
  359. C CONSTITUANT:
  360. C : air
  361. C : vapeur
  362. C : eau
  363. C : eau + vapeur
  364. C : air + eau + vapeur
  365. C--------------------------------------------------------
  366.  
  367. C masse volumique de air
  368. ROA2=(P2-PS2)/RA/T2
  369.  
  370. C masse volumique vapeur
  371. ROV2=ROVAP0(PS2,T2)
  372.  
  373. C--------------------------------------------------------
  374. C CALCUL DU COEFFICIENT DE COMPRESSIBILITE Z de vapeur
  375. C CONSTITUANT:
  376. C : vapeur
  377. C--------------------------------------------------------
  378. Z2 = ZVAP0(ROV2,T2)
  379.  
  380. C masse volumique eau liquide
  381. ROL2=ROLIQ0(P2,T2)
  382.  
  383.  
  384. C--------------------------------------------------------
  385. C Calcul du dQe
  386. C--------------------------------------------------------
  387.  
  388. C calcul increment de debit eau
  389. DQEE= -2.*HW*DX*XL*B*(T-TP)/XLAT
  390.  
  391. C calcul nouveau debit eau (liq+vap)
  392. QEE2=QEE1+DQEE
  393.  
  394. C calcul debit eau
  395. QE2=QEE2/E2/B2
  396.  
  397. C calcul du titre de vapeur sortie troncon
  398. Y2=(ROV2/ROA2)*(QA/QE2)
  399.  
  400. IF(KIMP.GE.3) THEN
  401. WRITE(6,*) ' ROA2 ROV2 Y2',ROA2,ROV2,Y2
  402. write(6,*) ' QE1 QE2 ',QE1,QE2
  403. ENDIF
  404.  
  405.  
  406. C--------------------------------------------------------
  407. C On raffine le maillage dans 2 cas
  408. C--------------------------------------------------------
  409.  
  410. IF ( ((QEE2/QEE1).LT.1.E-5).OR.(Y2.GE.1.)) THEN
  411. C WRITE(6,*) ' WALL FLUX ',QEE2/QEE1
  412. RINDEX=0.5
  413. RETURN
  414. ENDIF
  415.  
  416.  
  417. QV2 = QE2 * Y2
  418. QL2 = QE2 * (1.-Y2)
  419.  
  420.  
  421. C calcul la fraction volumique de gaz
  422. AL2=BALPHA(Y2,ROL2,ROV2)
  423.  
  424.  
  425. C calcul la masse volumique air + eau + vapeur
  426. RO2=AL2*(ROA2+ROV2)+(1.-AL2)*ROL2
  427.  
  428.  
  429.  
  430.  
  431.  
  432. C--------------------------------------------------------
  433. C CALCUL DE LA NOUVELLE VITESSE DE ECOULEMENT
  434. C--------------------------------------------------------
  435.  
  436. U2= QA/AL2/ROA2
  437.  
  438. C--------------------------------------------------------
  439. C DEFINITION DES NOUVELLES VALEURS MOYENNES
  440. C
  441. C : Pression
  442. C : Temperature
  443. C : pression vapeur
  444. C--------------------------------------------------------
  445.  
  446. P=(P1+P2)/2
  447. T=(T1+T2)/2
  448. PS=(PS1+PS2)/2
  449.  
  450. PHI = PS/P
  451.  
  452. C--------------------------------------------------------
  453. C CALCUL DES VITESSES MOYENNES SUR TRONCON
  454. C CONSTITUANT: melange
  455. C : masse volumique
  456. C : vitesse
  457. C : debit eau
  458. C : debit total
  459. C : titre
  460. C--------------------------------------------------------
  461.  
  462.  
  463. C calcul de la masse volumique moyenne melange
  464. RO=(RO1+RO2)/2
  465.  
  466. C calcul de la vitesse moyenne du melange
  467. U=Q/RO
  468.  
  469.  
  470. C calcul du debit eau moyen du melange
  471. QE=(QE1+QE2)/2.
  472.  
  473. C calcul du debit total moyen du melange
  474. QX=QA+QE
  475.  
  476. C calcul du titre moyen sur le troncon
  477. Y=(Y1+Y2)/2
  478.  
  479. C-----------------------------------------------------------------
  480. C CALCUL DU TAUX DE CONDENSATION EN FILM
  481. C-----------------------------------------------------------------
  482.  
  483. QW2=(QEE2-QEE1)/DX/XL/B
  484.  
  485.  
  486.  
  487. IF(KIMP.GE.3) write(6,*) ' P2 T2 Y2 ',P2,T2,Y2
  488.  
  489. C definition des conditions fin de boucle
  490.  
  491.  
  492. C critere sur P,T,Y
  493.  
  494. ERT=ABS((T2-TT2)/T2)
  495. ERP=ABS((P2-PP2)/P2)
  496. ERY=ABS((Y2-YY2)/Y2)
  497.  
  498. IF (((ERT.GT.1E-4).OR.(ERP.GT.1E-4).OR.(ERY.GT.1E-4))
  499. $ .AND.(NITER.LE.10)) GOTO 10
  500.  
  501. IF(KIMP.GE.2) THEN
  502. write(6,2300) X, P2/P0, T2-T0, Y2, RO2, QE1, QE2
  503. 2300 FORMAT((1X,'bcond2: X P2 T2 Y2 RO2 QE1 QE2 ',7E12.5))
  504. ENDIF
  505. IF(KIMP.GE.3) THEN
  506. WRITE(6,*) ' HM HW HW/HM ',HM,HW,HW/HM
  507. ENDIF
  508. PHI2=PHI1
  509.  
  510.  
  511.  
  512. C--------------------------------------------------------
  513. C BILAN ENERGIE SUR TRONCON
  514. C--------------------------------------------------------
  515.  
  516. C write(6,*) ' E1 E2 ',E1,E2
  517.  
  518. HH1=(QA*CA*T1)+(QV1*HVS0(PS1,T1))+(QL1*HLS0(P1,T1))
  519. HH1=HH1*E1*B1
  520. HH2=(QA*CA*T2)+(QV2*HVS0(PS2,T2))+(QL2*HLS0(P2,T2))
  521. HH2=HH2*E2*B2
  522. C puissance fluide
  523. C ! PF = HH2-HH1
  524. PF = (HH2-HH1) - (DQEE*HVS0(PS1,T1))
  525.  
  526. C puissance paroi
  527. PP = H *(T-TP)*2*DX*XL*B
  528.  
  529.  
  530. IF(KIMP.GE.2) THEN
  531. write(6,2200) X,HH1,HH2,PF,PP
  532. 2200 format(1X,'bcond2 X HH1 HH2 DH HDT ',5E12.5)
  533. ENDIF
  534. IF(KIMP.NE.0.AND.X.EQ.1.) THEN
  535. write(6,2110) X,RE,BK
  536. 2110 format(1X,'bcond2 X RE BKRO ',3E12.5)
  537. ENDIF
  538. DPF=PF/DX/XL/B
  539. DPP=PP/DX/XL/B
  540.  
  541. RETURN
  542. END
  543.  
  544.  
  545.  
  546.  
  547.  
  548.  
  549.  
  550.  
  551.  
  552.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales