MESA equations, time steps and meshes

1. MESA equations, time steps and meshes Michael Palmer

2. Overview of steps • Takes input files and also initializes the modules • These are then stored in memory • Prepares to take new time step but first checks to see if remeshing needs to be performed (more later) • Adjusts model to reflect mass loss or gain, abundance changes, calculates convective diffusion coefficients. Solves for this structure using Newton-Raphson solver • Next time step estimated (more later) • Output files generated (more later)

3. In the Beginning • Can used predefined models for solar metallicity in a variety of mass ranges and evolutionary stages • ZAMS - 32 masses 0.08 - 100Msun • Low mass PMS - 0.001 - 0.025Msun • WD - He cores 0.15 - 0.45Msun, CO cores 0.496 - 1.025Msun and ONe cores 1.259 - 1.376 Msun • Create PMS, specify mass, uniform composition, luminosity and temperature. Uses modules and NR-solver to find c for that model. • Works best for a mass range of 0.02 - 50Msun

4. Structure • 1-D spherically symmetric • Divided into cells, number depends on complexity • Numbered from outside in • Each cell has variables that are mass averaged or defined at the outer boundary • Tk, Xi,k, and k are mass averaged variables • Boundary variables mass interior to mk,rk,Lk, and vk • Inner boundary of innermost cell usually defined as centre of star, Rc,Lc and vc zero unless defined otherwise (special circumstances) • How does it handle equations of stellar evolution? Figure 9Paxton et al. 2010

5. Equations of Stellar Evolution • Ones we are used to • An example from the MESA code Paxton et al. 2010 (MESA Paper)

6. ESE continued • For enhanced stability pressure equation is divided by P_k = (Pk-1dmk + Pkdmk-1)/(dmk + dmk-1) • Boundary conditions constructed using Ps and Ts. The difference from surface to centre of first cell is found by • This gives us our first boundary conditions for both P1 and T1 Paxton et al. 2010 (MESA Paper)

7. Mesh Adjustment • Cells may be merged or split at the beginning of each time step, however, remeshing designed so most cells are left alone during remeshing • Two stages for remeshing, planning and adjustment stage

8. Planning • Decides which cells may be split or merged • Minimizes splits but maximizes mergers • Must ensure magnitude between two adjacent cells are below a certain threshold, ex ∆log(P) < P where the default for P is 1/30 • Local reductions in magnitudes will result in higher resolution in the desired region of the star. Example when ∆log(P) is small compared to ∆log(nuc) • Also increase near boundaries of convective zone, over a distance of the pressure scale height • Where ∆log(Xi) is large compared to ∆log(P) • Near locations where spatial gradients exist in the most abundant species

9. Adjustment • Executes remesh plan • Cells needing to be split have an interpolation in mass performed to obtain luminosity and enclosed volume at new cell boundary. From this a new cells density found by • New composition mass fractions calculated. For cell mergers simple, for cell divisions neighbouring cells used to linearly approximate mass fraction as a function of mass coordinate. Slopes adjusted so mass fractions sum to one everywhere and functions are then integrated over new cell mass to find mass fraction

10. Adjustment continued • Temperature calculations depend on degeneracy of electrons. For degenerate electrons, split cells take on their former temperature (from before they were split). While merging cells take a mass average of their temperatures • For non degenerate electrons, a parabola is formed for the specific internal energy of neighbouring and original cell. The parabola is integrated over the new cell to find total internal energy. EOS module is called with new parameters (composition and density) until a temperature that predicts this internal energy is found

11. Time Step • Time step based on digital control theory • Uses both previous and current time step to calculate the next • c unweighted average of the relative changes in log() log(T) and log(R) across all cells, also uses previous c • Target value of t = 10-4 • f(x) = 1 + 2tan-1[0.5(x-1)] • Next time step calculated by • Where t-1, t and t+1 are the previous, current and next time step respectively • Proposed step can be reduced based on tests • Can set both hard and soft limits, if change exceeds hard limit, solution is rejected and retries. If change exceeds soft limit next time step is reduced proportionally

12. Mass Adjustments • Mass adjustments, for either loss or gain, are done at the beginning of a new time step, before structure is solved • Several ways to account for dM/dt • Section 6.6 - To many to list • Write your own module for a mass loss or gain, example given is for mass loss of carbon stars • Total mass is M = Mc + Mm • Mc is core mass and Mm is modeled mass • Relative cell mass, dqk = dmk/Mm, and relative mass interior to cell qk = mk/Mk or 1-∑dqi ,i = 1,..,k-1, must be defined as this for accuracy close to outer envelope

13. Mass part 2 • Once change, M found, mass structure changed. This does not effect cell number. • Changes mass location of cells and revises composition based on location • Does not effect initial solutions, T,,.. • Mass structure divided into 3 regions • Inner, intermediate and outer • Inner region and outer boundaries defined by temperature, logT = 6 (inner) logT = 5 (outer) • Region expanded for stability if intermediate mass is not much larger than M. First expanded by moving outer boundary to surface and if mass in intermediate is still to small can move inner boundary in (to a point). Intermediate mass is not allowed to change by more than 10% or by more than a factor of 2

14. Mass finale • When regions defined, dqk rescaled by M / (M + M) in inner region • Outer region dqk left alone • Intermediate scaled so ∑dqk = 1 • Composition for last two are than updated

15. Convergence • yi is trial solution, F(yi) is residual, yi is the correction term and [dF/dy]i is the Jacobian • Jacobain gets its derivatives from modules like the eos module • Converges on a solution by building upon yi • Solution accepted when solution meets convergence criteria. However if convergence cannot be achieved • Retry => Uses a smaller time step to see if this results in convergence • Backup => If retry fails. Returns to previous model with smaller time step than was used to find it. • Keeps doing this until convergence happens or time step goes below specified threshold, sequence terminated then

16. Output • Two types • Photo and save • Photo - Binary file that has complete current state of star. There will be no changes in evolution if you restart from this • Save - Not as complete as binary • Two types - log and profiles • Log - contain properties like stellar age • Profiles - contains properties at specified time steps at each zone from centre to surface • Other options to, can produce images using PGPLOT or output into another format that a different code uses

17. Finally The End? • Please see resolution sensitivity (Not fair to Athira to take her time!) • Questions?