[mesa-users] optical_depth -- and a message for MESA BLACK BELT USERS

[Bill Paxton]

>>>   MESA USERS

Please consider that if you are publishing papers using mesa,
it might be a good idea to be able to poke around inside the code enough
to get an idea of what is going on.  If you can do that, you are a mesa black belt!
Are you a black belt yet?  Do you want to become one?
If so, take a moment to see if you could have answered this email yourself.   





On Mar 23, 2012, at 9:00 AM, jingluan at caltech.edu wrote:
> 
> At the kth shell where it satisfies s% R(k)/rsol==s% photosphere_r, I
> found that s% tau(k)=361 in a model by mesa2664;
> 
> Then I calculate a variable called tau_phot by summing s% opacity(j)*s%
> rho(j)*(s% R(j)-s% R(j+1)) over j=1,..,k; and I got tau_phot=5184;
> 
> It is pretty weird. What could be the problem? The model has not yet been
> converged well? But mesa does not report any ierr information though.


Hi Jing,

Good question!   Rather than just giving an answer, let's go through
the process of discovering how things work (which is what I just did!).

Question: why do I find weird value for tau before the model has converged.


1st step: grep to find where tau is set.
	cd star/private
	grep tau *.f -n > lst

	open the file lst and don't panic when you see a lot of lines
	we don't care about places that are reading tau;
	we want to find where it is being set.
	this looks like the best bet:

star_utils.f:566:            tau(k) = tau(k-1) + dtau

 	that is in a routine called "get_tau" which is a sign that we're on the right track.
	open star_utils and take a look at get_tau -- bingo! 


2nd step: find where and when get_tau is being called.

grep 'call get_tau' * -n	
report.f:203:         call get_tau(s, s% tau)

	open report and find the call on get_tau.  it is in a routine called do_report


3rd step: find where and when do_report is being called.

grep 'call do_report' * -n
evolve.f:285:         call do_report(s, ierr)
read_model.f:67:         call do_report(s, ierr) 
	
	we find 2 places: one in evolve and the other in read_model.
	since we're doing evolution rather than reading, we want the one in evolve.

	open evolve.f and find the call on do_report.
	turns out that it is called from do_evolve_step,
		after the model has converged (in do_struct_burn_mix) .
	

understanding the order of events in do_evolve_step is the hard part of this.
there's a lot of extra things in there for error checks and timekeeping, but
fundamentally here's the order in which things are done:

	prepare_for_new_step -- does remesh
	prepare_for_new_try -- save stuff in case need to do retry or backup later
	do_winds -- calculate mdot for this step
	do_adjust_mass -- add or remove material according to mdot
	do_set_vars -- set variables for new mesh and mass
	do_element_diffusion -- optional
	save_pre_hydro_values -- save starting values for solver
	do_struct_burn_mix -- solve for new model
	smooth_newly_non_conv -- optionally smooth abundances
	do_brunt_N2 -- optionally calculate brunt N^2 at each cell
	do_report -- calculate information for logs and profiles

One more subtle point that "black belt" mesa user's should know:
the diffusion coefficients for convective mixing, overshooting, semiconvection,
thermohaline mixing are calculated in the call on do_set_vars from
do_evolve_step and remain fixed during the rest of the calls from do_evolve_step;
in particular, they are not updated between iterations of the solver, so they
are not being treated in a fully implicit manner.  this is done for
numerical reasons only.

Okay --- now I expect that there will be lots of mesa users out there
in the future who'll be able to dig out the answers to questions like this.
Consider it part of your qualification for getting a mesa black belt.  ;-)

Cheers,
Bill




-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20120323/44621b7f/attachment.html>