[mesa-users] bad angular_momentum_j
Kenny Van
kvan at ualberta.ca
Tue Aug 9 12:37:48 EDT 2016
Hi Pablo,
I figured out what the issue was. I didn't realize that the different
timestep values were in different units so that caused issues.
Thanks for the help
On 8 August 2016 at 17:06, Kenny Van <kvan at ualberta.ca> wrote:
> Hi Pablo,
>
> I'll verify the problem I'm running into and get back to you.
>
> As for the snippet I added before, that was just an attempt to see if I
> could get MESA to decrease the timestep in the very first iteration. I
> include a similar code in extra_binary_check_model and just use a small
> initial timestep now.
>
> On 8 August 2016 at 17:02, Pablo Marchant <pamarca at gmail.com> wrote:
>
>> Kenny, if you want to debug this just do some small edits to
>> binary/private/binary_evolve.f90, in particular around these lines
>>
>> 257 ! get jdot and update orbital J
>> 258 b% jdot = get_jdot(b, b% mtransfer_rate, b% xfer_fraction)
>> 259 b% angular_momentum_j = b% angular_momentum_j + b% jdot*b%
>> s_donor% time_step*secyer
>> 260
>> 261 if (b% angular_momentum_j <= 0) then
>> 262 stop 'bad angular_momentum_j'
>> 263 end if
>>
>> Since the stop condition literally has a check for negative angular
>> momentum, I have a hard time believing what you say. You can put write
>> statements in there to verify you are getting the numbers you expect. If
>> this still gives you no answers, then if you want you can send a minimal
>> example that fails (truly minimal please), and I can give it a look.
>>
>> PS:
>>
>> Just saw in more detail this snippet you posted before
>>
>> write(*,*) "starting modified MB"
>> if (b% angular_momentum_j <= 0) then
>> write(*,*), "angular_momentum_j <= 0", b% angular_momentum_j
>> b% s1% dt_next = min(b% s1% dt * 0.50, b% s1% dt_next)
>> extras_binary_startup = retry
>> end if
>>
>> You said you add this to extras_binary_startup. That will do nothing,
>> extras_binary_startup is only ran at the beginning of the simulation, and
>> orbital angular momentum will be positive. You should be doing this in
>> extras_binary_check_model.
>>
>> On Mon, Aug 8, 2016 at 9:12 PM, Kenny Van <kvan at ualberta.ca> wrote:
>>
>>> Hi,
>>>
>>> So unfortunately after some testing I've started running into some
>>> problems that I cant quite solve. I've included a scaling factor in the
>>> timestep check to try and retry the run if we encounter a problem with
>>> angular momentum.
>>> The problem that I'm running into now is that the code triggers the "bad
>>> angular_momentum_j" portion even when the angular_momentum_j > 0. At every
>>> step I have the extras_binary_check_model output the
>>>
>>> b% angular_momentum_j,
>>>
>>> as well as
>>>
>>> b% angular_momentum_j + b% jdot * b% s1% dt_next
>>>
>>> so I have an idea of what the current angular momentum is as well as
>>> approx what it will be in the next timestep.
>>>
>>> At the point of the code breaking both of these values are positive and
>>> I encounter the termination code
>>>
>>> STOP bad angular_momentum_j.
>>>
>>> I'm not quite sure what I can do from here.
>>>
>>> On 4 August 2016 at 21:40, Pablo Marchant <pamarca at gmail.com> wrote:
>>>
>>>> You should be able to do that in extras_binary_startup. Anyhow, I'd
>>>> just go for the simple approach of setting a small initial timestep, if its
>>>> too small it will grow quite fast.
>>>>
>>>> On Fri, Aug 5, 2016 at 3:36 AM, Kenny Van <kvan at ualberta.ca> wrote:
>>>>
>>>>> Yes, but I was also hoping to create an adaptive way for MESA to
>>>>> adjust the first step if necessary
>>>>>
>>>>> On 4 August 2016 at 17:00, Pablo Marchant <pamarca at gmail.com> wrote:
>>>>>
>>>>>> Are you running into these issues at the first timestep? If so, you
>>>>>> need to ensure the first step is small as well.
>>>>>>
>>>>>> On Thu, Aug 4, 2016 at 11:02 PM, Kenny Van <kvan at ualberta.ca> wrote:
>>>>>>
>>>>>>> Hi Everyone,
>>>>>>>
>>>>>>> So I've checked the units for Jdot and they seem fine, Jorb/Jdot is
>>>>>>> on the order of 1e12. The issue I'm currently running into is that the code
>>>>>>> doesn't quite seem to be doing what I think it should be doing.
>>>>>>> For extra checks I also included:
>>>>>>>
>>>>>>> write(*,*) "starting modified MB"
>>>>>>> if (b% angular_momentum_j <= 0) then
>>>>>>> write(*,*), "angular_momentum_j <= 0", b%
>>>>>>> angular_momentum_j
>>>>>>> b% s1% dt_next = min(b% s1% dt * 0.50, b% s1% dt_next)
>>>>>>> extras_binary_startup = retry
>>>>>>> end if
>>>>>>>
>>>>>>> to extra_binary_startup. Just to check if my initial timestep is
>>>>>>> sufficiently small and the MESA code doesn't seem to be triggering this. I
>>>>>>> have a similar portion of code in extras_binary_check_model,
>>>>>>> but I haven't reached a point where the simulation breaks to
>>>>>>> determine if this works or not. Does the file binary_evolve.f90 trigger and
>>>>>>> instantly kill the MESA simulation before the code checks the
>>>>>>> run_binary_extras.f file?
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On 3 August 2016 at 21:25, Kenny Van <kvan at ualberta.ca> wrote:
>>>>>>>
>>>>>>>> Hey Pablo,
>>>>>>>>
>>>>>>>> Thanks for the advice, I'll try out the time step adjustment to see
>>>>>>>> how it goes. I'll also double check the units to be 100% sure that I
>>>>>>>> haven't crated any jdots that are a few orders too big.
>>>>>>>>
>>>>>>>>
>>>>>>>> On 3 August 2016 at 20:53, Pablo Marchant <pamarca at gmail.com>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> Also, be sure you are getting your units right. That's an easy way
>>>>>>>>> to create absurdly large jdots (been there). For the process you are
>>>>>>>>> considering, what's the decay timescale for Jorb (ie. Jorb/Jdot)? If that
>>>>>>>>> is absurdly tiny then I'd worry.
>>>>>>>>>
>>>>>>>>> On Thu, Aug 4, 2016 at 4:50 AM, Pablo Marchant <pamarca at gmail.com>
>>>>>>>>> wrote:
>>>>>>>>>
>>>>>>>>>> Kenny, I'd recommend you to manually adjust the timestep, this is
>>>>>>>>>> something we discussed here in mesa-users a while ago, copy from Rob
>>>>>>>>>> Farmer's answer to create a hard timestep limit
>>>>>>>>>>
>>>>>>>>>> (Untested)
>>>>>>>>>>
>>>>>>>>>> In extras_startup:
>>>>>>>>>> s%xtra1=0.d0
>>>>>>>>>>
>>>>>>>>>> extras_finish_step:
>>>>>>>>>> s%xtra1=s%log_R
>>>>>>>>>>
>>>>>>>>>> extras_check_model:
>>>>>>>>>>
>>>>>>>>>> if( abs( 10**s%xtra1_old - 10**s%log_R) > EPS) then
>>>>>>>>>> extras_check_model = retry
>>>>>>>>>> s% dt = s%dt * SOME_SCALE_FACTOR
>>>>>>>>>> end if
>>>>>>>>>>
>>>>>>>>>> And to create a soft timestep limit:
>>>>>>>>>>
>>>>>>>>>> (also untested)
>>>>>>>>>>
>>>>>>>>>> if( abs( 10**s%xtra1_old - 10**s%log_R) > EPS) then
>>>>>>>>>> s% dt_next = min(s% dt_next, s%dt * SOME_SCALE_FACTOR)
>>>>>>>>>> end if
>>>>>>>>>>
>>>>>>>>>> If the first timestep is the issue, then set a smaller initial
>>>>>>>>>> timestep. Check this from defaults/star_job.defaults
>>>>>>>>>>
>>>>>>>>>> 435 !### set_initial_dt
>>>>>>>>>> 436 !### years_for_initial_dt
>>>>>>>>>> 437 !### seconds_for_initial_dt
>>>>>>>>>> 438 ! if true, set initial timestep, dt, in years
>>>>>>>>>> 439
>>>>>>>>>> 440 set_initial_dt = .false.
>>>>>>>>>> 441 years_for_initial_dt = 0
>>>>>>>>>> 442 seconds_for_initial_dt = 0
>>>>>>>>>>
>>>>>>>>>> On Wed, Aug 3, 2016 at 4:56 PM, Kenny Van <kvan at ualberta.ca>
>>>>>>>>>> wrote:
>>>>>>>>>>
>>>>>>>>>>> Hi,
>>>>>>>>>>>
>>>>>>>>>>> I'm currently working on adapting additional magnetic braking
>>>>>>>>>>> terms into MESA using the run_binary_extras.f file. I'm currently running
>>>>>>>>>>> into an issue where the amount of angular momentum being removed is too
>>>>>>>>>>> great in a single timestep and causing the simulation to break. Looking at
>>>>>>>>>>> the code it seems like the MESA binary evolution dies immediately if the
>>>>>>>>>>> angular momentum loss is too great instead of retrying with a smaller
>>>>>>>>>>> timestep. Is there a way to get MESA to retry with a smaller timestep when
>>>>>>>>>>> it encounters this issue?
>>>>>>>>>>>
>>>>>>>>>>> Thanks
>>>>>>>>>>>
>>>>>>>>>>> ------------------------------------------------------------
>>>>>>>>>>> ------------------
>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> mesa-users mailing list
>>>>>>>>>>> mesa-users at lists.sourceforge.net
>>>>>>>>>>> https://lists.sourceforge.net/lists/listinfo/mesa-users
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> --
>>>>>>>>>> Pablo Marchant Campos
>>>>>>>>>> M.Sc on Astrophysics, Universidad Católica de Chile
>>>>>>>>>> PhD student, Argelander-Institut für Astronomie
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> --
>>>>>>>>> Pablo Marchant Campos
>>>>>>>>> M.Sc on Astrophysics, Universidad Católica de Chile
>>>>>>>>> PhD student, Argelander-Institut für Astronomie
>>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> Pablo Marchant Campos
>>>>>> M.Sc on Astrophysics, Universidad Católica de Chile
>>>>>> PhD student, Argelander-Institut für Astronomie
>>>>>>
>>>>>
>>>>>
>>>>
>>>>
>>>> --
>>>> Pablo Marchant Campos
>>>> M.Sc on Astrophysics, Universidad Católica de Chile
>>>> PhD student, Argelander-Institut für Astronomie
>>>>
>>>
>>>
>>
>>
>> --
>> Pablo Marchant Campos
>> M.Sc on Astrophysics, Universidad Católica de Chile
>> PhD student, Argelander-Institut für Astronomie
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20160809/31e282a9/attachment.html>
More information about the Mesa-users
mailing list