[mesa-users] approx13

Ken Shen kenshen at astro.berkeley.edu
Fri May 16 18:19:16 EDT 2014


Looks good, thanks!


On Fri, May 16, 2014 at 2:32 PM, Bill Paxton <paxton at kitp.ucsb.edu> wrote:

> Good job!   I think that's it.  That code wasn't ever intended for
> anything but approx20/21 --- so it never hit that particular bug.
>
> i'll fix it.
>
> thanks,
> bill
>
>
> On May 16, 2014, at 2:26 PM, Ken Shen wrote:
>
> > Hi Bill, just realized that my "understanding" last night wasn't quite
> complete (well, obviously).  Here are the lines in get1_derivs in
> net_derivs.f where I think the high temperature factors are getting used:
> >
> >          if (temp > 1d9) then
> >             has_neut = (i1 == neut .or. i2 == neut .or. i3 == neut .or. &
> >                o1 == neut .or. o2 == neut .or. o3 == neut)
> >             if (has_neut) then
> >                call limit_rate_at_very_high_T( &
> >                    ir, rates, rtab, n% factor_neut, n% d_factor_neut_dT,
> .false.)
> >             else if (prot /= 0) then
> >                has_prot = (i1 == prot .or. i2 == prot .or. i3 == prot
> .or. &
> >                   o1 == prot .or. o2 == prot .or. o3 == prot)
> >                if (has_prot) then
> >                   call limit_rate_at_very_high_T( &
> >                      ir, rates, rtab, n% factor_prot, n%
> d_factor_prot_dT, .false.)
> >                end if
> >             else
> >                has_h1 = (i1 == h1 .or. i2 == h1 .or. i3 == h1 .or. &
> >                   o1 == h1 .or. o2 == h1 .or. o3 == h1)
> >                if (has_h1) then
> >                   call limit_rate_at_very_high_T( &
> >                      ir, rates, rtab, n% factor_prot, n%
> d_factor_prot_dT, .false.)
> >                end if
> >             end if
> >          end if
> >
> > If neutrons are not in the isotopes list, "neut" is 0.  Then, has_neut
> should always be false, and factor_neut should never be used.  However,
> when neutrons are not in the list, factor_neut ends up being used for ALL
> reactions, the exact opposite of what I think you intended.
> >
> > I think what is happening is that i3 and o3 are set to 0 at the
> beginning of get1_derivs and stay that way (because no reactions are 3 body
> in these nets).  So when neut == 0, has_neut is always true, regardless of
> the reaction (which should have has_neut = false because there are no
> neutrons).
> >
> > I assume that's not what's supposed to happen in this piece of code,
> right?  Maybe there's supposed to be an "if (neut /= 0)", just like the "if
> (prot /= 0)" a little further down.  And maybe the same thing for h1
> (although I'm confused about prot vs. h1 further up in net_derivs.f).
> >
> > Or am I totally misinterpreting the code??
>
> it is now changed to this
>
>          if (temp > 1d9) then
>             has_neut = (neut /= 0) .and. (i1 == neut .or. i2 == neut .or.
> i3 == neut .or. &
>                o1 == neut .or. o2 == neut .or. o3 == neut)
>             if (has_neut) then
>                call limit_rate_at_very_high_T( &
>                    ir, rates, rtab, n% factor_neut, n% d_factor_neut_dT,
> .false.)
>             else if (prot /= 0) then
>                has_prot = (i1 == prot .or. i2 == prot .or. i3 == prot .or.
> &
>                   o1 == prot .or. o2 == prot .or. o3 == prot)
>                if (has_prot) then
>                   call limit_rate_at_very_high_T( &
>                      ir, rates, rtab, n% factor_prot, n% d_factor_prot_dT,
> .false.)
>                end if
>             else
>                has_h1 = (h1 /= 0) .and. (i1 == h1 .or. i2 == h1 .or. i3 ==
> h1 .or. &
>                   o1 == h1 .or. o2 == h1 .or. o3 == h1)
>                if (has_h1) then
>                   call limit_rate_at_very_high_T( &
>                      ir, rates, rtab, n% factor_prot, n% d_factor_prot_dT,
> .false.)
>                end if
>             end if
>          end if
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20140516/3d1ab45a/attachment.html>


More information about the Mesa-users mailing list