Scilab: Stoffwerte von feuchter Luft

// Funktionsbiliothek zur Berechnung der Stoffeigenschaften von feuchter Luft mit der Magnus-Formel
//Ausgabe der Funktion, die in der Bibliothek enthalten sind
printf('Funktionsbibliothek Stoffwerte Wasser nach der Magnusformel\n')
printf('p = SaettDruckWasserdampfMg(t)  Druck in Pa\n')
printf('rhoL = DichteWasserdampfGG(t) Dichte in kg/m³\n')
printf('x = WasseranteilLuft(p,phi,t)  in kg[H2O]/kg[tr. Luft] \n')
printf('h = h1x(p,x,t) (Enthalpie in J/kg, Flüsiges Wasser in kg[H20 fl.]/kg[tr. Luft])\n')
printf('dhdt = dhdT(p,x,t)  dh/dT in J/kgK\n')
printf('\n')

// Berechnung des Sättigungsdampfdrucks von Wasser nach der Magnusformel
function p = SaettDruckWasserdampfMg(t)
 p0 = 1917;         // Bezugsdruck bei 290 K
 T0 = 290;          // Bezugstemperatur in K
 rv = 2461000;      // Verdampfungswärme bei 290 K in J/kg
 Rw = 461.5;        // Gaskonstante in J/kgK
 T = t + 273.15;
 p = p0 * exp(((rv)/(Rw*T0))*((T-T0)/(T)));      // Druck wird in Pa ausgegeben
endfunction

// Berechnung der Dichte von Wasserdampf nach der allgemeinen Gasgleichung
function rhoL = DichteWasserdampfGG(t)
 RL = 461.5;        // Gaskonstante
 T = t + 237.15;
 p = SaettDruckWasserdampfMg(t);
 rhoL =(p)/(RL * T);
endfunction

// Berechnung der Absoluten Luftfeuchte anhand der relativen Luftfeuchte
function x = WasseranteilLuft(p,phi,t)
 ML = 0.02895;               // Molmasse der Luft in kg/mol
 MH2O = 0.01801;             // Molmasse des Wassers in kg/mol
 ps = SaettDruckWasserdampfMg(t);
 x = (MH2O/ML)*((phi*ps)/(p-phi*ps));  // Einheit in kg[H2O]/kg[tr. Luft]
endfunction

// Berechung der spezifischen Enthalpie von feuchter Luft
function h = h1x(p,x,t)           // Ausgabe: (Enthalpie in J/kg, Flüsiges Wasser in kg[H20 fl.]/kg[tr. Luft])
 cpl = 1004;                     //Spezifische Wärmekapazität der Luft in J/kgK
 cpd = 1860;                     //Spezifische Wärmekapazität des Wasserdampfs in J/kgK
 cw = 4189;                      //Spezifische Wärmekapazität des flüssigen Wassers in J/kgK
 ML = 0.02895;               // Molmasse der Luft in kg/mol
 MH2O = 0.01801;             // Molmasse des Wassers in kg/mol
 rv = 2500900;               // Verdampfungswärme bei 0°C in J/kg
 ps = SaettDruckWasserdampfMg(t);
 xs = (MH2O/ML)*((ps)/(p-ps));  // Einheit in kg[H2O]/kg[tr. Luft]
 if x <=xs                       // kein Nebel
   h(1,1) = cpl*t + x*(rv + cpd*t);   // Enthalpie der feuchten Luft in J/kg
   h(1,2) = 0;
 end
 //
 if x>xs                         // Nebelgebiet
   h(1,1) = cpl*t + xs*(rv + cpd*t)+ (x-xs)*cw*t;
   h(1,2) = x-xs;
 end
endfunction

// Berechnung der ersten Ableitung der Enthalpie nach der Temperatur
function dhdt = dhdT(p,x,t)
 dt = 1;                 // Temperaturdifferenz zur Berechnung von dh/dT
 h2 = h1x(p,x,t+dt);     //   Die Funktion h1x gibt einen Vektor aus, von dem aber nur die erste Spalte gebraucht wird.
 h1 = h1x(p,x,t);
 dhdt = (h2(1,1) - h1(1,1))/dt;
endfunction