Degree-1 Love Number Issues

[nlwagner]

Opening this bug report to document exploring why degree-1 Love numbers are often too large. Will add reports on test cases to see if this is an integration issue or not.

[jrenaud90]

Awesome, thanks for working on this.

Here is where the boundary condition is defined:

TidalPy/TidalPy/RadialSolver/solver.pyx

Line 306 in edab615

while num_ytypes > ytype_i:

And here is where it is applied:

TidalPy/TidalPy/RadialSolver/boundaries/boundaries.pyx

Line 47 in edab615

surface_matrix_ptr = &surface_matrix_solid[0][0]

You could make modifications to the first link if it is as simple as a factor off. The 3 BCs that are defined there are applied to y1, y3, and y6. If you need to apply a BC to a different y (maybe y5 as in LoadDef's table 5.2) then you will need to hack the code in the 2nd link. Notice the [0 * max_num_y + 1], [0 * max_num_y + 3], [0 * max_num_y + 5], etc? Those are pulling out the values of y2, y4, y6 respectively. So applying a BC to the other ys would require changing that index.

Note any changes you make will need to be recompiled by doing a python -m pip install -v . on your local, modified version of TidalPy.

[jrenaud90]

Putting a note here to remind us that if the trick of subtracting off k_1 from the other love numbers is indeed the correct approach that it should be easy to implement as a post-process. So that is a ToDo if we go that route.

[nlwagner]

truthfully through reading LoadDef's manual more and testing I think the k_1 trick is enough to get more accurate h1s. However! Ive found something kind of fascinating. In the attached plot I have h1 as a function of forcing frequency. (note these are all Load Love numbers and elastic)

I figured this sort of insight might be useful for you. The slight decrease in h1 after about 100 hours is strange to me though, do you have any insight? k2 and h2 don't have this decrease around the same value:

[nlwagner]

Note: I was getting stable h1s with a simpler interior model. Increasing complexity makes h1 want to misbehave. I'll work on quantifying this more in the future.

[nlwagner]

Here's a more complete plot before I sign off for the day:

[jrenaud90]

That is really interesting. Thanks for sharing - would you be willing to share your inputs to radial solver? I can poke around to see if I can learn anything more.

My one guess is that you are hitting some sort of resonance with the rheology. A good check is if you change the viscosity (not sure if you are using a constant model or some $\eta(r)$ function. But either way if you take it and bump it up or down $eta(r)_n = 1000 \times \eta(r)_o$ (or divide by 1000). And see if that changes anything? I am surprised that would change l=1 and not l>1 though...

[nlwagner]

Oops sorry just saw this. Definitely. Once I push my workspace to my forked repository I'll let you know.

Additionally, I've noticed singularities occuring at really long periods (>1000 hours) that appear/disappear with varying integration tolerances. Seen below.

If we tighten the tolerances we see:

So probably no surprise we see numerical issues at really long periods I'm just trying to double check degree 1 isn't actually dependent on frequency in an elastic case :)

[nlwagner]

Within my forked repository "mars", the file that generated the above plots is deg1deg2_stability_complexmodel.py

Note I had to add 1500 slices for each layer so it takes a while to run.