psychrometrics / psychrolib

📚 Library of psychrometric functions to calculate 🌡️ thermodynamic properties of air for Python, C, C#, Fortran, R, JavaScript and VBA/Excel
MIT License
222 stars 61 forks source link

Find the intersection of relative humidity and enthalpy lines #70

Open heymikid opened 3 years ago

heymikid commented 3 years ago

I am enjoying this library. Thank you for making it available.

I had to find the intersection point between RH lines and moist air enthalpy lines. I was using the js version and ended up writing this outside the library, following the bisection search code found in GetTWetBulbFromHumRatio.

You may consider adding this (and another one like it - to find the intersection of TWetBulb with RH) to the library itself if you wish.

function GetTDryBulbFromEnthalpyRelHum(MoistAirEnthalpy, RelHum, Pressure) {
  /*
  We want to determine the psychrometric properties from the moist air enthalpy 
  and relative humidity. These two lines (RH(T), h(T)) cross at one point. We find TDryBulb of this point.
  */
  var T_TOLERANCE = 1E-10;     

  var TDryBulbInf, TDryBulbSup, TDryBulb;
  var HumRatio, currRelHum;
  var index = 1;

  //initial guesses
  TDryBulbInf = psychrolib.GetTDryBulbFromEnthalpyAndHumRatio(MoistAirEnthalpy, 0);
  TDryBulbSup = psychrolib.GetTWetBulbFromRelHum(TDryBulbInf, 0, Pressure) - 2; //wet bulb T is above the enthalpy line, not good enough for high target RH. -2deg is enough to cover it
  TDryBulb = (TDryBulbSup + TDryBulbInf) / 2.;

  //bisection loop
  while ((TDryBulbInf - TDryBulbSup) > T_TOLERANCE) {
    //compute the relative humidity for the current guess of TDryBulb (first get Humidity Ratio)
    HumRatio = psychrolib.GetHumRatioFromEnthalpyAndTDryBulb(MoistAirEnthalpy, TDryBulb);
    currRelHum = psychrolib.GetRelHumFromHumRatio(TDryBulb,HumRatio, Pressure);

    //get new bounds
    if (currRelHum < RelHum)
      TDryBulbInf = TDryBulb;
    else
      TDryBulbSup = TDryBulb;

    //new guess for TDryBulb
    TDryBulb = (TDryBulbSup + TDryBulbInf) / 2.;

    if (index > MAX_ITER_COUNT)
      throw new Error("Convergence not reached. Stopping.");

    index++;
  }
  return TDryBulb;
}
dmey commented 3 years ago

Thanks @heymikid -- would you be willing to submit a PR for this for most/all supported languages? I am a bit short of time at the moment...