[mesa-users] a little bug

Ehsan Moravveji e.moravveji at gmail.com
Fri Aug 26 15:40:50 EDT 2016


> On 26 Aug 2016, at 20:28, Bill Paxton <paxton at kitp.ucsb.edu> wrote:
> 
> 
> On Aug 26, 2016, at 11:06 AM, Ehsan Moravveji wrote:
> 
>> Thanks Bill for the encouragement. 
>> 
>> Before you close the mlt_info.f90 file on your editor, you may perhaps consider the following two cases as well:
>> 
>> - At lines 983-4, the two derivatives of grada are missing:
>> 
>>                s% d_gradT_dlnR(k) = one_m_f*s% d_gradr_dlnR(k)
>>                s% d_gradT_dL(k) = one_m_f*s% d_gradr_dL(k)
>> 
>> I thought, perhaps, the more consistent differential of gradT requires propagating that of grada too:
>> 
>>                s% d_gradT_dlnR(k) = f * d_grada_dlnR + one_m_f*s% d_gradr_dlnR(k)
>>                s% d_gradT_dL(k) = f* d_grada_dL + one_m_f*s% d_gradr_dL(k)
>> 
>> Am I right or wrong here?
> 
> Check on how grada is computed -- and see what you can find about d_grada_dlnR and d_grada_dL.
> 
You’re absolutely right; I overlooked the fact that grada has a vanishing derivative w.r.t. any other variable, except itself.
Thanks again for this pointer. Indeed, the d_gradT_dlnR and d_gradT_dL are correctly implemented.

Have a nice weekend,
Ehsan.

> 
> 
> 
> 
>> 
>> - Next stop is at line 922:
>> 
>>                   s% d_gradT_dlnTm1(k) = &
>>                      f*beta*s% d_eos_dlnT(i_grad_ad, k-1) + one_m_f*s% d_gradr_dlndm1(k)
>> 
>> would better be:
>> 
>>                   s% d_gradT_dlnTm1(k) = &
>>                      f*beta*s% d_eos_dlnT(i_grad_ad, k-1) + one_m_f*s% d_gradr_dlnTm1(k)
>> Perhaps, another example of how cut-paste can get us wrong. 
> 
> Yup.  that's 2 stars for you in one day.  Good job!
> 
> thanks,
> Bill
> 
> 
> 
> 
> 
> 
> 
> 
>> 
>> Thanks again,
>> Ehsan.
>> 
>> 
>>> On 26 Aug 2016, at 18:41, Bill Paxton <paxton at kitp.ucsb.edu <mailto:paxton at kitp.ucsb.edu>> wrote:
>>> 
>>> 
>>> On Aug 26, 2016, at 6:49 AM, Ehsan Moravveji wrote:
>>> 
>>>> Hi Bill,
>>>> 
>>>> I hope this email finds you well.
>>>> I just found a little bug in mlt_info.f90, line 978:
>>>>                 s% d_gradT_dlnTm1(k) = s% d_gradr_dlndm1(k)
>>>> should be
>>>>                 s% d_gradT_dlnTm1(k) = s% d_gradr_dlnTm1(k)
>>>> 
>>>> MLT is important, and interesting for me ;-)
>>>> 
>>>> Kind regards,
>>>> Ehsan.
>>>> 
>>> 
>>> 
>>> Hi Ehasn,
>>> 
>>> Thanks!  I've fixed it now.  The dangers of cut-and-paste editing is that sometimes things get forgotten -- like this one.
>>> 
>>> I'm Cc'ing mesa-users in hopes that your effort will inspire others to find and report bugs!  ;)
>>> 
>>> Cheers,
>>> Bill
>>> 
>>> 
>>> ------------------------------------------------------------------------------
>>> _______________________________________________
>>> mesa-users mailing list
>>> mesa-users at lists.sourceforge.net <mailto:mesa-users at lists.sourceforge.net>
>>> https://lists.sourceforge.net/lists/listinfo/mesa-users <https://lists.sourceforge.net/lists/listinfo/mesa-users>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20160826/3917286f/attachment.html>


More information about the Mesa-users mailing list