Add associated Legendre function calculations #505
Conversation
Allocate the output array out of the loop. This way we can do it once only when evaluating the spherical harmonics. Pre-compute the square roots of integers that we use in the loops. Doesn't work well if the square roots are calculated at the module level (probably a numba thing).
Results are bad beyond that even with the scaling
|
@santisoler and @MarkWieczorek this is ready now. I added the Holmes and Featherstone scaling and tested until degree 2800. I was surprised that calculating the derivatives on rescaled functions works well. That's nice because it means we can keep them separate and use the derivative functions to calculate higher order derivatives (which would appear in gravity gradients). If you could both have a quick look I'd really appreciate it. And thanks for all the information and code examples @MarkWieczorek! It really helped. |
There was a problem hiding this comment.
Great work @leouieda! I think it will be great to have the spherical harmonics in Harmonica (redundancy required 😁), and this is step forward that.
I left some minor comments here. Let me know what do you think.
|
Thanks, @santisoler! What do you think? |
|
Looking great @leouieda! I'm updating the branch so we can merge it. |
The numba-based functions calculate the associated Legendre functions (unnormalized, Schmidt normalized, and fully normalized) and their derivatives$\dfrac{\partial P_n^m}{\partial \theta}(\cos\theta)$ . The derivatives functions can be used to calculate higher order derivatives as well. Values are tested against analytical solutions for the first 4 degrees as well as an identity for the Schmidt functions and the Legendre equation for the Schmidt and fully normalized functions. The values are returned as 2D numpy arrays. This is a bit of a waste of memory but it's much easier to use than storing these in sparse matrices or 1D arrays. But it should be acceptable since we're unlikely to ever go over ~2700 degree so the waste is in the order of a few Mb at worse.
Relevant issues/PRs:
Related to #504