Issue 165: Added calculate_deltaT function into spa.py #240
Conversation
|
Thanks for the contribution, @VolkrB. It appears to me that this code only works with scalar input and that it would not work with array or series input. However, I can't actually tell for certain because there are no tests. So, can you please also add some unit tests for this code? I recommend using numpy's where function to support both scalar and array conditionals. ref #165. |
|
Ok I'll work on that. Thanks! |
|
I worked on this some today, it can now take scalar, series and array input, but I just have a few questions, if you don't mind, before I begin working on the unit tests:
Thanks, (also this was my first time playing around with the numpy where function. Pretty cool!) |
|
Sounds great, though you'd need to push the changes if you want us to see them at this point.
Are you planning on integrating this function with the spa calculators or leaving it to the user to separately call the function and then supply the result to a spa calculator? Automatic would be nice, but only if it doesn't significantly slow down the computation. It's also slightly complicated to make it work automatically with the numba-accelerated code. We can certainly just add the function for now and address the integration later on, too. |
|
I pushed some changes, but first I'll address your comments:
My changes:
|
|
Interesting approach, though I think we'll want to do something else. numpy.vectorize is basically short hand for a for loop, so I'm guessing that the code will be unusably slow for people that don't use numba (which is most people). How long does it take to do a calculation with 1 minute resolution for 1 year? If your code only makes the entire ephemeris calculation a tiny bit slower, then it's ok. We also really want to avoid duplicating code like this. Numba supports np.where, so I think you can remove the second function if you change your if statements to np.where statements. The tests look good, though please use the python-standard underscore_case instead of camelCase. |
|
Sounds good, thanks for the feedback! I'll start working on this soon. I am curious, why won't most people use numba? |
|
Yes, it appears that the (Month + 0.5) is a bug, it should be (Month - 0.5). I'll fix that in the MATLAB version. Thanks! |
|
Would it be worth making some kind of brief assessment of the impact of On 2016-09-26 20:59, DanRiley wrote:
Photovoltaic Performance Labs Germany +49-761-8973-5603 (Europe) |
|
@VolkrB I don't know of any other libraries that require or are even enhanced by numba, so I think it's fairly unknown and maybe untrusted. I think it's mostly still used by people for their own code. It's also hard to install if you're not using conda (at least it used to be), and not everyone uses conda. So, I think sticking with our current approach to numba is a good idea: use numba to make it easy to multithread code that's already reasonably fast. I'd be happy to reconsider that if there's a problem that would benefit from a different approach, but I don't think this one qualifies. @DanRiley thanks for confirming. @adriesse I don't think the deltat calculation is important for most applications, and I'm guessing that the importance of this error in the deltat calculation is even smaller. I'd still be interested in seeing it quantified, though! |
|
The maximum difference between the correct and incorrect version over the years -1999 to 3000 is about 1.5 seconds (at year -500). For the period of importance for PV (1960 to 2050), the error in delta T is 0.1 seconds. The error in sun hour angle is roughly 1/240th of a degree per second, so we're looking at hour angle errors of 1/2400th of a degree for the period of importance for most PV applications. |
|
Thanks! That puts it into perspective! On 2016-09-26 23:49, DanRiley wrote:
Photovoltaic Performance Labs Germany +49-761-8973-5603 (Europe) |
|
Took a step back and made the following changes from my original pull request:
for a 1E6 elements array: |
|
Thanks, I like the direction. It's probably more readable if you make a new deltat assignment for each where statement. The first where call needs to set the default value. deltat = np.where(cond1, values1, year*0)
deltat = np.where(cond2, values2, deltat)
...Not sure if this will be a little slower or a lot slower. For context, how long does the get_solarposition with method='nrel_numpy' take on your machine for a similar number of points? |
|
Back on track! I''ll work on making it more readable and testing soon. For now, here's the %timeit for my potato: |
|
Adding deltat assignments didn't slow it down at all, it was actually 20ms faster for a 1E6 element array, so I went ahead and pushed the change |
There was a problem hiding this comment.
Great, this is looking good. Are you going to hook this up to the solar position calculators in this PR or later? You'll need more tests when you do that.
Can you run flake8 on the module? It will complain about line lengths, spaces, etc. so that I don't have to. Don't worry about fixing things that you didn't touch.
|
|
||
| -20+32*((y-1820)/100)**2, deltat) | ||
|
|
||
| deltat = np.asscalar(deltat) if np.isscalar(year) & np.isscalar(month) else deltat |
There was a problem hiding this comment.
Is this line necessary? I thought that np.where would preserve the type.
There was a problem hiding this comment.
I think the output is always either ndarray or tuple of ndarrays.
There was a problem hiding this comment.
hm, I was mistaken. I am guessing that pvlib has a test or two that unintentionally passes because you can assert e.g. 1 == np.array(1) or np.testing.assert_allclose(1., np.array(1.)) without an exception.
|
|
Should Also, should the year and month variables input be called I guess if the |
|
It seems to me like most of the tests in solarposition can be added by simply removing the delta_t keyword. Is this what you imagined for the testing? I pushed an example of what i mean in test_solarposition |
|
Probably make a new test or tests with _no_deltat in the name so that it's more clear what we're testing. Or you can try to get more clever with pytest's parametrize if you want to avoid copy paste. |
|
I wrote the test using pytest's parametrize decorator for the get_solarposition function. Let me know if this is what you imagined. I was debating whether the test should be put with the spa functions and/or I should add more test cases. |
|
Works for me. You could add another test for a time at which the previous default would fail, but I don't know if that's really necessary. The failing tests in irradiance.py will pass if you add an atol=1e-4 to the assert statements. |
|
Since deltat can now be numpy array, there is an issue with the numba compiled solar_position function.
gives ValueError. This makes sense, because delta_t is now an array. However, the fix for this seems not so straight forward to me. |
|
Yeah, that's a little tricky. You'll need to remove delta_t from loc_args, chunk it with the input data and output data, and make it an argument of solar_position_loop. Not tested, but something like this should work... @jcompile('void(float64[:], float64[:], float64[:], float64[:,:])', nopython=True,
nogil=True)
def solar_position_loop(unixtime, delta_t, loc_args, out):
"""Loop through the time array and calculate the solar position"""
lat = loc_args[0]
lon = loc_args[1]
elev = loc_args[2]
pressure = loc_args[3]
temp = loc_args[4]
atmos_refract = loc_args[5]
sst = loc_args[6]
esd = loc_args[7]
...
def solar_position_numba...
loc_args = np.array([lat, lon, elev, pressure, temp,
atmos_refract, sst, esd])
...
# split the input and output arrays into numthreads chunks
split0 = np.array_split(unixtime, numthreads)
split2 = np.array_split(result, numthreads, axis=1)
if np.isscalar(delta_t):
chunks = [[a0, delta_t, loc_args, split2[i]] for i, a0 in enumerate(split0)]
else:
chunks = [[a0, delta_t[i], loc_args, split2[i]] for i, a0 in enumerate(split0)]The docstring for spa.solar_position should also be updated. |
|
I forgot that you originally wanted to make the delta_t calculation compatible with numba. The approach I outlined in my previous comment doesn't quite allow for that. Instead, you could add numba decorators to the calculate_deltaT function and call this function within the loop of solar_position_loop and within solar_position_numpy. You'd need to enclose those calls in an if statement that checked for if the input delta_t is None. I'd also consider adding an xfail to the numba test, merging the PR, marking it as a known incompatibility, and addressing it in a separate PR. |
|
I've added a docstring to the |
|
I realized that this will break anyone's code (including mine!) that uses
I recommend 2, with the appropriate notes in the documentation. Please also add a note to the 0.4.2 whatsnew document (you might need to first rebase/merge the latest commits) Sorry for dragging this out and thanks for your patience! |
|
No problem! I'll work on getting this done asap. |
|
I've sprinkled some docs, added a note in what's new, and replaced Not sure what's up with test_forecast |
|
|
||
|
|
||
| Enhancements | ||
| ~~~~~~~~~~~~ | ||
|
|
||
| * Adding a complete_irradiance method to the ModelChain to make it possible to calculate missing irradiation data from the existing columns [beta] (:issue:`239`) | ||
|
|
||
| * Added calculate_deltat method to the spa module to calculate the time difference between terrestrial time and UT1. Specfiying a scalar is sufficient for most calculations. (:issue:`165`) |
| @@ -61,6 +61,7 @@ def test_extraradiation(input, expected, method): | |||
|
|
|||
|
|
|||
| @requires_numba | |||
| @pytest.mark.xfail(raises=ValueError, reason = 'spa.calculate_deltat not implemented for numba yet') | |||
There was a problem hiding this comment.
We can now get rid of the xfail on test_extraradiation_nrel_numba, right? It would be good to add a more localized test in test_spa.py or test_solarposition.py that does fail under the expected conditions (decorated with xfail, of course).
There was a problem hiding this comment.
you are correct, don't need that xfail anymore
|
My approach is as follows:
|
| "delta_t, method, expected", [ | ||
| (None, 'nrel_numpy', expected_solpos_multi()), | ||
| (67.0, 'nrel_numpy', expected_solpos_multi()), | ||
| (None, 'nrel_numba', expected_solpos_multi()), |
There was a problem hiding this comment.
Nice. Only thing is that you need to decorate just this one parameter set so that we don't get XPASS results for the others. I think there's an example in the pytest docs.
There was a problem hiding this comment.
yea i was wondering about that. thanks for the heads up
|
This should fix it. |
|
Thanks @VolkrB! |
|
No problem! I am going to start looking at making the numba calculation working. |
Hey!
I'd like to get involved and figured this would be a good place for me to start. I have appended a delta t function (similar to the one that was recently added into the matlab version) to the spa module. Anything you want me to change, add or improve, I will happily do!
Volker