* fichier : condens.dgibi ************************************************************************ ************************************************************************ GRAPH = FAUX ; ************************************************************************ * condens.dgibi : bas Mach + condensation en paroi * ************************************************************************ * MODELE MISTRA AVEC CONDENSEURS * * ------------------------------ * * INJECTION DE VAPEUR DANS UNE ENCEINTE FERMEE CONTENANT DE L'AIR * * TEMPERATURE DE PAROI IMPOSEE (PAROI VERTICALE) NON CALCULEE * * MODELE DE CONDENSATION : CHILTON-COLBURN + ANALOGIE HMT * * ECHANGES PAR CONVECTION ET CONDENSATION * * EQUATION D'ENERGIE EN TEMPERATURE * * ECOULEMENTS TURBULENTS * * MAILLAGE 2D AXI FORMULATION EF QUA8 * * P. CORNET SEMT/TTMF 23 DECEMBRE 1998 * ************************************************************************ * ************************************************************************ * PROCEDURE FILTREKE * * VERSION non encore évoluée à ce jour * ************************************************************************ * * 4 arguments : U0 L0 NU UN * * Filtre sur K et Epsilon * - Echelle de vitesse (K**0.5) inférieure à une fraction (alfk) * de Uref (vitesse caractéristique) (alfk=1 pour l'instant) * Uref=max(UN,U0) * - K > K0 * - Epsilon tel que l'echelle de longueur reste inférieure * à (L0/a) où L0 = diamètre enceinte et a=f(Re) * => Nut < Uref*L0/a rvp=rv.'PRESSION' ; iarg=rx.'IARG' ; si( non ( ega iarg 4)) ; finsi ; U1=rv.'INCO'.(rx.'ARG1') ; sinon ; U1=rx.'ARG1' ; sinon ; finsi ; finsi ; L0=rv.'INCO'.(rx.'ARG2') ; sinon ; L0=rx.'ARG2' ; sinon ; finsi ; finsi ; NU=rv.'INCO'.(rx.'ARG3') ; sinon ; NU=rx.'ARG3' ; sinon ; finsi ; finsi ; UN=rv.'INCO'.(rx.'ARG4') ; sinon ; UN=rx.'ARG4' ; sinon ; finsi ; finsi ; si( non ( ega nic 2)) ; finsi ; en=rv.'INCO'.nom2 ; kn=rv.'INCO'.nom1 ; Rec=100.; k0 = 1.e-10 ; cnu=0.09; mdu=mdu ** 0.5 ; rv.'INCO'.nom2=en ; rv.'INCO'.nom1=kn ; finproc ; ************************************************************************ * PROCEDURE CALCUL1 * ************************************************************************ DEBPROC CALCUL1 ; $MT = RV.'DOMAINE' ; * on récupère le pas de temps (tn - tn-1) DT = RV.PASDETPS.'DELTAT-1' ; SI (RV.PASDETPS.'NUPASDT' EGA 1) ; RV.PASDETPS.'DELTAT' = 1.0E-5 ; DT = 0. ; FINSI ; * filtre sur la fraction massique de vapeur et calcul fraction air YH2O = RV.INCO.'YH2O' ; RV.INCO.'YAIR' = YAIR ; * calcul de Cv, Cp et R du mélange RV.INCO.'CV' = CV ; RV.INCO.'CP' = CP ; RV.INCO.'RN' = RR ; * calcul de la densité totale moyenne 0D (kg/m3) Q0 = RV.INCO.'DEBIT' ; QC = RV.INCO.'QC' ; DRO = DT*(Q0+QC)/VOLT ; ROT0D = RV.INCO.'ROT0D' + DRO ; RV.INCO.'ROT0D' = ROT0D ; * calcul de la densité de vapeur moyenne 0D (kg/m3) QV = Q0*YH2Ojet ; DROV = DT*(QV+QC)/VOLT ; ROV0D = RV.INCO.'ROV0D' + DROV ; RV.INCO.'ROV0D' = ROV0D ; * calcul de l'énergie emmenée par la condensation (W) ENC = RV.'DI' '*' ENC ; * calcul de l'énergie emmenée par la convection (W) * calcul de l'énergie moyenne 0D (J/m3) DES = (Q0*EMjet) + Etran - ENC - EHI ; DES = DES*DT/VOLT ; ES0D = RV.INCO.'ES0D' + DES ; RV.INCO.'ES0D' = ES0D ; * Cp et Cv 0D du mélange YV0D = ROV0D/ROT0D ; YA0D = 1.0 - YV0D ; CP0D = (YV0D*CPH2O) + (YA0D*CPAIR) ; CV0D = (YV0D*CVH2O) + (YA0D*CVAIR) ; * calcul de la pression moyenne 0D GAM0D = CP0D/CV0D ; RV.INCO.'PR0D' = (GAM0D - 1.0)*ES0D ; * calcul de la masse volumique du melange DRO1 = RV.INCO.'ROT0D' - ROM ; RV.INCO.'RHO1' = RHO ; * contrainte sur rhoh2o pour conservation 0D DROV = RV.INCO.'ROV0D' - ROVM ; DYV = DROV/ROM ; RV.INCO.'YH2O' = YH2OF ; * contrainte sur l'énergie e=rho.Cv.T pour conservation 0D DEN = RV.INCO.'ES0D' - ENM ; DTN = DEN/RCVM ; * calcul du terme (rho-rhom)*g de la QDM * calcul la vitesse à partir de la qdm GN = RV.INCO.'GN' ; RV.INCO.'UN' = un ; * calcul de nut et uet (pour FPT) *-------------------------------------------------* * calcul du flux condensé Jv par Chilton-Colburn * * et la corrélation de convection naturelle * * Sh = kL/Dv = 0.13 (Gr Sc)**1/3 * * Jv = k ro (Yv - Yvsat) en kg/m²s * *-------------------------------------------------* * pression de saturation à l'interface * pression de vapeur pres de la paroi * masse volumique à l'interface * fraction massique vapeur à l'interface RV.INCO.'YVS' = YVI ; * coefficient de diffusion de la vapeur dans l'air (m2/s) AA = 2.2 / RV.INCO.'PM' ; * coefficient d'échange de masse k (m/s) BB = 0.13 * (9.81 * RV.INCO.'ROT0D' / MU)**0.33333 ; DRO = ABS DRO ; RV.INCO.'K' = KK ; * fonction indicatrice de la condensation IND = DPV MASQUE 'EGSUP' 0. ; * coefficient d'échange pour yh2o (kg/m2s) RV.INCO.'KS' = KRO ; * calcul de la masse de vapeur condensée QC (kg/s) QC = RV.'DI' '*' FCOND ; RV.INCO.'FCOND' = FCOND ; RV.INCO.'QC' = QC ; * suivi masse condensée RV.INCO.'MCOND' = RV.INCO.'MCOND' + (QC*DT) ; * suivi des inconnues OD et messages DD = RV.PASDETPS.'NUPASDT' ; NN = DD/NSAUV ; LO = (DD-(NSAUV*NN)) EGA 0 ; SI ( LO ) ; MESSAGE '============================================================' ; * suivi pression, masse vol tot, temperature, enthalpie OD, condensation R0D = (YV0D*RH2O) + (YA0D*RAIR) ; PVAP = RV.INCO.'PR0D'*YV0D*RH2O/R0D ; PAIR = RV.INCO.'PR0D'*YA0D*RAIR/R0D ; T0D = ES0D/ROT0D/CV0D - 273.15 ; H0D = GAM0D * ES0D / ROT0D + LAT ; QC1 = -1.* QC ; en2 = ES0D*VOLT ; FINSI ; FINPROC ; ************************************************************************ * PROCEDURE CALCUL2 * ************************************************************************ DEBPROC CALCUL2 ; $MT = RV.'DOMAINE' ; * on récupère le pas de temps Dt = RV.PASDETPS.'DELTAT' ; Dt = Dt * RV.'ALFA' ; * calcul de l'intégrale se trouvant dans le terme Dp/Dt. * calcul de dP/dt * Q0 = 'DBIT' RV.INCO.'GN' $BAS ; dPdt = ZZT + RV.INCO.'DEBIT' + RV.INCO.'QC' ; dPdt = dPdt*(RV.INCO.'PM')/VOLT/RV.INCO.'ROT0D' ; * terme source de l'équation d'énergie (terme en dP/dt) * terme source équation de Poisson Z1 = dPdt/RV.INCO.'PM' ; ZZS = ZZS * (-1.0) ; * calcul de l'évolution de la pression RV.INCO.'PM' = RV.INCO.'PM' + (dPdt*dt) ; * affichage informations toutes les NSAUV iterations ************************************************************************ DD = RV.PASDETPS.'NUPASDT' ; NN = DD/NSAUV ; LO = (DD-(NSAUV*NN)) EGA 0 ; SI ( LO ) ; PMas = MASINI + ( RV.INCO.'DEBIT' * RV.PASDETPS.'TPS') + RV.INCO.'MCOND' ; RHO = RV.INCO.'RHO1' ; RV.INCO.'TEMPS' = RV.INCO.'TEMPS' ET PT ; RV.INCO.'PRESS' = RV.INCO.'PRESS' ET PP ; RV.INCO.'MTOTA' = RV.INCO.'MTOTA' ET PMTO ; RV.INCO.'MTHEO' = RV.INCO.'MTHEO' ET PMTH ; MESSAGE 'MASSES (kg): H2O =' MASH2O 'AIR =' MASAIR 'TOTALE = ' (MASH2O + MASAIR) ; RV.PASDETPS.'TPS'* YH2OJet) + RV.INCO.'MCOND'); RV.INCO.'PMH2O' = RV.INCO.'PMH2O' ET PMH2O ; RV.INCO.'PMAIR' = RV.INCO.'PMAIR' ET PMAIR ; RV.INCO.'PMH2OT' = RV.INCO.'PMH2OT' ET PMH2OT ; RV.INCO.'PMAIRT' = RV.INCO.'PMAIRT' ET PMAIRT ; FINSI ; FINPROC ; ************************************************************************ * GENERATION DU MAILLAGE * ************************************************************************ * * cuve cylindrique de rayon interieur RA * de hauteur HA * de rayon d'injection RJ. *------------------------------------------ 'OPTION' 'MODE' 'AXIS' ; * Dimensions de la cuve et de la zone d'injection HA = 6.2 ; RA = 1.9 ; RJ = 0.1 ; * Nombre d'elements injection et cotés * en MACRO, N réel = n x 2 nj = 1 ; nr = 10 ; nh = 20 ; * Definition des points utiles pour la creation du maillage * --------------------------------------------------------- * * B0| BJ B1 * .----------------| * | | * . | * | | * . | * | | * . | * | | * . | * | | * . | * | | * . | * | | * . | * | | * . | * | | * .----------------| * A0| |AJ A1 * * ---------------------------------------------------------------------- A0 = 0. 0. ; AJ = RJ 0. ; A1 = RA 0. ; B0 = 0. HA ; BJ = RJ HA ; B1 = RA HA ; * Définition des segments A01 = A0J ET AJ1 ; B01 = B0J ET BJ1 ; * Maillage total 'ORIENTER' MT ; * Lignes pour conditions aux limites et post-traitement AXE = AB0 ; BAS = A01 ; INJ = A0J ; BAS2 = AJ1 ; HAUT = B01 ; PCON = AB1 ; * Définition des tables domaine * Maillage final MT = $MT.MAILLAGE ; * TRAC (BAS ET AXE ET PCON ET HAUT) ; * TRAC (MT ET INJ ET PCON ET AXE) ; * TRAC (PTOT ET AXE) ; * Volume, nbre d'éléments, surface condensante ************************************************************************ * DONNEES PHYSIQUES * ************************************************************************ *-------------------------------------- * masses molaires et constantes des gaz *-------------------------------------- MH2O = 18.0E-3 ; MAIR = 28.8E-3 ; Rg = 8.313 ; RH2O = Rg/MH2O ; RAIR = Rg/MAIR ; *------------------------------------------------ * fractions massiques des differents constituants *------------------------------------------------ YAIR = 0.664; YH2O = 1-YAIR ; *--------------------------------------------- * constante du gaz *--------------------------------------------- R = (RH2O*YH2O)+(RAIR*YAIR) ; *------------------------------------------------------------ * pression, temperature et densite du melange (loi d'etat GP) *------------------------------------------------------------ Pbar = 1.5 ; T°C = 100. ; P = Pbar /1.013 ; P = P * 1.0E5 ; T = T°C + 273.15 ; PSI = R*T ; RHO = P/PSI ; *------------------------------------------------------------- * Masse volumique des constituants *------------------------------------------------------------- RHOAIR= P*YAIR/(T*R); RHOH2O= RHOAIR*YH2O/YAIR; *--------------------- * pressions partielles *--------------------- PH2O = RHOH2O*RH2O*T ; PAIR = RHOAIR*RAIR*T ; *---------------------------------------------------------- * fractions molaires des differents constituants du melange *---------------------------------------------------------- XH2O = PH2O/P ; XAIR = PAIR/P ; *--------------------------------------------- * Masse molaire du mélange *--------------------------------------------- M = (XH2O*MH2O) + (XAIR*MAIR) ; *---------------------------------------------------- * Cp et Cv des différents constituants et du mélange *---------------------------------------------------- CPH2O = 1715.632 + (0.552805*T) ; CPN2 = 1006.15 + (0.1387166*T) ; CPO2 = 907.580 + (0.1420522*T) ; CPAIR = (0.2*CPO2) + (0.8*CPN2) ; CVH2O = CPH2O - RH2O ; CVAIR = CPAIR - RAIR ; CP = (YH2O*CPH2O) + (YAIR*CPAIR) ; CV = (YH2O*CVH2O) + (YAIR*CVAIR) ; gamma = CP/(CP-R) ; *------------------------ * Caractéristiques du jet * (Qjet en kg/s) *------------------------ XH2OJet = 1.0 ; XAIRJet = 1.0 - XH2OJet ; YAIRJet = 0. ; YH2OJet = 1-YAIRJet; MJet = (XH2OJet*MH2O) + (XAIRJet*MAIR) ; RJet = Rg/MJet ; CVJet = (YH2OJet*CVH2O) + (YAIRJet*CVAIR) ; SJet = PI*RJ*RJ ; *°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°° * ° * ATTENTION : ° * ° * Pour que le débit effectif soit 'Q0' kg/s ° * ° * il faut que Qjet imposé soit < Q0 ° * ° * ==> imposer un certain Qjet par tatonnement et vérifier ensuite ° * avec l'opérateur DBIT que l'on a bien le débit voulu ° * ° *°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°° * nr = 10 , QJet = 0.069124 => Q0 effectif = 0.149999079999999979 kg/s * TJ°C = 180. ; QJet = 0.069124 ; GJet = QJet/SJet ; * Q0 = QJet ; Q0 = 0.149999079999999979 ; TJet = TJ°C + 273.15 ; PSIJet = Rjet*Tjet ; RhoJet = P/PSIJet ; EMJet = CVJet*TJet ; EJet = Q0*EMJet ; UJet = GJet/RhoJet ; * K et epsilon KJet = 4.3E-3*UJet*UJet ; EPJet = 4.3E-4*(UJet**3.)/RJ ; * énergie de transvasement PUjSj Etran = Q0*Rjet*Tjet ; *------------------------------------------------- * energie, enthalpie massique et masses du melange *------------------------------------------------- ESini = RHO*CV*T ; ETini = ESini*VOLT ; Pini = (gamma-1.)*ESini ; LAT = 2.3E6 ; MASINI = RHO*VOLT ; Hini = gamma*ETini/MASINI + LAT ; MH2OINI = MASINI*YH2O ; MAIRINI = MASINI*YAIR ; *-------------------------------------------------- * vitesses de référence (tracés, filtre, fct paroi) *-------------------------------------------------- ampl = 0.25/Ujet ; U0 = 1.0*UJet ; L0 = 0.7 ; UET0 = Ujet/10.0 ; YP = 0.05 ; YPT= 0.05 ; *------------------------------------------------------------------------ * Prandtl, Schmidt, viscosité (loi de Sutherland), conductivité thermique * coefficients de diffusion moléculaire, coeff. d'échange par convection *------------------------------------------------------------------------ Pr = 0.7 ; Prt = 1.0 ; Sct = 1.0 ; MUvap = 1.6E-5 ; MU = 1.716E-5*((T/273.0)**1.5)*(273.0+110.5)/(T+110.5) ; NU = MU/RHO ; Lambda = Mu*CP/Pr ; alpha = Lambda/RHO/CP ; DH2O = 2.55E-5 ; DH2 = 7.12E-5 ; DO2 = 2.06E-5 ; H1 = 10.0 ; H2 = H1/CP/RHO ; *----------------------------------- * température initiale du condenseur *----------------------------------- TP°C = 40. ; Tcond = TP°C + 273.15 ; *----------------------------------- * Reynolds et Richardson d'injection *----------------------------------- VD = Q0 / (RHOjet*Sjet) ; RI0 = (RHO - RHOjet)*2.0*RJ*9.81 / (RHO*VD*VD) ; RE0 = RHOjet*VD*2.0*RJ / MUvap ; mess '______________________________________________________________'; mess ' ' ; 'MH2Oini =' MH2OINI ; mess '______________________________________________________________'; * opti donn 5 ; *=====================================================================* * TABLE EQEX CONTENANT LES EQUATIONS A RESOUDRE : * * 0) FILTRE K-E * * 1) QUANTITE DE MOUVEMENT --- OPERATEUR NSKE * * 2) EQUATION DE L'ENERGIE --- OPERATEUR TSCAL * * 3) TRANSPORT DES ESPECES --- OPERATEUR TSCAL * * AINSI QUE 2 PROCEDURES (CALCUL1 ET CALCUL2) POUR L'ALGORITHME * * 'FAIBLE MACH' * * 4) CONDITIONS AUX LIMITES * *=====================================================================* 'ZONE' $MT 'OPER' 'CALCUL1' 'ZONE' $MT 'OPER' 'FILTREKE' U0 L0 NU 'UN' 'INCO' 'KN' 'EN' 'ZONE' $PTOT 'OPER' 'FPU' NU 'UET' YP 'INCO' 'UN' 'KN' 'EN' 'OPTI' 'CONS' 'SUPG' 'ZONE' $MT 'OPER' 'NSKE' 'RHOC' MU 'MUT' 'UN' 'ROG' 'INCO' 'GN' 'KN' 'EN' ; RV = 'EQEX' RV 'OPTI' 'NOCONS' 'SUPG' 'ZONE' $MT 'OPER' 'TSCAL' alpha 'UN' 'SE' 'NUT' Prt 'INCO' 'TN' 'ZONE' $PCON 'OPER' 'ECHI' 'KY' 'YVI' 'INCO' 'YH2O' 'OPTI' 'NOCONS' 'SUPG' 'ZONE' $MT 'OPER' 'TSCAL' DH2O 'UN' 0. 'NUT' SCT 'INCO' 'YH2O' 'ZONE' $MT 'OPER' 'CALCUL2' ; *===================== conditions aux limites ========================== 'GN' 'VIMP' (BAS2 et HAUT) 0. 'GN' 'UIMP' INJ 0. 'GN' 'VIMP' INJ GJet 'KN' 'TIMP' INJ KJet 'EN' 'TIMP' INJ EPJet ; RV = 'EQEX' RV 'CLIM' 'TN' 'TIMP' INJ TJet 'YH2O' 'TIMP' INJ YH2OJet ; *=========================================================== * TABLE EQPR POUR L'EQUATION DE POISSON (SOLVEUR ELLIPTIQUE) *=========================================================== 'OPER' 'PRESSION' 'SP' 'PIMP' 0. ; *========================================================== * TABLE INCO CONTENANT LES INCONNUES ET DONNEES DU PROBLEME *========================================================== RV.'PRESSION' = RVP ; RV.'NOMVI' = 'GN' ; RV.'INCO' = 'TABLE' 'INCO' ; RVP.'INCO' = RV.'INCO' ; $FLUID = $MT ; $PAROI = $PTOT ; * inconnues des équations différentielles * matrices masse pour calcul des moyennes aux noeuds * variables 0D : masse totale, masse vapeur, energie, pression RV.INCO.'ROT0D' = RHO ; RV.INCO.'ROV0D' = RHOH2O ; RV.INCO.'ES0D' = ESini ; RV.INCO.'PR0D' = P ; * debit de masse injecté, masse condensée RV.INCO.'DEBIT' = Q0 ; RV.INCO.'MASCD' = 0. ; * inconnues auxiliaires * variables modèle de turbulence * masse volumique * fractions massiques * Cv, R, et Cp du mélange * à la paroi coef d'ech, frac vap, temp, flux de masse RV.INCO.'QC' = 0. ; RV.INCO.'MCOND' = 0. ; * termes sources (Poisson/énergie/QDM) * pression RV.INCO.'PM' = P ; * suivi des valeurs moyennes calculées en 2D * suivi des variables purement 0D *===================================================================* * EXECUTION ET SAUVEGARDE * *===================================================================* * on gardera une valeur tous les NSAUV pas de temps NSAUV = 20 ; RV.'FIDT' = NSAUV ; exec rv ; *===================================================================* * POST TRAITEMENT * *===================================================================* SI GRAPH ; MFLUID = MT ; 'TRAC' gn1 mt cnt 'TRAC' RV.'INCO'.'RHO' MFLUID cnt 'TRAC' RV.'INCO'.'YH2O' MFLUID cnt 'TITR' 'CHAMP DE LA FRACTION MASSIQUE DE VAPEUR YH2O' ; nut = RV.INCO.'NUT' ; 'TRAC' NUTNU MFLUID cnt RV.INCO.'PRESS' ; RV.INCO.'PT0D' ; RV.INCO.'PV0D' ; RV.INCO.'PA0D' ; TAB1 = TABLE ; TAB1.1 = 'REGU MARQ CROI ' ; TAB1.3 = 'REGU MARQ CARR ' ; TAB1 . 'TITRE' = TABLE ; RV.INCO.'T0D' ; RV.INCO.'H0D' ; RV.INCO.'FCDST' ; RV.INCO.'MTHEO' ; RV.INCO.'MTOTA' ; TAB1 = TABLE ; TAB1.1 = 'REGU MARQ CROI' ; TAB1 . 'TITRE' = TABLE ; RV.INCO.'PMH2OT' ; RV.INCO.'PMH2O' ; TAB1 = TABLE ; TAB1.1 = 'REGU MARQ CROI ' ; TAB1 . 'TITRE' = TABLE ; RV.INCO.'PMAIRT' ; RV.INCO.'PMAIR' ; TAB1 = TABLE ; TAB1.1 = 'REGU MARQ CROI ' ; TAB1 . 'TITRE' = TABLE ; RV.INCO.'ETOT' ; RV.INCO.'ET0D' ; TAB1 = TABLE ; TAB1.1 = 'REGU MARQ CROI ' ; TAB1.2 = 'REGU MARQ CARR ' ; TAB1 . 'TITRE' = TABLE ; FINSI ; *===================================================================* * TESTS DE BON FONCTIONNEMENT (à 5 %) * *===================================================================* DELTAP = (P - RV.INCO.'PM' - 1673.)/1673. ; MH2O = RV.INCO.'ROV0D' * VOLT ; DELTAM = (MH2OINI - MH2O - 0.586)/0.586 ; SI ( (ABS DELTAP) > 0.05 ) ; ERREUR 5 ; FINSI ; SI ( (ABS DELTAM) > 0.05 ) ; ERREUR 5 ; FINSI ; FIN ;
© Cast3M 2003 - Tous droits réservés.
Mentions légales