[gradsusr] odd behavior of pow() function

[Jeff Duda]

Hello,
I am trying to compute some thermodynamic parameters, shown below.

 'define e = rhprs*6.11*exp((2.5e6/461.5)*((1/273)-(1/tmpprs)))'
 'define r = 622*e/(lev-e)'
 'define condtemp = 55 + 1 / (1/(tmpprs-55) - log(rhprs/100)/2840)'
 'define thetaDL =
tmpprs*pow(1000/(lev-e),0.2854)*pow(tmpprs/condtemp,0.0028*r)'
 'define thetae =
thetaDL*exp(((3.036/condtemp)-0.00178)*r*(1+0.448+r/1000))'
 'define thetaw = 45.114 - 51.489*pow(thetae/273.15,-3.504)'
 'define LSI = max(thetaw,lev=1000,lev=500) - max(thetaw,lev=1000,lev=900)'

I get a lot of NaNs.  They are coming from this line:
 'define thetaDL =
tmpprs*pow(1000/(lev-e),0.2854)*pow(tmpprs/condtemp,0.0028*r)'

And in fact, from this piece:
pow(1000/(lev-e),0.2854)

I know this because I tested this for a model grid column (from 1000 - 500
mb).
ga-> d 1000/(lev-e)
Printing Grid -- 13 Values -- Undef = -9.99e+08
-3.55595 -4.04332 -4.69364 -9.98891 39.3747 -5.3341 -11.5849 4.19346
2.07427 2.42358 2.11214 2.1287 2.24423
ga-> d pow(1000/(lev-e),0.2854)
Printing Grid -- 13 Values -- Undef = -9.99e+08
-nan -nan -nan -nan 2.85282 -nan -nan 1.5055
1.2315 1.28743 1.23787 1.24063 1.25949

For whatever reason, the pow() function doesn't like negative bases,
because for example, I can use a calculator to compute -3.55595^0.2854 as
-1.43629.  Why, then, does the pow() function give a -nan as a result?

Jeff Duda

-- 
Jeff Duda
Graduate research assistant
University of Oklahoma School of Meteorology
Center for Analysis and Prediction of Storms
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://gradsusr.org/pipermail/gradsusr/attachments/20130118/c6a63449/attachment-0003.html