Hi Rich and the rest of the GYRE community,
I am having some issues with non-adiabatic analysis of some high-Teff pre-WD models produced by MESA. I have produced pre-WD models with a range of masses and I have been probing the stability of these models, focusing on the radial fundamental mode ( l=0, n_p = 1).
My analysis works fine in most cases, however I am finding some issues at high effective temperatures. In some cases, while the adiabatic analysis finds a monotonic sequence of radial modes starting at n_p = 1, the non-adiabatic results fail to find the fundamental mode.
I include 3 models as examples. The '290p622' files relate to a 0.295 Msun model with a Teff of 39600 K where the and the '295p694' files relate to a 0.29 Msun with a Teff of 41900 K and the '410p431' files are from a 0.41 Msun model with a Teff of 63500 K.
Example 1: GYRE terminal output of '290p622' gives:
Root bracketing
Time elapsed : 4.502 s
Root Solving
l m n_pg n_p n_g Re(omega) Im(omega) chi n_iter
0 0 1 1 0 0.38398863E+01 0.00000000E+00 0.2357E-11 5
0 0 2 2 0 0.40989930E+01 0.00000000E+00 0.3482E-12 5
0 0 3 3 0 0.54916991E+01 0.00000000E+00 0.1214E-11 5
0 0 4 4 0 0.67685270E+01 0.00000000E+00 0.1466E-11 5
0 0 5 5 0 0.80640478E+01 0.00000000E+00 0.2751E-11 5
0 0 6 6 0 0.92936948E+01 0.00000000E+00 0.2538E-11 5
Time elapsed : 0.059 s
Starting search (non-adiabatic)
Root Solving
l m n_pg n_p n_g Re(omega) Im(omega) chi n_iter
0 0 1 1 0 0.39887494E+01 0.39071418E-01 0.1549E-13 8
0 0 2 2 0 0.42228618E+01 0.24874928E+00 0.2481E-13 9
0 0 3 3 0 0.55254721E+01 0.84734945E-03 0.3422E-12 9
0 0 4 4 0 0.67998720E+01 -0.52523593E-01 0.1926E-12 6
0 0 5 5 0 0.79904143E+01 -0.77464999E-01 0.1700E-12 8
0 0 6 6 0 0.89757654E+01 -0.11089162E+00 0.6503E-13 9
Time elapsed : 0.528 s
In this case everything seems to function quite well and the result is as I expected.
Example 2: When I move to an adjacent model in my grid with a reasonably comparable effective temperature, the terminal output of '295p694' gives:
Root bracketing
Time elapsed : 5.284 s
Root Solving
l m n_pg n_p n_g Re(omega) Im(omega) chi n_iter
0 0 1 1 0 0.36373966E+01 0.00000000E+00 0.3061E-11 6
0 0 2 2 0 0.41554992E+01 0.00000000E+00 0.1570E-11 5
0 0 3 3 0 0.56193145E+01 0.00000000E+00 0.4839E-12 5
0 0 4 4 0 0.68593011E+01 0.00000000E+00 0.5886E-12 4
0 0 5 5 0 0.81117628E+01 0.00000000E+00 0.1229E-11 4
0 0 6 6 0 0.93183111E+01 0.00000000E+00 0.2885E-11 5
Time elapsed : 0.063 s
Starting search (non-adiabatic)
Root Solving
l m n_pg n_p n_g Re(omega) Im(omega) chi n_iter
0 0 2 2 0 0.41426836E+01 0.25466404E-01 0.2217E-14 11
0 0 2 2 0 0.40676369E+01 0.42280368E+00 0.3338E-12 8
0 0 3 3 0 0.56488349E+01 0.19009108E-02 0.2359E-11 7
0 0 4 4 0 0.68724787E+01 -0.51255914E-01 0.5178E-12 7
0 0 5 5 0 0.79815857E+01 -0.81906322E-01 0.3438E-12 7
0 0 6 6 0 0.90179117E+01 -0.10594145E+00 0.1795E-12 10
Time elapsed : 0.556 s
In this second instance, the non-adiabatic analysis fails to find the fundamental mode, instead collapsing to twice find the first overtone. I have tried adding extra mesh points through the resampling parameters but see little/no change in the outcome. Looking at both stellar models, there doesn't appear to be a significant physical difference between them.
Example 3: I have not had time to manually review each of the analyses carried out but one other case that I pulled out of may dataset to analyse was '410p431'. Here the terminal output shows:
Root bracketing
Time elapsed : 4.423 s
Root Solving
l m n_pg n_p n_g Re(omega) Im(omega) chi n_iter
0 0 1 1 0 0.22680694E+01 0.00000000E+00 0.1002E-10 10
0 0 2 2 0 0.54627480E+01 0.00000000E+00 0.1978E-11 7
0 0 3 3 0 0.66095621E+01 0.00000000E+00 0.4970E-12 5
0 0 4 4 0 0.74153814E+01 0.00000000E+00 0.1856E-11 5
0 0 5 5 0 0.91974173E+01 0.00000000E+00 0.2463E-11 5
Time elapsed : 0.053 s
Starting search (non-adiabatic)
Root Solving
l m n_pg n_p n_g Re(omega) Im(omega) chi n_iter
0 0 2 3 1 0.54630543E+01 -0.11899732E-03 0.2052E-17 1>
0 0 2 3 1 0.54630543E+01 -0.11899732E-03 0.1073E-10 6
0 0 4 4 0 0.73842745E+01 -0.28920292E-02 0.4230E-14 8
0 0 4 4 0 0.73842745E+01 -0.28920292E-02 0.1357E-11 6
0 0 5 5 0 0.90493436E+01 0.65355728E-02 0.3007E-12 8
Time elapsed : 0.387 s
In this case, the first mode found has a mode with mixed p-mode and g-mode components (n_p = 3, n_g = 1), again failing to locate the fundamental radial mode.
I also include in my zip file 2 plots of the logTeff-log(g) plane to show a) all of the models in my grid and b) the models where non-adiabatic GYRE finds the fundamental mode (n_p = 1). This seems to be a particular issue with my high effective temperature models, as none of them have an n_p = 1 mode found in the non-adiabatic analysis.
I also attach an example GYRE inlist from my grid. When I was using the COLLOC_GL2 scheme, the first non adiabatic mode found in some of these problem cases mentioned above would have a negative frequency, so I preferred using MAGNUS for these models.
So I suppose my questions are:
1. Does anyone have ideas why I would be failing to locate the radial fundamental mode in these models?
2. Is this an issue of numerics or physics?
3. What would be causing some of these models to have a 'first' p-mode found which has a mixed g-mode component?
In the zip file I attach the MESA .data and .data.GYRE files, along with the GYRE summary.txt files for the three models mentioned above.
Many thanks in advance for any assistance you can provide,
Conor