[Mesa-users] Warnings with additional energy sink and rel_E_err
Adrian Ayala
aayala at iaa.es
Tue May 12 15:21:10 EDT 2020
On 12/05/2020 19:43, Adrian Ayala wrote:
>
> Hi,
>
> I tried to compute a 1 solar mass model at r12778, implementing an
> extension in run_star_extras, to account for an additional sink in the
> star energy budget.
>
> I got a value of rel_E_err of 10^-1, which is, in my views, too much.
>
> The point is that the error remains, which the same order of
> magnitude, even when I turned off the additional energy sink, running
> a simple 1 solar mass, solar metallicity model.
>
> Please, find the inlists attached.
>
> Best regards,
>
> Adrián
>
> On 12/05/2020 17:22, Bill Paxton wrote:
>>
>>> On May 12, 2020, at 1:24 AM, Ian Foley <ifoley2008 at gmail.com
>>> <mailto:ifoley2008 at gmail.com>> wrote:
>>>
>>> Hi Bill and Josiah,
>>>
>>> Working through your suggestions.
>>>
>>> Two runs with gold tolerances
>>>
>>> 1M2gold.txt - stopped after 2000 models. Has nearly 400 retries and
>>> 11 backups but still a long way from reaching the main sequence. My
>>> guess would be that I'll need 10 million models to reach a wd.
>>>
>>> 1M2goldmod.txt. Modified gold
>>>
>>> gold_tol_residual_norm1 = 1d-11
>>> gold_tol_max_residual1 = 1d-9
>>> gold_iter_for_resid_tol2 = 5
>>> gold_tol_residual_norm2 = 1d-8
>>> gold_tol_max_residual2 = 1d-6
>>> gold_iter_for_resid_tol3 = 10
>>> gold_tol_residual_norm3 = 1d99 ! default 1d-6
>>> gold_tol_max_residual3 = 1d99 ! default 1d-4
>>> gold_newton_iterations_limit = 11 ! default 14 ! for timestep
>>> control
>>>
>>> In this case, (if I've understood correctly), if beyond 10
>>> iterations a max_correction will be applied to provide a solution.
>>> So time step is not impacted so the previous performance is
>>> retained. We hit the man sequence after about 300 models. Problem is
>>> that there is one bad lg_LH at model 299. Most of the time the
>>> correction applied does not have a significant impact on the energy
>>> error. Well within the default 0.02 which
>>> triggers the default warning message.
>>
>> ok. gold_tolerances cause more iterations in an attempt to get a
>> better residual. but as you can see the default settings allow the
>> solver to give up on residuals after 10 iterations and just consider
>> the size of corrections. so you must have a case where the
>> corrections are tiny, but they cannot yield small residuals. can
>> you think of why that might happen? consider finding x such that
>> f(x) = 0 where f is making a nearly discontinuous jump from large
>> negative to large positive. tiny changes in x around the root lead
>> to large jumps in f(x). nonlinearity is nasty for a root solver
>> that depends on nicely behaved functions with not-too-extreme jumps.
>> sometimes the extreme nonlinearity is physical; but it can also be
>> a numerical problem. that might turn out to be the issue here since
>> it shows up when using rates that aren’t often included.
>>>
>>> Applying
>>> delta_lgL_H_limit = 0.01d0
>>> delta_lgL_H_hard_limit = 0.05d0
>>>
>>> makes no difference. Still get exactly the same value for lg_LH at
>>> the same model number. The large change in lg_LH does not generate
>>> retries.
>>
>> ? hmm. that is a problem. changes in lg_LH greater than the hard
>> limit should force a retry. if that isn’t happening, please save a
>> model just before an example. send that along with the modified
>> inlists so i can reproduce the problem easily.
>>>
>>> My suspicion is that the error occurs when a max_correction is being
>>> applied to get an accepted solution. And usually the max_correction
>>> works fine so we get reasonable solutions without getting a time
>>> step adjustment and without getting too big an energy error.
>>>
>>> btw, I do have 2 runs with the o18_and_ne22.net
>>> <http://o18_and_ne22.net/> with a mesh coeff of 1.0 for 1M stars
>>> The first, using gold tolerances, took 12033 models to wd and
>>> finished with log_rel_run_E_err = -3.5875688657807161
>>> The second, without gold tolerances, took 11991 models to wd and
>>> finished with log_rel_run_E_err = -3.7343232181208545
>>> This made me question whether gold tolerances really are the way to go.
>>
>> gold tolerances can be helpful in many cases, but they are certainly
>> not going to be a miracle cure for all bugs! the experiment of
>> running with them shows that this problem can’t be fixed by having
>> the newton solver just try harder. that’s valuable information.
>> next comes the hard work to find the source of the problem.
>>
>> b
>>
>>
>>>
>>> I will continue to work on the other suggestions.
>>>
>>> Kind regards
>>> ian
>>>
>>>
>>>
>>> On Tue, 12 May 2020 at 14:20, Bill Paxton <paxton at kitp.ucsb.edu
>>> <mailto:paxton at kitp.ucsb.edu>> wrote:
>>>
>>>
>>>> On May 11, 2020, at 8:40 PM, Ian Foley via Mesa-users
>>>> <mesa-users at lists.mesastar.org
>>>> <mailto:mesa-users at lists.mesastar.org>> wrote:
>>>>
>>>> Ok, I will try to do that. Thanks
>>>>
>>>> Personally, I was examining the accumulation of energy energy
>>>> errors and was puzzled by the occasional huge rel_E_err.
>>>
>>> You are right to be puzzled - in fact, you should be horrified!
>>> The system has definitely gone astray at step 274. The
>>> following are from the log you sent (thank you for including
>>> that - it’s a big help).
>>>
>>> step lg_Tcntr Teff lg_LH lg_Lnuc Mass
>>> H_rich H_cntr N_cntr Y_surf X_avg eta_cntr
>>> zones retry
>>> lg_dt_yr lg_Dcntr lg_R lg_L3a lg_Lneu
>>> lg_Mdot He_core He_cntr O_cntr Z_surf Y_avg
>>> gam_cntr iters bckup
>>> age_yr lg_Pcntr lg_L lg_gsurf lg_Psurf
>>> lg_Dsurf C_core C_cntr Ne_cntr Z_cntr Z_avg
>>> v_div_cs dt_limit
>>> __________________________________________________________________________________________________________________________________________________
>>>
>>> 273 7.149445 5660.570 -0.102621 -0.102620 0.999678
>>> 0.999678 0.588152 0.005009 0.280000 0.687990 -1.638059
>>> 421 2
>>> 7.544841 1.992413 -0.033735 -46.026721 -1.760965
>>> -12.689535 0.000000 0.391282 0.009290 0.020000 0.291896
>>> 0.094468 6 0
>>> 1.6342E+09 17.232539 -0.102547 4.505219 4.819642
>>> -6.735807 0.000000 0.000016 0.001942 2.057E-02 2.011E-02
>>> 0.000E+00 retry
>>> rel_E_err -1.1024069047620799D-10
>>> log_rel_run_E_err -2.6819105959472362
>>> rel_run_E_err -0.0020801248584473
>>>
>>> 274 7.149838 5661.517 0.549531 0.549531 0.999671
>>> 0.999671 0.585504 0.005009 0.280000 0.687719 -1.636798
>>> 421 2
>>> 7.544841 1.994223 -0.033378 -45.996348 -1.759985
>>> -12.688180 0.000000 0.393930 0.009288 0.020000 0.292167
>>> 0.094680 6 0
>>> 1.6692E+09 17.233776 -0.101543 4.504503 4.819064
>>> -6.736461 0.000000 0.000016 0.001942 2.057E-02 2.011E-02
>>> 0.000E+00 hold
>>> rel_E_err -3.8973049465573766
>>> log_rel_run_E_err 0.5905325552405367
>>> rel_run_E_err 3.8952250568823850
>>>
>>> rel_E_err of 10^-10 is fine, so step 273 is ok. but 274 is
>>> terrible. let’s change the inlist controls so that it will not
>>> allow 274 to be accepted with that large rel_E_err.
>>>
>>> one thing to do is turn use_gold_tolerances back on. (i noticed
>>> that you’ve turned it off in your inlist for some reason.)
>>> with gold_tolerances enabled, the newton solver will keep trying
>>> for more iterations to get a better solution with smaller
>>> residuals, and that might be enough to drive down the rel_E_err.
>>>
>>> also, you should try using the timestep limits for rel_E_err in
>>> order to force a retry when a step is accepted but gives a large
>>> rel_E_err. try something like the following:
>>>
>>> limit_for_abs_rel_E_err = 1d-7
>>> hard_limit_for_abs_rel_E_err = 1d-4
>>>
>>> also note that the lg_LH jumped in a single step from -0.1 to
>>> 0.5. you should try limiting the size of that change also.
>>> You have a timestep limit for lgL_He. Add one for lgL_H. And
>>> add hard limits too to force retries with smaller timesteps if
>>> get really large jumps.
>>>
>>> delta_lgL_H_limit = 0.01d0
>>> delta_lgL_H_hard_limit = 0.05d0
>>>
>>> Please try these changes and let me know what happens. You
>>> should either get a nice clean run keeping rel_E_err and
>>> log_rel_run_E_err small, or it will die with tiny timesteps. If
>>> that bad result happens, then it will be time to dig deeper to
>>> see where the problem is arising. It might be related to the
>>> larger nets you are using - especially if the run is ok with the
>>> smaller more commonly used ones. So you should certainly check
>>> if the problem goes away with a smaller net such as
>>> pp_and_cno_extras. I’m not saying you shouldn’t be able to
>>> use the larger nets; just that they are not used very much and
>>> there may be bugs lurking in them that no one has reported
>>> before now (e.g., rarely used reactions may be producing bad
>>> rates and/or partials of the rates with respect to T or Rho -
>>> either of those could make it impossible for the newton solver
>>> to get residuals down).
>>>
>>> Thanks for reporting these puzzles. Not all of them will get
>>> satisfactory solutions, but looking at them may reveal things
>>> that we can do to improve mesa. So I appreciate your feedback.
>>>
>>> Cheers,
>>> Bill
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>> I was trying to track the cause thinking it was related to the
>>>> solution of the equations in star_newton so I'm not sure if the
>>>> problem occurs from the lg_LH value or from the equation
>>>> solution. Normally, I only dump the terminal every 10 models,
>>>> so could easily have missed seeing the funny lg_LH values.
>>>>
>>>> Kind regards
>>>> Ian
>>>>
>>>>
>>>>
>>>> On Tue, 12 May 2020 at 13:20, Josiah Schwab <jwschwab at ucsc.edu
>>>> <mailto:jwschwab at ucsc.edu>> wrote:
>>>>
>>>> Hi Ian,
>>>>
>>>> Thanks. I can reproduce this. We'll look into it and get
>>>> back to you.
>>>>
>>>> In the future, it would be helpful if you could take the
>>>> time to reduce the problem further. It is easiest if you
>>>> send a minimal work directory that we can compile and run
>>>> and immediately see the problem. For example, save a model
>>>> 5 or 10 steps before the problem, start the run from that,
>>>> and only include the inlists needed for that model. I've
>>>> already done this here, but keep that in mind.
>>>>
>>>> Josiah
>>>>
>>>> _______________________________________________
>>>> mesa-users at lists.mesastar.org
>>>> <mailto:mesa-users at lists.mesastar.org>
>>>> https://lists.mesastar.org/mailman/listinfo/mesa-users
>>>>
>>>
>>> <1M2goldmod_299.txt><1M2gold.txt>
>>
>>
>> _______________________________________________
>> mesa-users at lists.mesastar.org
>> https://lists.mesastar.org/mailman/listinfo/mesa-users
>>
>
> _______________________________________________
> mesa-users at lists.mesastar.org
> https://lists.mesastar.org/mailman/listinfo/mesa-users
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20200512/ef7f882a/attachment.htm>
-------------- next part --------------
! inlist to evolve a 0.82 solar mass star, Reimers 0.4 and Y 0.248
! For the sake of future readers of this file (yourself included),
! ONLY include the controls you are actually using. DO NOT include
! all of the other controls that simply have their default values.
&star_job
! begin with a pre-main sequence model
create_pre_main_sequence_model = .true.
! save a model at the end of the run
save_model_when_terminate = .false.
save_model_filename = 'M1.22solar_Reimers04_Y0273_HP_Chi0'
! display on-screen plots
pgstar_flag = .true.
! Set of Castelli 2003 bolometric corrections in UBV (Johnson), RI (Cousins) and JHKL (Bessell)
! Can be used at the same time as the lcb98cor.dat file
! Filter names BC_v U-B B-V V-I V-R V-J V-K V-L H-K
color_num_files=2
color_file_names(2)='vmag_tables_Adrian.dat'
color_num_colors(2)=9
/ !end of star_job namelist
&controls
! starting specifications
initial_mass = 1.0d0
initial_y = 0.2676d0
initial_z = 0.0125d0
! use_Type2_opacities = .true.
! Zbase = 0.001
! Use neutrino hook
use_other_neu = .true.
! Activate semi convection while mixing
use_Ledoux_criterion = .true.
alpha_semiconvection = 1
semiconvection_option = 'Langer_85'
thermohaline_coeff = 1
thermohaline_option = 'Brown_Garaud_Stellmach_13'
! mass gain or loss
mdot_omega_power = 0.43
cool_wind_RGB_scheme = 'Reimers'
Reimers_scaling_factor = 0.2d0
cool_wind_AGB_scheme = 'Blocker'
Blocker_scaling_factor = 0.2d0
RGB_to_AGB_wind_switch = 1d-4
! cool_wind_full_on_T = 8d3
! cool_wind_full_off_T = 1d4
! relax_mass = .true. ! gradually change total mass by a wind ! Comentado porque la actual version no reconoce el parametro
! relax_initial_mass = .true. ! Comentado porque la actual version no reconoce el parametro
! MLT
mixing_length_alpha = 1.9658
MLT_option = 'Henyey'
! atmosphere
atm_table = 'photosphere'
! predictive mixing (to avoid breathing pulses and allow core expansion)
predictive_mix(1) = .true.
predictive_zone_type(1) = 'any'
predictive_zone_loc(1) = 'any'
predictive_bdy_loc(1) = 'any'
predictive_superad_thresh(1) = 0.01
predictive_avoid_reversal(1) = 'he4'
! stop passed horizontal branch
HB_limit = 0.0000001
write_pulse_data_with_profile = .true.
pulse_data_format = 'GYRE'
/ ! end of controls namelist
-------------- next part --------------
&pgstar
! MESA uses PGPLOT for live plotting and gives the user a tremendous
! amount of control of the presentation of the information.
! show HR diagram
! this plots the history of L,Teff over many timesteps
HR_win_flag = .true.
! set static plot bounds
HR_logT_min = 3.0
HR_logT_max = 4.6
HR_logL_min = -1.0
HR_logL_max = 4.0
! set window size (aspect_ratio = height/width)
HR_win_width = 6
HR_win_aspect_ratio = 1.0
! show temperature/density profile
! this plots the internal structure at single timestep
TRho_Profile_win_flag = .true.
! add legend explaining colors
show_TRho_Profile_legend = .true.
! display numerical info about the star
show_TRho_Profile_text_info = .true.
! set window size (aspect_ratio = height/width)
TRho_Profile_win_width = 8
TRho_Profile_win_aspect_ratio = 0.75
Mode_Prop_win_flag = .true.
Mode_Prop_win_width = 6
Mode_Prop_win_aspect_ratio = 0.75
Mode_Prop_xleft = 0.15
Mode_Prop_xright = 0.85
Mode_Prop_ybot = 0.15
Mode_Prop_ytop = 0.85
Mode_Prop_txt_scale = 1.0
Mode_Prop_title = 'Mode_Prop'
Mode_Prop_nu_max_obs = -999
/ ! end of pgstar namelist
-------------- next part --------------
! this is the master inlist that MESA reads when it starts.
! This file tells MESA to go look elsewhere for its configuration
! info. This makes changing between different inlists easier, by
! allowing you to easily change the name of the file that gets read.
&star_job
read_extra_star_job_inlist1 = .true.
extra_star_job_inlist1_name = 'inlist_project'
/ ! end of star_job namelist
&controls
read_extra_controls_inlist1 = .true.
extra_controls_inlist1_name = 'inlist_project'
/ ! end of controls namelist
&pgstar
read_extra_pgstar_inlist1 = .true.
extra_pgstar_inlist1_name = 'inlist_pgstar'
/ ! end of pgstar namelist
More information about the Mesa-users
mailing list