[Mesa-users] Resolving TP-AGB

Thu Mar 26 14:50:57 UTC 2026

Hi Théo,

There is a rich literature on the difficulty and ultimatel inability to evolve TP-AGB stars realistically with 1D stellar evolution codes, going back half a century. Alan Sweigart was on of the early contributors to this field.

With a 1D stellar evolution code we need to keep in mind that it is genuinely blind to 3D macro physics. The 1D version of the conservation laws can’t see 3D macro physics. TP-AGB stars with their fluffy envelopes, with mass out into the shallow wings o the gravitational potential, and with their global convection morphologies and pulsations are dominated by 3D macro physics.

A 1D code not converging may not be due to a bug or due to a weak numerical scheme or even due to a poor choice of insist parameters. A 1D code not converging may just be saying to you: “Don’t ask me to do the impossible”

Best, Falk.

—
Dr Falk Herwig, Professor

Dept of Physics & Astronomy
University of Victoria
https://web.uvic.ca/~fherwig

On Mar 26, 2026, at 04:51, Théo-Courty via Mesa-users <mesa-users at lists.mesastar.org> wrote:

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/6c1a354a/attachment.htm>