[mesa-users] optical_depth -- and a message for MESA BLACK BELT USERS
jingluan at caltech.edu
jingluan at caltech.edu
Fri Mar 23 18:44:25 EDT 2012
> Hi,
>
> Let's see if there are any black belts out there who can help on this.
>
> But first you need to provide enough info to make it possible for them.
> Make a saved model a few steps before the bad one.
Thank you, Bill :-)
I start with WD0.6-Hrich-purecooling.mod, which is attached. I evolved a 7
Msun star from pre-main sequence. When the helium core is 0.5Msun, I
started losing mass on a thermal timescale until the total mass is
0.6Msun. Then I terminated the mass-loss and evolves the 0.6Msun star
until no nuclear burning is going on, i.e., it is purely cooling down
white dwarf which is saved as this model I attached here.
> Make a modified inlist that starts by loading that newly saved model.
inlist_project is also attached. I change the net to 17_to_fe56.net which
I also used to make the starting model from pre-main sequence star.
> Confirm that the problem still happens when you run with the new inlist.
> Send the starting model and the modified inlist.
> If possible, provide a plot to show where things are bad.
> explain what is wrong in enough detail so someone else can see the
> problem.
What I check is quite simple. I first find the shell number kphot for the
photosphere. I expect s% tau(kphot)~0.66. But s% tau(kphot)~360.
then I thought maybe s% tau(:) is not the optical depth. Then I add some
lines in src/run_star_extras.f(also attached) to calculate the optical
depth (tau_phot):
kphot=0
tau_phot=0.
do k=1, s% nz
if (s% R(k)/rsol<= s% photosphere_r) then
kphot=k
tau_phot=tau_phot+s% opacity(k)*s% rho(k)*(s% R(k)-s% R(k+1))
exit
end if
end do
And it turns out that tau_phot is an even crazily larger number than s%
tau(kphot).
Then I thought maybe mesa had not yet converged but it did not report any
err or warning information which made me worried that "Do I need to check
manually every variable in s even if mesa seems to run happily?"
Thank you very much!
Sincerely,
Jing
> Include the log file for that model.
> include version number info for mesa, the compiler, and the operating
> system.
>
> -Bill
>
>
>
>
> On Mar 23, 2012, at 10:18 AM, jingluan at caltech.edu wrote:
>
>> Dear Bill,
>>
>> Thank you so much for guiding me. Although I have not yet figured out
>> what
>> specific problem causes the strange tau(:), I will try to follow yours
>> to
>> see.
>>
>> I still have one question. If this weird tau were caused by
>> non-convergence, I would NOT be able to find it out if I did not
>> accidentally check the tau values. MESA does not report err information
>> like 'convergence problem', normally I would accept this model as a
>> converged one and trust the outputs by mesa, e.g. the figures plotted by
>> pgplot. But actually there is something wrong in the model (e.g. the tau
>> in this case) which is NOT revealed. I feel that in practice people do
>> not
>> check every variable's value in s to see whether it makes sense if MESA
>> does not report error.
>>
>> Thank you so much!
>>
>> Sincerely,
>> Jing
>>
>>
>>>
>>>
>>>
>>>
>>>>>> 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 --------------
A non-text attachment was scrubbed...
Name: WD0.6-Hrich-purecooling.mod
Type: audio/x-mod
Size: 1111578 bytes
Desc: not available
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20120323/6740a2c2/attachment.bin>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: inlist_project
Type: application/octet-stream
Size: 12527 bytes
Desc: not available
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20120323/6740a2c2/attachment.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: run_star_extras.f
Type: text/x-fortran
Size: 150620 bytes
Desc: not available
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20120323/6740a2c2/attachment-0001.bin>
More information about the Mesa-users
mailing list