TypeError being raised from some state calculations

[TimTimMadden]

Hi,

Firstly, thanks for this well written and concise library - I'm attempting to simulate adiabatic cooling of air and have found this to be very helpful!

I noticed a small bug when iterating through a bunch of scenarios (some purely theoretical states).

Firstly, there is sensible behaviour around RH == 0 behaviour, which generates a math error. As it's practically impossible to have this state, I'm happy with this!

However, if something causes a state calculation to have a WBT < 0 C, or causes the lower wet bulb temp window edge to dip below zero, a type error is thrown because of a None value in the depths of the calculation:

>>> si.state('DBT', 283.15, 'WBT', 273, 101325)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/srv/anaconda2/envs/devicegenerator/lib/python3.6/site-packages/psypy/psySI.py", line 
213, in state
    W=__W_DBT_WBT_P(DBT, WBT, P)
  File "/srv/anaconda2/envs/devicegenerator/lib/python3.6/site-packages/psypy/psySI.py", line 
342, in __W_DBT_WBT_P
    return ((2501-2.326*WBT)*__W_DBT_RH_P(WBT+273.15,1,P)-1.006*(DBT-WBT))/\
TypeError: unsupported operand type(s) for *: 'float' and 'NoneType'

and

>>> si.state('DBT', 303.15, 'RH', 0.01, 101325)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/srv/anaconda2/envs/devicegenerator/lib/python3.6/site-packages/psypy/psySI.py", line 
222, in state
    WBT=__WBT_DBT_W_P(DBT, W, P)
  File "/srv/anaconda2/envs/devicegenerator/lib/python3.6/site-packages/psypy/psySI.py", line 
352, in __WBT_DBT_W_P
    Ws=__W_DBT_WBT_P(DBT, WBT, P)
  File "/srv/anaconda2/envs/devicegenerator/lib/python3.6/site-packages/psypy/psySI.py", line 
342, in __W_DBT_WBT_P
    return ((2501-2.326*WBT)*__W_DBT_RH_P(WBT+273.15,1,P)-1.006*(DBT-WBT))/\
TypeError: unsupported operand type(s) for *: 'float' and 'NoneType'

This is happening in the following class (psySI.py, lines 337 to 343):

# ASHRAE 2009 Chapter 1 Equation 35
def __W_DBT_WBT_P(DBT, WBT, P):
    if __valid_DBT(DBT):
        DBT=DBT-273.15
        WBT=WBT-273.15
        return ((2501-2.326*WBT)*__W_DBT_RH_P(DBT+273.15,1,P)-1.006*(DBT-WBT))/\
               (2501+1.86*DBT-4.186*WBT)

which calls (psySI.py, lines 326 to 330):

# ASHRAE 2009 Chapter 1 Equation 22 and Equation 24
def __W_DBT_RH_P(DBT, RH, P):
    if __valid_DBT(DBT):
        Pw=RH*__Pws(DBT)
        return 0.621945*Pw/(P-Pw)

Note that in the return statement of the first class, the Wet Bulb Temperature is being passed to a fn named for the Dry Bulb Temperature, which is then checked against the valid DBT check method. The need for a DBT check is obvious as behaviour of air changes for DBT < 0 C, but unfortunately this causes nothing to be returned from these two classes in the instance that we are passing a -ve C wet bulb temperature.

The same happens here when called from a second class (psySI.py, lines 179 - 187):

# ASHRAE 2009 Chapter 1 Equation 6
def __Pws(DBT):
    if __valid_DBT(DBT):
        C8=-5.8002206*10**3
        C9=1.3914993
        C10=-4.8640239*10**-2
        C11=4.1764768*10**-5
        C12=-1.4452093*10**-8
        C13=6.5459673
        return m.exp(C8/DBT+C9+C10*DBT+C11*DBT**2+C12*DBT**3+C13*m.log(DBT))

Unfortunately I don't have access to the 2009 edition of the ASHRAE handbook - the copy I have here is the older 1993 version with rather different equation numbers!

I originally attempted to pass dry bulb temp in line 342 instead, but that caused my unit tests (sanity checked against values read of the ASHRAE psychrometric chart) to fail and be obviously incorrect. Instead I added optional params to these two broken classes, like so:

lines 179 - 187:

# ASHRAE 2009 Chapter 1 Equation 6
def __Pws(DBT, check_dry_bulb=True):
    if __valid_DBT(DBT) or not check_dry_bulb:
        C8=-5.8002206*10**3
        C9=1.3914993
        C10=-4.8640239*10**-2
        C11=4.1764768*10**-5
        C12=-1.4452093*10**-8
        C13=6.5459673
        return m.exp(C8/DBT+C9+C10*DBT+C11*DBT**2+C12*DBT**3+C13*m.log(DBT))

lines 326 - 330:

# ASHRAE 2009 Chapter 1 Equation 22 and Equation 24
def __W_DBT_RH_P(DBT, RH, P, check_for_dry_bulb=True):
    if __valid_DBT(DBT) or not check_dry_bulb:
        Pw=RH*__Pws(DBT)
        return 0.621945*Pw/(P-Pw)

line 342 - 343:
```python
        return ((2501-2.326*WBT)*__W_DBT_RH_P(WBT+273.15,1,P, False)-1.006*(DBT-WBT))/\
              (2501+1.86*DBT-4.186*WBT)

This allows for sensible handling of low RH and negative wet bulb values, and all unit tests still pass:

>>> si.state('DBT', 303.15, 'RH', 0.01, 101325)
[30.0, 30.846649010311683, 0.01, 0.8591488970583029, 0.00026073568926458275, 10.812545040052214]

and

>>> si.state('DBT', 303.15, 'WBT', 271.35, 101325)
[5.0, 6.467328987664403, 0.1068161111995133, 0.7886921777450676, 0.000572572595970363, -1.8000000000000114]

I'm writing this story to help other users, if you're happy with this fix I can add a pull request with some logging changes and these improvements implemented.

Thanks,

Tim

[TimTimMadden]

updating comment above with pythonic highlighting in markdown

[longapalooza]

Thanks for pointing this out. I'll get around to creating a fix for it eventually but for the time being, your solution will work great