[Mesa-users] Resolving TP-AGB

Thu Mar 26 12:17:43 UTC 2026

Hi Theo,
The late TP-AGB is a very tricky evolutionary stage, especially at higher
masses.
I'd recommend having a look at Rees&Izzard 2024
https://ui.adsabs.harvard.edu/abs/2024MNRAS.531.4033R/abstract , as they go
through the main instabilities TP-AGB models experience, and offer possible
solutions as mesa subroutines.
Given  the max residuals  at the surface and the envelope mass, you might
be encountering HRI .

However, their subroutines might not (and at higher masses certainly
won't)  solve your problem, merely delaying it to a later point in the
TP-AGB.
If your interest lies *in* the TP-AGB stage, the default approach is to
carry calculations on for as long as possible and go through
post-processing. If you're interested in the WD stage,  the only
workaround is to remove the envelope before the instabilities crash the
code.

Good luck!
Ana Antonini

Théo-Courty via Mesa-users <mesa-users at lists.mesastar.org> escreveu
(quinta, 26/03/2026 à(s) 08:52):

> Dear Mesa users,
>
>
> I’m trying to evolve stars with mass varying between 1-5 Msun from ZAMS to
> white dwarf cooling track with version r25.12.1. However, I encounter
> numerical difficulties during the late TP-AGB phase for a 5Msun model with
> Z=0.012 (remaining mass ~1.25 Msun). The first thermal pulses go smoothly,
> but the late evolution struggles and the timestep becomes small and goes
> into a loop of increase/decrease between -3 < log dt <  -2.5. The limiting
> factor for dt is labelled as max increase until a retry and loops like
> this. I tried to lower max_timestep_factor to 1.01 but it does not solve
> the issue. I observe a lot of retries during that part of evolution until
> crash due to min_timestep_limit.
>
>
> I am fairly new to mesa and not very experienced with debugging. Here are
> typical examples of error messages I get:
>
>
> avg resid  0.223E-06     max resid dv_dt     1  0.35367E-03 mix type xx000
> avg corr  0.474E-03      max corr lnd   451 -0.21100E-01 mix type 00000
> avg+max corr+resid -- give up
>
>
> hydro_mtx: logT too small       26764         407           3
> 6.7237179921904777D-01    5.3757388036508846D+00
>
>
> hydro_mtx: logRho too large       26738         421           5
> 5.6268674362592407D+01   -7.8503397605837941D+00
>
>
> I tried to include okay_to_reduce_gradT_excess = .true. which I found used
> in some inlists online, but with this parameter activated the star shows
> non-physical behavior (log L > 10^5-6 with few small timesteps) so I keep
> it to false.
>
>
> Please find below the enlists I use which are based on make_co_wd
> test_suite which I modified to include element diffusion during MS and
> thermal pulses to remove the enveloppe instead of remove_env routine.
> Thanks a lot!
>
>
> Best regards,
>
>
> Théo
>
>
> *Common inlist for shared parameters*
>
> &star_job
>
>       show_log_description_at_start = .false.
>
>
>
>
>       change_net = .true.
>
>       new_net_name = 'co_burn_extras.net'
>
>
>
>       ! age start at ZAMS
>
>       set_initial_age = .false.
>
>       set_initial_model_number = .false.
>
>
>       num_special_rate_factors = 2
>
>       reaction_for_special_factor(1) = 'r_c12_ag_o16'
>
>       special_rate_factor(1) = 1
>
>       filename_of_special_rate(1) = 'c12ag_deboer_sigma_0p0_2000_Tgrid.dat'
>
>
>       reaction_for_special_factor(2) = 'r_he4_he4_he4_to_c12'
>
>       special_rate_factor(2) = 1
>
>       filename_of_special_rate(2) = 'r_he4_he4_he4_to_c12_cf88.txt'
>
>
> / ! end of star_job namelist
>
>
> &eos
>
>
>
> / ! end of eos namelist
>
>
> &kap
>
>       Zbase = 1.2d-02
>
>
>       kap_file_prefix = 'gs98'
>
>       use_Type2_opacities = .true.
>
>
> / ! end of kap namelist
>
>
> &controls
>
>
> initial_mass = 5.00d+00
>
>
>
> ! winds
>
>       cool_wind_full_on_T = 9.99d9
>
>       hot_wind_full_on_T = 1d10
>
>       cool_wind_RGB_scheme = 'Reimers'
>
>       cool_wind_AGB_scheme = 'Blocker'
>
>       RGB_to_AGB_wind_switch = 1d-4
>
>       Reimers_scaling_factor = 0.1d0
>
>       Blocker_scaling_factor = 0.2d0
>
>
>
> ! when to stop
>
>       max_age = 1.4d+10
>
>
> ! mesh
>
>       max_allowed_nz = 20000 ! preventing 'fail to adjust mesh'
>
>
> ! diffusion
>
>       diffusion_v_max = 1d-5
>
>       do_element_diffusion = .false. ! Activated during MS and settling
> only
>
>       diffusion_use_full_net = .true.
>
>
>       ! Solver controls.
>
>       diffusion_use_cgs_solver = .true.
>
>       diffusion_use_isolve = .true.
>
>       diffusion_rtol_for_isolve = 1d-4
>
>       diffusion_atol_for_isolve = 1d-5
>
>       diffusion_min_X_hard_limit = -1d-4
>
>       diffusion_maxsteps_for_isolve = 1000
>
>       diffusion_isolve_solver = 'ros2_solver'
>
>
>       ! Timestep controls to prevent steps that are difficult for
> diffusion.
>
>       diffusion_steps_limit = 20
>
>       diffusion_steps_hard_limit = 150
>
>       diffusion_iters_limit = 50
>
>       diffusion_iters_hard_limit = 100
>
>
> ! convection
>
>       mlt_option = 'Henyey'
>
>       mixing_length_alpha = 1.8
>
>       use_Ledoux_criterion = .true.
>
>       thermohaline_coeff = 1d2
>
>
> ! timesteps
>
>       max_years_for_timestep = 1d7
>
>       time_delta_coeff = 1.0d0
>
>       use_gold2_tolerances = .true.
>
>
>
>       varcontrol_target = 1d-3
>
>       dX_nuc_drop_limit = 1d-2
>
>       delta_HR_limit = 0.1
>
>
> ! Atm
>
>       atm_option = 'T_tau'
>
>       atm_T_tau_relation = 'Eddington'
>
>       atm_T_tau_opacity = 'iterated'
>
>       atm_T_tau_max_iters = 100
>
>
> ! Overshoot
>
>       overshoot_scheme(1) = 'exponential'
>
>       overshoot_zone_type(1) = 'any'
>
>       overshoot_zone_loc(1) = 'any'
>
>       overshoot_bdy_loc(1) = 'any'
>
>       overshoot_f(1) = 0.014
>
>       overshoot_f0(1) = 0.004
>
>
> / ! end of controls namelist
>
>
> *AGB inlist*
>
>
> &star_job
>
>
>       load_saved_model = .true.
>
>       load_model_filename = 'end_he_core_burn.mod'
>
>       save_model_when_terminate = .true.
>
>       save_model_filename = 'end_AGB.mod'
>
>       required_termination_code_string = ''
>
>       set_initial_model_number = .true.
>
>       initial_model_number = 0
>
>       pgstar_flag = .true.
>
>
>
> / ! end of star_job namelist
>
>
> &eos
>
>
> / ! end of eos namelist
>
>
> &kap
>
>
> / ! end of kap namelist
>
>
> &controls
>
>       log_directory = 'LOGS_AGB'
>
>
> ! energy conservation
>
>       energy_eqn_option = 'eps_grav'
>
>       use_time_centered_eps_grav = .true.
>
>
> ! when to stop
>
>       envelope_mass_limit = 1e-3
>
>
> ! winds
>
>       max_wind = 1d-2
>
>
> ! solver
>
>       retry_hold = 0
>
>       neg_mass_fraction_hold = 3
>
>
>
> ! mesh
>
>       mesh_delta_coeff = 1.0
>
>
> ! convection & dragging
>
>       num_cells_for_smooth_gradL_composition_term = 0
>
>       drag_coefficient = 1d0
>
>       min_q_for_drag = 0.98d0
>
>       alpha_semiconvection = 0.1
>
>
> ! timesteps
>
>
>       varcontrol_target = 1d-4
>
>       delta_lgL_He_limit = 0.005
>
>       lgL_He_burn_min = 2.0
>
>       dX_nuc_drop_limit = 1d-3
>
>       delta_HR_limit = 0.1
>
>       max_timestep_factor = 1.01d0
>
>
>
>       delta_lgTeff_limit = 0.005
>
>       delta_lgTeff_hard_limit = 0.01
>
>       delta_lgL_limit = 0.02
>
>       delta_lgL_hard_limit = 0.05
>
>
>
> ! output
>
>       num_trace_history_values = 2
>
>       trace_history_value_name(1) = 'rel_E_err'
>
>       trace_history_value_name(2) = 'log_rel_run_E_err'
>
>
>       !report_solver_progress = .true. ! set true to see info about solver
> iterations
>
>       !report_ierr = .true. ! if true, produce terminal output when have
> some internal error
>
>       !stop_for_bad_nums = .true.
>
>
> / ! end of controls namelist
>
>
>
> _______________________________________________
> 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/20260326/06e0bd34/attachment.htm>