[Mesa-users] using an index I found in another array

Pablo Marchant pamarca at gmail.com
Fri Feb 23 09:41:15 EST 2018

Hi Aldana, various comments to your code snippet

>          if(avg_charge_H(j) .eq. 0.5d0) then

It's almost impossible that the avg_charge will be exactly 0.5. You likely
want to check for a range.

>             loc=j

and because the previous if will never work, loc will never be assigned a
value, it's value will just be random

>                m_50_new= s% m(k .eq. loc) *!s% m(loc) *

Index arrays cannot be logicals, they have to be integers. The reason your
code fails is because loc is never set, in your error message it is a very
large negative number. That is not a valid array index, it has nothing to
do with the content of s% m

               A= 1.3d13 !13.6 ev/hydrogen mass where the energy in erg/gr
>                E_recomb= s% X(loc)*delta_m*A !erg
>                s% extra_heat(k) = E_recomb/ s% nz !erg/g/s

This is wrong, you are computing an energy in ergs, and then dividing it by
the number of zones to get the desired specific energy injection rate. That
does not have the right units, and is completely resolution dependent. Even
if you fix this, I'd urge you to critically evaluate the code you have
written. Even if you fix the units I don't think it would do what you want.

But more importantly, why are you including the recombination energy? MESA
already has this in its EOS (for clarification, see Section 8 of


Pablo Marchant Campos
M.Sc on Astrophysics, Universidad Católica de Chile
PhD on Astrophysics, Argelander-Institut für Astronomie, Universität Bonn
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20180223/8beb9f8f/attachment.html>

More information about the Mesa-users mailing list