[mesa-users] mesa version 4566

[Bill Paxton]

You may have noticed that the entire month of September went by with no new mesa releases.
I'm sure you are getting anxious by now.  So here you go -- a bunch of new toys in mesa version 4566.


opacities 

	ideally, we'd just compute the opacity from the temperature, density, and the full set of isotope abundances.
	we don't do that yet because it is still too computationally expensive.  hopefully that will change before long.
	so for now we are stuck with collapsing the full abundance information down to a small set of parameters.
	We have two types of opacity tables for two choices of abundance parameters (this is the OPAL standard).
	Type1 opacities are for a fixed "normal" mixture of metals so tables are given for different (X,Z) -- 2 parameters.
	Type2 opacities allow for the normal mixture to be enhanced by extra C and O.
		the tables are given for (X,Zbase,dXC,dXO) -- 4 parameters. note that Z = Zbase + dXC + dXO
	Of course there are complications.  The Type2 tables have a less dense and less extensive set of X and Zbase
		values than the Type1 tables have for X and Z.  So there are combinations of (X,Zbase,dXC,dXO)
		that need to be redirected from Type2 to Type1.   We are now doing that automatically using several new controls.
	
	here are the new controls along with default values
         ! switch to Type1 if X too large
         kap_Type2_full_off_X = 0.70001d0 ! Type2 full off for X >= this
         kap_Type2_full_on_X = 0.68d0 ! Type2 full on for X <= this
         
         ! switch to Type1 if Zbase too large
         kap_Type2_full_off_Zbase = 0.100d0 ! Type2 full off for Zbase >= this
         kap_Type2_full_on_Zbase = 0.095d0 ! Type2 full on for Zbase <= this
         
         ! switch to Type1 if dXC+dXO too small
         kap_Type2_full_off_sum_dXCO = 0.0d0 ! Type2 full off for dXC+dXO <= this
         kap_Type2_full_on_sum_dXCO = 0.01d0 ! Type2 full on for dXC+dXO >= this

	Using these limits and the given values of X, Zbase, and dXC+dXO, we determine whether
	the opacity will be calculated as pure Type2, pure Type1, or some blend of the two for transition values.

	We calculate Zbase, dXC, and dXO automatically from (Z,XC,XO,base_fC,base_fO)
		where XC = mass fraction of carbon, XO = mass fraction of oxygen,
			base_fC = fraction XC/Z for unenhanced material,
			base_fO = fraction XO/Z for unenhanced material.
		the equations relating these are
			Z =  Zbase + dXC + dXO
			XC = dXC + Zbase*base_fC
			XO = dXO + Zbase*base_fO
		we switch to Type1 if we get a negative value for any of Zbase, dXC, or dXO.

		the user can set base_fC and base_fO.   here are the default values.
         		base_fC = 0.173312d0 ! base value for C fraction of metals
         		base_fO = 0.482398d0 ! base value for O fraction of metals

	by default, mesa/star still just uses Type1 opacities.  That's fine for cases like solar models that
	only do H burn and don't have any regions of enhanced C or O.  But if you are getting into He burning
	you will have enhanced C and O, so you should probably be using Type2 (with the usual caveat
	that things may break when you make even a reasonable change like going from Type1 to Type2).

	to enable the new Type2 scheme for mesa/star, set use_Type2_opacities = .true. in your &controls inlist.
	(NOTE: this flag was previously called use_CO_enhanced_opacities; I've given it a new name to
		make sure people who might care are forced to think about what has changed.)


newton damping to reduce problems with negative mass fractions

	each iteration of the newton solver for star uses a linear approximation to create a vector of corrections to the model.
	the linear solution knows nothing about the requirement that mass fractions stay non-negative, and we sometimes
	get proposed corrections that if applied in full would give such results.  previous versions of star rejected proposed
	solutions that gave negative abundances exceeding a limit (min_xa_hard_limit).  for smaller violations of non-negativity,
        star used the "clipping" scheme in which all mass fractions are clipped to the range 0 to 1 and then the sum is
	renormalized to 1.  The obvious problem with the clipping scheme is that it violates conservation of species.

	to reduce the occurrence of such non-conservation, star now uses the "newton damping" scheme.  This checks for
	proposed corrections that would produce negative abundances and multiplies the entire correction vector by a 
	factor less than 1 to fix the problem.   in many cases, this is enough to get the newton iterations back on track.
	however, in other cases, it may force so many small corrections that the newton solver gives up forcing a backup.
	that is considered to be a reasonable price to pay to get improved conservation of species and more accurate solutions.


element diffusion

	there are several new controls to let you clump the outer cells into a single cell for the purpose of element diffusion.
	this helps with both performance and accuracy.   here are the new items in star/defaults/controls_defaults
            
         diffusion_min_dq_at_surface = 1d-9
            ! treat at least this much at surface as a single cell for purposes of diffusion
            
         diffusion_min_T_at_surface = 1d4 ! this should be large enough to ensure hydrogen ionization
            ! treat cells cells at surface with T < this as a single cell for purposes of diffusion
            
         diffusion_min_dq_ratio_at_surface = 50
            ! combine cells at surface until have total mass >= this factor times the next cell below them
            ! this helps with surface boundary condition for diffusion by putting large cell at surface
	
	as you might recall, element diffusion is solved separately from the rest of the stellar equations.
	in the past, star has done diffusion first and then the rest of the structure and abundances.
	we now allow the user to partition the diffusion into a part done before and a part done after the other equations.
         
         diffusion_frac_before = 1d0 
            ! do diffusion for this fraction of the timestep before doing solve for structure+compositon
            ! do rest after.


create initial model

	thanks to Phil Arras we have a new alternative to create pre-main-sequence that works vastly better for very low masses.
	Phil has plans to make it even better, but you might want to start using it now for making planets.
	here are the new controls in &star_job
   	
       create_initial_model = .false. ! this is an alternative to create_pre_main_sequence_model
      ! creates an adiabatic, contracting model for given mass and radius
      ! assumes no nuclear burning and constant entropy. ignores radiation pressure.
      
      radius_in_cm_for_create_initial_model = 0 ! in cm
      mass_in_gm_for_create_initial_model = 0 ! in grams
      
      initial_model_relax_num_steps = 10 
         ! let initial model settle in for this many steps before changing anything else


we've recently done a lot of renaming to make things more consistent.  I missed "terminal_cnt".  It is now called
"terminal_interval" to match profile_interval and history_interval.


new wind option: 'van Loon' -- see van Loon et al. 2005, A&A, 438, 273
option to use different schemes for low T in "Dutch" wind.
	Dutch_wind_lowT_scheme = 'de Jager' ! alternative is 'van Loon'


I've removed the "supersonic_wind" option because I realized the implementation was bogus.  


I've removed approx22 and approx24 nets because of problems with compound reactions for cno falsely activated at high T.  Please use approx21 instead.


we have a new way to calculate Brunt N^2 thanks to Mike Montgomery.   enable it by setting use_brunt_basic2 = .true.  this is similar to the Brassard scheme,
but generalized to work for arbitrary abundance gradients rather than just Helium.  more details will come in the future.


many improvements to star/astero_adipls since last release.  thanks to Mia Lundkvist for being the brave pioneer user!  


there are several small but important fixes the mesa version of ADIPLS.  many thanks to Jørgen Christensen-Dalsgaard.


the bcyclic matrix solver now supports sparse blocks using the KLU sparse solver.  
this comes into play when using very large nets with mesa/star (> 100 isotopes).


we have a new "relax" option -- relax_to_xaccrete -- that gradually transforms the star to a given uniform composition


new "when to stop" controls based on total mass of Helium
         star_He_mass_min_limit = 0
            ! stop when star helium mass in Msun units is < this.  <= 0 means no limit.
         star_He_mass_max_limit = 0
            ! stop when star helium mass in Msun units is > this.  <= 0 means no limit.


fixed bugs in the calculation of the turbulent pressure term in the dlnP/dm equation. (see conv_dP_term_factor in controls.defaults)


new "hook" for alternative implementation of non-nuclear neutrino losses.  see star/others/other_neu.f


newly functional routine for spin-orbit synchronization in binaries -- see star/binary_systems/tide_routines.inc


new members of test_suite

	make_massive_with_uniform_composition -- uses relax_to_xaccrete

	make_low_mid_mass_with_uniform_composition -- uses relax_to_xaccrete

	make_0.6m_co_wd -- a bit different recipe than the one in make_co_wd which makes a more massive wd

	wd_diffusion -- does wd cooling with element diffusion enabled

	synch_spin_to_orbit -- based on rlo, donor spin is adjusted using other_torque to synchronize with binary orbit.
		this uses the synch_spin_to_orbit routine in star/binary_systems/tide_routines to calculate jdot vectors at each step

	axion_cooling -- demonstrates use of other_neu hook. based on work using MESA by Friedland, Giannotti, and Wise 


there are probably some other things I'm missing here.  let me know if you find something that seems important.
and as usual, you can expect a flurry of updates with small corrections to follow this one.
so check the webpage to find the current version number if you wait a few days before checking out the new version.


And finally, big thank you's to the folks like Josiah Schwab, Max Katz, and others who are stepping forward to respond to mesa-users questions.
I really appreciate it.  

Cheers,
Bill