MESAHub / mesa

Modules for Experiments in Stellar Astrophysics
http://mesastar.org
GNU Lesser General Public License v2.1
144 stars 38 forks source link

NaN from Skye EOS #86

Closed jschwab closed 4 years ago

jschwab commented 4 years ago

I saw the following NaN come out of the eos during a run with r14459. Looks like Skye conditions.

                                    solver_call_number         322
                                  gold tolerances level           2
                       correction tolerances: norm, max    3.0000000000000001D-05    3.0000000000000001D-03
                    tol1 residual tolerances: norm, max    9.9999999999999994D-12    1.0000000000000001D-09
                                               s% cp(k)        1893                       NaN
                                           s% csound(k)        1893                       NaN
                                              s% lnP(k)        1893                       NaN
                                              s% gam(k)        1893    9.8556840800187140D-01
                                                s% P(k)        1893                       NaN
                                             s% Pgas(k)        1893                       NaN
                                              s% rho(k)        1893    1.5317863773612282D+06
                                                s% T(k)        1893    6.0395571548145366D+08
                                                 logRho        1893    6.1851982028931074D+00
                                                   logT        1893    8.7810050955298653D+00
                                                   abar        1893    1.4883701474075236D+01
                                                   zbar        1893    7.4332806659512825D+00

                                           xa(j,k) neut           1           5        1893    1.0157613202179674D-04
                                             xa(j,k) h1           2           6        1893    1.2770953517588663D-21
                                             xa(j,k) h2           3           7        1893    2.3604348415364603D-30
                                            xa(j,k) he3           4           8        1893    2.0129928398075110D-34
                                            xa(j,k) he4           5           9        1893    7.8683159649280527D-15
                                            xa(j,k) li7           6          10        1893    2.7478137711801651D-23
                                            xa(j,k) be7           7          11        1893    5.6793721690784418D-39
                                             xa(j,k) b8           8          12        1893    3.0391754977941294D-62
                                            xa(j,k) c12           9          13        1893    2.3588720496100168D-01
                                            xa(j,k) c13          10          14        1893    1.9286328335483586D-05
                                            xa(j,k) n13          11          15        1893    3.8705961435540981D-11
                                            xa(j,k) n14          12          16        1893    8.3052223263579574D-07
                                            xa(j,k) n15          13          17        1893    1.3982164869287179D-09
                                            xa(j,k) o16          14          18        1893    7.4534628207214471D-01
                                            xa(j,k) o17          15          19        1893    1.4269407278953196D-05
                                            xa(j,k) o18          16          20        1893    4.9072403538797979D-08
                                            xa(j,k) f19          17          21        1893    7.7702849661993753D-11
                                           xa(j,k) ne20          18          22        1893    4.1923538660657463D-03
                                           xa(j,k) ne21          19          23        1893    8.7683713724065543D-06
                                           xa(j,k) ne22          20          24        1893    9.4293397671965686D-03
                                           xa(j,k) na21          21          25        1893    2.6164461245748743D-17
                                           xa(j,k) na22          22          26        1893    2.0418185326468750D-11
                                           xa(j,k) na23          23          27        1893    2.5198430502819408D-03
                                           xa(j,k) mg24          24          28        1893    4.2891444212455319D-04
                                           xa(j,k) mg25          25          29        1893    3.5174721789680607D-06
                                           xa(j,k) mg26          26          30        1893    1.1351776483079530D-04
                                           xa(j,k) al25          27          31        1893    8.4478859234650238D-16
                                           xa(j,k) al26          28          32        1893    6.5594972501604542D-15
                                           xa(j,k) al27          29          33        1893    1.9342452354714360D-03
                                              sum(xa) =    1.0000000000000002D+00
                                                Pgas =                        NaN
                                             logPgas =                        NaN
                                                 rho =     1.5317863773612282D+06
                                                   T =     6.0395571548145366D+08
                                                logQ =     6.2318801183337769D-01

                                                   P =                        NaN
                                                Prad =     3.3554452152910034D+20
                                                logP =                        NaN
                                                logS =                        NaN
                                                logE =                        NaN
                                              energy =                        NaN
                                               grada =                        NaN
                                                  cv =                        NaN
                                                  cp =                        NaN
                                              gamma1 =                        NaN
                                              gamma3 =                        NaN
                                                 eta =     3.1282095632650346D+00
                                                 gam =     9.8556840800187140D-01
                                                  mu =     1.7648768093497738D+00
                                          log_free_e =    -3.0153034496077669D-01
                                              chiRho =                        NaN
                                                chiT =                        NaN
                  s% eos_rq% logRho_min_OPAL_SCVH_limit   -1.4298999999999999D+01

 do_eos_for_cell k, nz        1893        1918
 lnPgas_flag F
                                              logRho =     6.1851982028931074D+00
                                                logT =     8.7810050955298653D+00
                                                   z =     9.9999999999999212D-01
                                                   x =     1.2770953541193011D-21
                                                abar =     1.4883701474075236D+01
                                                zbar =     7.4332806659512825D+00

 s% eos_rq% tiny_fuzz   9.9999999999999995E-007
 s% eos_rq% logRho_min_for_all_OPAL  -9.0999999999999996     
 s% eos_rq% logRho_min_for_any_OPAL  -9.1999999999999993     
 s% eos_rq% logQ_max_OPAL_SCVH   5.2999999999999998     
adamjermyn commented 4 years ago

I believe this is caused by the solid phase mixing corrections, which produce NaN whenever any species with abundance above 1d-7 has zero charge (neutrons). This issue was in PC too but I think neutrons usually came up in matter that's warm enough to be liquid. In Skye we need to be able to evaluate the solid phase corrections at all temperatures and densities for phase determination purposes, so this is a problem.

To resolve this I've told the solid mixing routines to skip species with zero charge.

@jschwab Can you confirm and close the issue if this fixes it for you?

jschwab commented 4 years ago

Thanks Adam! I confirm this resolves the issue.