Télécharger bsur2.eso

Retour à la liste

Numérotation des lignes :

  1. C BSUR2 SOURCE CB215821 16/04/21 21:15:25 8920
  2. SUBROUTINE BSUR2 (X,DX,XL,RUG,XW,XN,TN,EN,BN,KIMP,PSLIM,REL,
  3. & RINDEX,P1,T1,QAE,QEE1,PHI1,P2,T2,U2,QEE2,PHI2,QW2,RE,H,
  4. & PSQ,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 version gaz reel
  14. C QAE: QA.E.B constant air flowrate
  15. C QEEi: QEi.E.B 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 T
  21. C
  22. IF(KIMP.GE.1) THEN
  23. write(6,*) 'entree bsur2'
  24. ENDIF
  25.  
  26. HW=0.D0
  27. RINDEX=1.0
  28. DQEE = 0.D0
  29.  
  30. C IMPLEMENTATION DES PROPRIETES PHYSIQUES
  31. CALL BPHYS(T0,P0,RA,RV,CA,XCV,XCL,XLAT0,XROL,XKL,XKT,XREL)
  32.  
  33. C TEMPERATURE PAROI A ENTREE (KELVIN)
  34. CALL BTPAR(XN,TN,X,TP,NPP,ITP)
  35. TP=TP+T0
  36.  
  37. C EPAISSEUR DE FISSURE EN ENTREE
  38. CALL BTPAR(XN,EN,X,E1,NPP,ITP)
  39. C etendue a l'entree
  40. CALL BTPAR(XN,BN,X,B1,NPP,ITP)
  41.  
  42. C--------------------------------------------------------
  43. C DEFINITION DES DIFFERENTS DEBITS
  44. C--------------------------------------------------------
  45. C debit air
  46. QA=QAE/E1/B1
  47.  
  48. C debit eau (vapeur+liquide)
  49. QE1=QEE1/E1/B1
  50.  
  51. C debit total
  52. Q=QA+QE1
  53. C
  54. C DEFINITION DES MASSES VOLUMIQUES
  55. C CONSTITUANT:
  56. C : air
  57. C : vapeur
  58.  
  59. PV1=PHI1*P1
  60.  
  61. C masse volumique de air
  62. ROA1=(P1-PV1)/RA/T1
  63. C WRITE(6,*) 'T1,P1,PV1,RA',T1,P1,PV1,RA
  64.  
  65. C masse volumique vapeur
  66. ROV1=ROVAP0(PV1,T1)
  67.  
  68. C masse volumique du melange
  69. RO1 = ROA1 + ROV1
  70. C write(6,*) ' ROA ROV RO ',ROA1,ROV1,RO1
  71.  
  72. C--------------------------------------------------------
  73. C CALCUL DE LA VITESSE DE ECOULEMENT EN ENTREE
  74. C--------------------------------------------------------
  75.  
  76. C calcul de la vitesse
  77. U1=QA/ROA1
  78. C write(6,*) ' U1 ',U1
  79.  
  80.  
  81. C--------------------------------------------------------
  82. C CALCUL DU COEFFICIENT DE COMPRESSIBILITE Z de vapeur
  83. C CONSTITUANT:
  84. C : vapeur
  85. C--------------------------------------------------------
  86. C calcul du ceofficient de compressibilite Z vapeur
  87. ZV = ZVAP0(ROV1,T1)
  88. C WRITE(6,*) 'Z',ZV
  89.  
  90.  
  91. C--------------------------------------------------------
  92. C CALCUL DES CONDUCTIVITES THERMIQUES
  93. C CONSTITUANT:
  94. C : air
  95. C--------------------------------------------------------
  96.  
  97. C calcul la conductivite de air
  98. XLA = BLA(T1,T0)
  99.  
  100.  
  101.  
  102. C--------------------------------------------------------
  103. C CALCUL DES CHALEURS SPECIFIQUES
  104. C CONSTITUANT:
  105. C : vapeur
  106. C----------------------------------------------------
  107.  
  108. CV1 = DHVDT0(PV1,T1)
  109. C WRITE(6,*) 'cv1',CV1
  110.  
  111.  
  112.  
  113. C--------------------------------------------------------
  114. C CALCUL DES VISCOSITES DYNAMIQUES
  115. C CONSTITUANT:
  116. C : air + vapeur
  117. C--------------------------------------------------------
  118.  
  119. XMU=BMUG(T1,PHI1,T0)
  120. C WRITE(6,*) 'xmu',XMU
  121.  
  122.  
  123. C--------------------------------------------------------
  124. C CALCUL DU NOMBRE DE PRANDT DE ECOULEMENT
  125. C CONSTITUANT:
  126. C : air + vapeur
  127. C--------------------------------------------------------
  128.  
  129. PR=XMU*CV1/BLA(T1,T0)
  130. C WRITE(6,*) 'PR',PR
  131. C--------------------------------------------------------
  132. C CALCUL DES ENTHALPIES SPECIFIQUES
  133. C CONSTITUANT:
  134. C : eau liquide
  135. C : vapeur
  136. C--------------------------------------------------------
  137. HVS01 = HVS0(PV1,T1)
  138. HLS01 = HLS0(P1,T1)
  139. C WRITE(6,*) 'HV',HVS01
  140. C WRITE(6,*) 'HL',HLS01
  141.  
  142. C--------------------------------------------------------
  143. C CALCUL CHALEUR LATENTE DU MELANGE
  144. C (utile pour la condensation en film)
  145. C CONSTITUANT:
  146. C : eau+vapeur
  147. C--------------------------------------------------------
  148.  
  149. XLAT = (HVS01-HLS01)
  150. C WRITE(6,*) 'L',XLAT
  151.  
  152. C--------------------------------------------------------
  153. C ACTUALISATION DES VARIABLES
  154. C CONSTITUANT:
  155. C : vitesse U
  156. C : pression totale P
  157. C : temperature T
  158. C----------------------------------------------------
  159.  
  160. E=E1
  161. B=B1
  162. U=U1
  163. T=T1
  164. P=P1
  165. RO=RO1
  166.  
  167.  
  168. IF ((X+DX).GT.1.) DX=1-X
  169. X=X+DX
  170.  
  171. DU=0
  172. T2=0
  173. P2=0
  174.  
  175. CALL BTPAR(XN,EN,X,E2,NPP,ITP)
  176. CALL BTPAR(XN,BN,X,B2,NPP,ITP)
  177.  
  178. NITER=0
  179.  
  180. 10 CONTINUE
  181.  
  182. NITER=NITER+1
  183. IF(KIMP.GE.1) THEN
  184. write(*,*) 'NITER= ',NITER
  185. ENDIF
  186.  
  187. TT2=T2
  188. PP2=P2
  189.  
  190. C--------------------------------------------------------
  191. C CALCUL DU NOMBRE DE REYNOLDS
  192. C CONSTITUANT:
  193. C : melange
  194. C-------------------------------------------------------
  195. RE=Q*2*E/XMU
  196. IF(KIMP.GE.1) THEN
  197. WRITE(6,*) 'bsur2: RE',RE
  198. WRITE(6,*) 'bsur2: Q,E,B,XMU',Q,E,B,XMU
  199. ENDIF
  200. C--------------------------------------------------------
  201. C CALCUL DE LA PERTE DE CHARGE DP DUE AU FROTTEMENT
  202. C CONSTITUANT:
  203. C : Melange
  204. C--------------------------------------------------------
  205.  
  206. DPP = XL*Q*Q/4/E
  207. BK = BKFRO(RE,REL,XKL,XKT,2*E,RUG,RECU,XKUL,XKUT1,XKUT2,XKUT3,
  208. *XKUT4)
  209. DP=DX*BK*DPP/RO
  210. C WRITE(6,*) 'XL Q E ', XL,Q,E
  211. C WRITE(6,*) 'DX RO DP', DX,RO,DP
  212.  
  213. P2=P1-DP
  214. IF (P2.LT.PSLIM) THEN
  215. PSQ = -1.
  216. RETURN
  217. ENDIF
  218.  
  219.  
  220. C--------------------------------------------------------
  221. C CALCUL DU TERME SOURCE S POUR LE CALCUL DE T
  222. C CONSTITUANT:
  223. C : melange
  224. C----------------------------------------------------
  225.  
  226. S=0.D0
  227.  
  228. C--------------------------------------------------------
  229. C CALCUL DU COEFFICIENT ECHANGE THERMIQUE H
  230. C CONSTITUANT:
  231. C : Melange
  232. C----------------------------------------------------
  233.  
  234. HM = BHECH(RE,REL,PR,XLA,E)
  235. C write(6,*) RE,REL,PR
  236. C write(6,*) XLA,E
  237. C write(6,*) ' HM ',HM
  238.  
  239. C if vapour at the inlet
  240.  
  241. IF(PHI1.GT.1.D-2) THEN
  242. C
  243. C wall condensation
  244. D=2.5D-5
  245. A=2.3D-5
  246.  
  247. C--------------------------------------------------------
  248. C CALCUL PRESSION DE SATURATION
  249. C CONSTITUANT:
  250. C : vapeur
  251. C----------------------------------------------------
  252.  
  253. PP=PSATT0(TP)
  254. C WRITE(6,*) 'TP',TP
  255. C--------------------------------------------------------
  256. C CALCUL PRESSION PARTIELLE DE VAPEUR
  257. C CONSTITUANT:
  258. C : vapeur
  259. C----------------------------------------------------
  260. C WRITE(6,*) 'ROV1',ROV1
  261. PV=PVAP0(ROV1,T1)
  262. C WRITE(6,*) 'PV',PV
  263.  
  264.  
  265. C--------------------------------------------------------
  266. C CONDENSATION EN FILM
  267. C--------------------------------------------------------
  268.  
  269.  
  270. RAPP=(P1-PP)/(P1-PV)
  271.  
  272. IF (RAPP.GT.(1.0001)) THEN
  273. XJ=XW*(HM/CV1)*(D/A)**(2./3.)*LOG(RAPP)
  274. ELSE
  275. XJ=0.D0
  276. ENDIF
  277.  
  278. DQEE=-2.*XJ*DX*XL*B
  279.  
  280. SW=(XLAT/Q/CV1/E)*DQEE
  281.  
  282. DELTAT = T-TP
  283. IF (ABS(DELTAT).GT.1D-3) THEN
  284. HW=XJ*XLAT/(T-TP)
  285. C write(6,*) ' T TP HW ',T,TP,HW
  286. ENDIF
  287.  
  288. QEE2=QEE1+DQEE
  289. QE2=QEE2/E2/B2
  290. C write(6,*) ' QE1 QE2 ',QE1,QE2
  291.  
  292. IF ((QEE2/QEE1).LT.(-1.D-5)) THEN
  293. C write(6,*) ' relative water flowrate ',QE2/QE1
  294. RINDEX = 0.5
  295. RETURN
  296. C PHI2=0.D0
  297. C QEE2=0.D0
  298. C QE2=0.D0
  299. ELSE
  300. A=RA*QA/RV/(QE2*ZV)
  301. PHI2=1./(1.+A)
  302. ENDIF
  303.  
  304.  
  305. * no vapour at the inlet
  306. ELSE
  307. XLAT=(HVS01-HLS01)
  308. HW=0.D0
  309. PHI2=0.D0
  310. QEE2=0.D0
  311. QE2=0.D0
  312. ENDIF
  313.  
  314.  
  315. C write(6,*) ' QE ',QE1,QE2
  316. C--------------------------------------------------------
  317. C CALCUL DES COEFFICIENTS ALPHAi
  318. C CONSTITUANT:
  319. C : on se sert des fonctions de VARI
  320. C--------------------------------------------------------
  321.  
  322.  
  323. C definition des alphai pour resolution systeme
  324. alpha1=DHVDT0(PV1,T1)
  325. alpha2=DHVDP0(PV1,T1)
  326.  
  327. alpha6=DZVDP0(PV1,T1)
  328. alpha7=DZVDT0(PV1,T1)
  329. alpha8=DRVDP0(PV1,T1)
  330. alpha9=DRVDT0(PV1,T1)
  331.  
  332. C WRITE(6,*) 'ALPHA',alpha1
  333. C WRITE(6,*) 'ALPHA',alpha2
  334. C WRITE(6,*) 'ALPHA',alpha3
  335. C WRITE(6,*) 'ALPHA',alpha4
  336. C WRITE(6,*) 'ALPHA',alpha5
  337. C WRITE(6,*) 'ALPHA',alpha6
  338. C WRITE(6,*) 'ALPHA',alpha7
  339. C WRITE(6,*) 'ALPHA',alpha8
  340. C WRITE(6,*) 'ALPHA',alpha9
  341.  
  342. C--------------------------------------------------------
  343. C DEFINITION DU Qxx,Qy
  344. C CONSTITUANT:
  345. C : on se sert des fonctions de VARI
  346. C--------------------------------------------------------
  347.  
  348. Qxx= QAE*CA+QEE1*alpha1
  349. s1=ZV*((alpha9*T1)+ROV1*T1*alpha7)
  350. s2= RA*T1*(alpha8*ZV+ROV1*alpha6)
  351. s3= 1-s2
  352. s4=s1/s3
  353. Qy= Qxx + QEE1*alpha2*RA*s4
  354.  
  355. C--------------------------------------------------------
  356. C CALCUL DE LA TEMPERATURE T2 EN SORTIE DU TRONCON
  357. C----------------------------------------------------
  358.  
  359. HEFF = 2 * HM
  360. C bt3 traite l echange de chaleur par unite de largeur
  361. Qy = Qy/B
  362. T2 = BT3(T1,TP,S,HEFF,Qy,DX,XL)
  363. IF(KIMP.GE.1) THEN
  364. WRITE(6,*) 'bsur2 : T2= ',T2
  365. ENDIF
  366.  
  367. C--------------------------------------------------------
  368. C CALCUL DE LA VARIATION DE TEMPERATURE SUR TRONCON
  369. C--------------------------------------------------------
  370.  
  371.  
  372. DT = T2-T1
  373.  
  374. C--------------------------------------------------------
  375. C CALCUL DE LA VITESSE DE ECOULEMENT U2 SORTIE TRONCON
  376. C---------------------------------------------------------
  377.  
  378. C WRITE(6,*) 'P2 RA T2 ',P2,RA,T2
  379.  
  380. C initialisation avec la valeur a l'entree du troncon
  381.  
  382. ZTRAN = ZV
  383. NZ = 0
  384. 20 CONTINUE
  385. NZ = NZ + 1
  386.  
  387. AA = (QEE2/QAE)* (RV*ZTRAN/RA)
  388. PHI2 = AA/(1.+AA)
  389. PV2 = PHI2 * P2
  390. ROV2 = ROVAP0(PV2,T2)
  391. Z2 = ZVAP0(ROV2,T2)
  392. C write(6,*) ' AA PHI2 ROV2 ',AA,PHI2,ROV2
  393. DELTAZ = ABS(Z2 - ZTRAN)
  394. C WRITE(6,*) ' NZ ZVAP ', NZ,Z2
  395. ZTRAN=Z2
  396. IF ((DELTAZ.GT.1D-5).AND.(NZ.LE.20)) GOTO 20
  397.  
  398. IF (KIMP.GE.2) WRITE(6,*) ' NZ ZVAP ', NZ,Z2
  399.  
  400. ROA2 = (1.-PHI2)*P2/RA/T2
  401. RO2 = ROA2 + ROV2
  402.  
  403. U2 = QA/ROA2
  404.  
  405. C WRITE(6,*) ' ROA ROV RO Z PV ',ROA2,ROV2,RO2,Z2,PV2
  406. C WRITE(6,*) ' PHI U ',PHI2,U2
  407.  
  408.  
  409. IF(KIMP.GE.2) THEN
  410. write(6,1200) X,(P2/P0),(T2-T0),PV2,Z2
  411. 1200 FORMAT(1X,' bsur2 X P2 T2 PV2 Z2 ',5E12.5)
  412. write(6,1300) QE1,QE2
  413. 1300 FORMAT(1X,' QE1 QE2 ',2E12.5)
  414. ENDIF
  415.  
  416.  
  417. C---------------------------------------------------------
  418. C CALCUL DES VALEURS MOYENNES DE P,T,U,RO
  419. C---------------------------------------------------------
  420.  
  421. T=(T1+T2)/2
  422. P=(P1+P2)/2
  423. RO= (RO1+RO2)/2
  424. C write(6,*) ' RO1 RO2 ',RO1,RO2
  425. U=Q/RO
  426.  
  427. C write(6,*) ' T P RO U ',T,P,RO,U
  428.  
  429. QE = (QE1 + QE2)/2
  430. Q = QA + QE
  431.  
  432. QW2=(QEE2-QEE1)/DX/XL/B
  433.  
  434. ERT=ABS((T2-TT2)/T2)
  435. ERP=ABS((P2-PP2)/P2)
  436.  
  437. IF(((ERT.GT.1E-4).OR.(ERP.GT.1E-4))
  438. $ .AND.(NITER.LE.10)) GOTO 10
  439.  
  440.  
  441. IF(KIMP.GE.2) THEN
  442. write(6,*) ' HM HW HW/HM ',HM,HW,HW/HM
  443. ENDIF
  444.  
  445. C calcul des puissances echangees
  446.  
  447. C RAPH=(CV1*T+XLAT0)/XLAT
  448. C H=HM+HW*RAPH
  449. H = HM
  450.  
  451.  
  452. C-------------------------------------------
  453. C BILAN ENERGIE TOTALE
  454. C-------------------------------------------
  455.  
  456. HH1=(QAE*CA*T1) + (QEE1*HVS0(PV1,T1))
  457. HH2=(QAE*CA*T2) + (QEE2*HVS0(PV2,T2))
  458. C puissance fluide
  459. C ! PF = HH2-HH1
  460. PF = (HH2-HH1) - (DQEE*HVS0(PV1,T1))
  461.  
  462. C puissance paroi
  463. PP = H*(T-TP)*2*DX*XL*B
  464.  
  465.  
  466. IF(KIMP.GE.2) THEN
  467. write(6,1000) X,HH1,HH2,PF,PP
  468. 1000 FORMAT(1X,' bsur2 X HH1 HH2 DH HDT ',5E12.5)
  469. ENDIF
  470. IF(KIMP.NE.0.AND.X.EQ.1.) THEN
  471. write(6,2110) X,RE,BK
  472. 2110 format(1X,'bsur2 X RE BKRO ',3E12.5)
  473. ENDIF
  474. DPF=PF/DX/XL/B
  475. DPP=PP/DX/XL/B
  476.  
  477. RETURN
  478. END
  479.  
  480.  
  481.  
  482.  
  483.  
  484.  
  485.  
  486.  
  487.  

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