CenterForTheBuiltEnvironment / pythermalcomfort
 1 9 ```import math ``` 2 9 ```from pythermalcomfort.utilities import * ``` 3 4 9 ```c_to_k = 273.15 ``` 5 9 ```cp_vapour = 1805.0 ``` 6 9 ```cp_water = 4186 ``` 7 9 ```cp_air = 1004 ``` 8 9 ```h_fg = 2501000 ``` 9 9 ```r_air = 287.055 ``` 10 11 12 9 ```def f_svv(w, h, d): ``` 13 ``` """ Calculates the sky-vault view fraction ``` 14 15 ``` Parameters ``` 16 ``` ---------- ``` 17 ``` w : float ``` 18 ``` width of the window, [m] ``` 19 ``` h : float ``` 20 ``` height of the window, [m] ``` 21 ``` d : float ``` 22 ``` distance between the occupant and the window, [m] ``` 23 24 ``` Returns ``` 25 ``` ------- ``` 26 ``` f_svv : float ``` 27 ``` sky-vault view fraction ranges between 0 and 1 ``` 28 ``` """ ``` 29 30 9 ``` return ( ``` 31 ``` math.degrees(math.atan(h / (2 * d))) ``` 32 ``` * math.degrees(math.atan(w / (2 * d))) ``` 33 ``` / 16200 ``` 34 ``` ) ``` 35 36 37 9 ```def p_sat_torr(tdb): ``` 38 ``` """ Estimates the saturation vapor pressure in [torr] ``` 39 40 ``` Parameters ``` 41 ``` ---------- ``` 42 ``` tdb : float ``` 43 ``` dry bulb air temperature, [C] ``` 44 45 ``` Returns ``` 46 ``` ------- ``` 47 ``` p_sat : float ``` 48 ``` saturation vapor pressure [torr] ``` 49 ``` """ ``` 50 9 ``` return math.exp(18.6686 - 4030.183 / (tdb + 235.0)) ``` 51 52 53 9 ```def v_relative(v, met): ``` 54 ``` """ Estimates the relative air velocity which combines the average air velocity of ``` 55 ``` the space plus the relative air velocity caused by the body movement. ``` 56 57 ``` Parameters ``` 58 ``` ---------- ``` 59 ``` v : float ``` 60 ``` air velocity measured by the sensor, [m/s] ``` 61 ``` met : float ``` 62 ``` metabolic rate, [met] ``` 63 64 ``` Returns ``` 65 ``` ------- ``` 66 ``` vr : float ``` 67 ``` relative air velocity, [m/s] ``` 68 ``` """ ``` 69 70 9 ``` if met > 1: ``` 71 0 ``` return round(v + 0.3 * (met - 1), 3) ``` 72 ``` else: ``` 73 0 ``` return v ``` 74 75 76 9 ```def clo_dynamic(clo, met, standard="ASHRAE"): ``` 77 ``` """ Estimates the dynamic clothing insulation of a moving occupant. The activity as ``` 78 ``` well as the air speed modify the insulation characteristics of the clothing and the ``` 79 ``` adjacent air layer. Consequently the ISO 7730 states that the clothing insulation ``` 80 ``` shall be corrected [2]_. The ASHRAE 55 Standard, instead, only corrects for the effect ``` 81 ``` of the body movement, and states that the correction is permitted but not required. ``` 82 83 ``` Parameters ``` 84 ``` ---------- ``` 85 ``` clo : float ``` 86 ``` clothing insulation, [clo] ``` 87 ``` met : float ``` 88 ``` metabolic rate, [met] ``` 89 ``` standard: str (default="ASHRAE") ``` 90 ``` - If "ASHRAE", uses Equation provided in Section 5.2.2.2 of ASHRAE 55 2017 ``` 91 92 ``` Returns ``` 93 ``` ------- ``` 94 ``` clo : float ``` 95 ``` dynamic clothing insulation, [clo] ``` 96 ``` """ ``` 97 98 9 ``` if standard.lower() not in ["ashrae"]: ``` 99 0 ``` raise ValueError( ``` 100 ``` "PMV calculations can only be performed in compliance with ISO or ASHRAE " ``` 101 ``` "Standards" ``` 102 ``` ) ``` 103 104 9 ``` if 1.2 < met < 2: ``` 105 0 ``` return round(clo * (0.6 + 0.4 / met), 3) ``` 106 ``` else: ``` 107 9 ``` return clo ``` 108 109 110 9 ```def running_mean_outdoor_temperature(temp_array, alpha=0.8, units="SI"): ``` 111 ``` """ Estimates the running mean temperature ``` 112 113 ``` Parameters ``` 114 ``` ---------- ``` 115 ``` temp_array: list ``` 116 ``` array containing the mean daily temperature in descending order (i.e. from ``` 117 ``` newest/yesterday to oldest) :math:`[\Theta_{day-1}, \Theta_{day-2}, \dots , ``` 118 ``` \Theta_{day-n}]`. ``` 119 ``` Where :math:`\Theta_{day-1}` is yesterday's daily mean temperature. The EN ``` 120 ``` 16798-1 2019 [3]_ states that n should be equal to 7 ``` 121 ``` alpha : float ``` 122 ``` constant between 0 and 1. The EN 16798-1 2019 [3]_ recommends a value of 0.8, ``` 123 ``` while the ASHRAE 55 2017 recommends to choose values between 0.9 and 0.6, ``` 124 ``` corresponding to a slow- and fast- response running mean, respectively. ``` 125 ``` Adaptive comfort theory suggests that a slow-response running mean (alpha = ``` 126 ``` 0.9) could be more appropriate for climates in which synoptic-scale (day-to- ``` 127 ``` day) temperature dynamics are relatively minor, such as the humid tropics. ``` 128 ``` units: str default="SI" ``` 129 ``` select the SI (International System of Units) or the IP (Imperial Units) system. ``` 130 131 ``` Returns ``` 132 ``` ------- ``` 133 ``` t_rm : float ``` 134 ``` running mean outdoor temperature ``` 135 ``` """ ``` 136 137 9 ``` if units.lower() == "ip": ``` 138 9 ``` for ix, x in enumerate(temp_array): ``` 139 9 ``` temp_array[ix] = units_converter(tdb=temp_array[ix])[0] ``` 140 141 9 ``` coeff = [alpha ** ix for ix, x in enumerate(temp_array)] ``` 142 9 ``` t_rm = sum([a * b for a, b in zip(coeff, temp_array)]) / sum(coeff) ``` 143 144 9 ``` if units.lower() == "ip": ``` 145 9 ``` t_rm = units_converter(tmp=t_rm, from_units="si")[0] ``` 146 147 9 ``` return round(t_rm, 1) ``` 148 149 150 9 ```def units_converter(from_units="ip", **kwargs): ``` 151 ``` """ Converts IP values to SI units ``` 152 153 ``` Parameters ``` 154 ``` ---------- ``` 155 ``` from_units: str ``` 156 ``` specify system to convert from ``` 157 ``` **kwargs : [t, v] ``` 158 159 ``` Returns ``` 160 ``` ------- ``` 161 ``` converted values in SI units ``` 162 ``` """ ``` 163 9 ``` results = list() ``` 164 9 ``` if from_units == "ip": ``` 165 9 ``` for key, value in kwargs.items(): ``` 166 9 ``` if "tmp" in key or key == "tr" or key == "tdb": ``` 167 9 ``` results.append((value - 32) * 5 / 9) ``` 168 9 ``` if key in ["v", "vr", "vel"]: ``` 169 9 ``` results.append(value / 3.281) ``` 170 9 ``` if key == "area": ``` 171 9 ``` results.append(value / 10.764) ``` 172 9 ``` if key == "pressure": ``` 173 9 ``` results.append(value * 101325) ``` 174 175 9 ``` elif from_units == "si": ``` 176 9 ``` for key, value in kwargs.items(): ``` 177 9 ``` if "tmp" in key or key == "tr" or key == "tdb": ``` 178 9 ``` results.append((value * 9 / 5) + 32) ``` 179 9 ``` if key in ["v", "vr", "vel"]: ``` 180 0 ``` results.append(value * 3.281) ``` 181 9 ``` if key == "area": ``` 182 0 ``` results.append(value * 10.764) ``` 183 9 ``` if key == "pressure": ``` 184 0 ``` results.append(value / 101325) ``` 185 186 9 ``` return results ``` 187 188 189 9 ```def t_o(tdb, tr, v): ``` 190 ``` """ Calculates operative temperature in accordance with ISO 7726:1998 [5]_ ``` 191 192 ``` Parameters ``` 193 ``` ---------- ``` 194 ``` tdb: float ``` 195 ``` air temperature, [°C] ``` 196 ``` tr: float ``` 197 ``` mean radiant temperature temperature, [°C] ``` 198 ``` v: float ``` 199 ``` air velocity, [m/s] ``` 200 201 ``` Returns ``` 202 ``` ------- ``` 203 ``` to: float ``` 204 ``` operative temperature, [°C] ``` 205 ``` """ ``` 206 207 9 ``` return (tdb * math.sqrt(10 * v) + tr) / (1 + math.sqrt(10 * v)) ``` 208 209 210 9 ```def enthalpy(tdb, hr): ``` 211 ``` """ Calculates air enthalpy ``` 212 213 ``` Parameters ``` 214 ``` ---------- ``` 215 ``` tdb: float ``` 216 ``` air temperature, [°C] ``` 217 ``` hr: float ``` 218 ``` humidity ratio, [kg water/kg dry air] ``` 219 220 ``` Returns ``` 221 ``` ------- ``` 222 ``` enthalpy: float ``` 223 ``` enthalpy [J/kg dry air] ``` 224 ``` """ ``` 225 226 9 ``` h_dry_air = cp_air * tdb ``` 227 9 ``` h_sat_vap = h_fg + cp_vapour * tdb ``` 228 9 ``` h = h_dry_air + hr * h_sat_vap ``` 229 230 9 ``` return round(h, 2) ``` 231 232 233 9 ```def p_sat(tdb): ``` 234 ``` """ Calculates vapour pressure of water at different temperatures ``` 235 236 ``` Parameters ``` 237 ``` ---------- ``` 238 ``` tdb: float ``` 239 ``` air temperature, [°C] ``` 240 241 ``` Returns ``` 242 ``` ------- ``` 243 ``` p_sat: float ``` 244 ``` operative temperature, [Pa] ``` 245 ``` """ ``` 246 247 9 ``` ta_k = tdb + c_to_k ``` 248 9 ``` c1 = -5674.5359 ``` 249 9 ``` c2 = 6.3925247 ``` 250 9 ``` c3 = -0.9677843 * math.pow(10, -2) ``` 251 9 ``` c4 = 0.62215701 * math.pow(10, -6) ``` 252 9 ``` c5 = 0.20747825 * math.pow(10, -8) ``` 253 9 ``` c6 = -0.9484024 * math.pow(10, -12) ``` 254 9 ``` c7 = 4.1635019 ``` 255 9 ``` c8 = -5800.2206 ``` 256 9 ``` c9 = 1.3914993 ``` 257 9 ``` c10 = -0.048640239 ``` 258 9 ``` c11 = 0.41764768 * math.pow(10, -4) ``` 259 9 ``` c12 = -0.14452093 * math.pow(10, -7) ``` 260 9 ``` c13 = 6.5459673 ``` 261 262 9 ``` if ta_k < c_to_k: ``` 263 0 ``` pascals = math.exp( ``` 264 ``` c1 / ta_k ``` 265 ``` + c2 ``` 266 ``` + ta_k * (c3 + ta_k * (c4 + ta_k * (c5 + c6 * ta_k))) ``` 267 ``` + c7 * math.log(ta_k) ``` 268 ``` ) ``` 269 ``` else: ``` 270 9 ``` pascals = math.exp( ``` 271 ``` c8 / ta_k ``` 272 ``` + c9 ``` 273 ``` + ta_k * (c10 + ta_k * (c11 + ta_k * c12)) ``` 274 ``` + c13 * math.log(ta_k) ``` 275 ``` ) ``` 276 277 9 ``` return round(pascals, 1) ``` 278 279 280 9 ```def psy_ta_rh(tdb, rh, patm=101325): ``` 281 ``` """ Calculates psychrometric values of air based on dry bulb air temperature and ``` 282 ``` relative humidity. ``` 283 ``` For more accurate results we recommend the use of the the Python package ``` 284 ``` `psychrolib`_. ``` 285 286 ``` .. _psychrolib: https://pypi.org/project/PsychroLib/ ``` 287 288 ``` Parameters ``` 289 ``` ---------- ``` 290 ``` tdb: float ``` 291 ``` air temperature, [°C] ``` 292 ``` rh: float ``` 293 ``` relative humidity, [%] ``` 294 ``` patm: float ``` 295 ``` atmospheric pressure, [Pa] ``` 296 297 ``` Returns ``` 298 ``` ------- ``` 299 ``` p_vap: float ``` 300 ``` partial pressure of water vapor in moist air, [Pa] ``` 301 ``` hr: float ``` 302 ``` humidity ratio, [kg water/kg dry air] ``` 303 ``` t_wb: float ``` 304 ``` wet bulb temperature, [°C] ``` 305 ``` t_dp: float ``` 306 ``` dew point temperature, [°C] ``` 307 ``` h: float ``` 308 ``` enthalpy [J/kg dry air] ``` 309 ``` """ ``` 310 9 ``` psat = p_sat(tdb) ``` 311 9 ``` pvap = rh / 100 * psat ``` 312 9 ``` hr = 0.62198 * pvap / (patm - pvap) ``` 313 9 ``` tdp = t_dp(tdb, rh) ``` 314 9 ``` twb = t_wb(tdb, rh) ``` 315 9 ``` h = enthalpy(tdb, hr) ``` 316 317 9 ``` return {"p_sat": psat, "p_vap": pvap, "hr": hr, "t_wb": twb, "t_dp": tdp, "h": h} ``` 318 319 320 9 ```def t_wb(tdb, rh): ``` 321 ``` """ Calculates the wet-bulb temperature using the Stull equation [6]_ ``` 322 323 ``` Parameters ``` 324 ``` ---------- ``` 325 ``` tdb: float ``` 326 ``` air temperature, [°C] ``` 327 ``` rh: float ``` 328 ``` relative humidity, [%] ``` 329 330 ``` Returns ``` 331 ``` ------- ``` 332 ``` tdb: float ``` 333 ``` wet-bulb temperature, [°C] ``` 334 ``` """ ``` 335 9 ``` twb = round( ``` 336 ``` tdb * math.atan(0.151977 * (rh + 8.313659) ** (1 / 2)) ``` 337 ``` + math.atan(tdb + rh) ``` 338 ``` - math.atan(rh - 1.676331) ``` 339 ``` + 0.00391838 * rh ** (3 / 2) * math.atan(0.023101 * rh) ``` 340 ``` - 4.686035, ``` 341 ``` 1, ``` 342 ``` ) ``` 343 9 ``` return twb ``` 344 345 346 9 ```def t_dp(tdb, rh): ``` 347 ``` """ Calculates the dew point temperature. ``` 348 349 ``` Parameters ``` 350 ``` ---------- ``` 351 ``` tdb: float ``` 352 ``` dry-bulb air temperature, [°C] ``` 353 ``` rh: float ``` 354 ``` relative humidity, [%] ``` 355 356 ``` Returns ``` 357 ``` ------- ``` 358 ``` t_dp: float ``` 359 ``` dew point temperature, [°C] ``` 360 ``` """ ``` 361 362 9 ``` c = 257.14 ``` 363 9 ``` b = 18.678 ``` 364 9 ``` d = 234.5 ``` 365 366 9 ``` gamma_m = math.log(rh / 100 * math.exp((b - tdb / d) * (tdb / (c + tdb)))) ``` 367 368 9 ``` return round(c * gamma_m / (b - gamma_m), 1) ``` 369 370 371 9 ```def t_mrt(tg, tdb, v, d=0.15, emissivity=0.9): ``` 372 ``` """ Converts globe temperature reading into mean radiant temperature in accordance ``` 373 ``` with ISO 7726:1998 [5]_ ``` 374 375 ``` Parameters ``` 376 ``` ---------- ``` 377 ``` tg: float ``` 378 ``` globe temperature, [°C] ``` 379 ``` tdb: float ``` 380 ``` air temperature, [°C] ``` 381 ``` v: float ``` 382 ``` air velocity, [m/s] ``` 383 ``` d: float ``` 384 ``` diameter of the globe, [m] default 0.15 m ``` 385 ``` emissivity: float ``` 386 ``` emissivity of the globe temperature sensor, default 0.9 ``` 387 388 ``` Returns ``` 389 ``` ------- ``` 390 ``` tr: float ``` 391 ``` mean radiant temperature, [°C] ``` 392 ``` """ ``` 393 9 ``` tg += c_to_k ``` 394 9 ``` tdb += c_to_k ``` 395 396 ``` # calculate heat transfer coefficient ``` 397 9 ``` h_n = 1.4 * (abs(tg - tdb) / d) ** 0.25 # natural convection ``` 398 9 ``` h_f = 6.3 * v ** 0.6 / d ** 0.4 # forced convection ``` 399 400 ``` # get the biggest between the tow coefficients ``` 401 9 ``` h = max(h_f, h_n) ``` 402 9 ``` print(h_n, h_f, h) ``` 403 404 9 ``` tr = (tg ** 4 + h * (tg - tdb) / (emissivity * (5.67 * 10 ** -8))) ** 0.25 - c_to_k ``` 405 406 9 ``` return round(tr, 1) ```

Read our documentation on viewing source code .