[Mesa-users] convergence problem

[sheyda najafi]

Dear MESA users

I have changed extra_dpdm to add a new acceleration term to HSE
equation in MESA code.
(p(i-1)- p(i))/dm_bar(i)=( -standard_cgrav*m(r))/(4*pi* s% r(i)**4) -
beta*(clight**2) *rho_face(i) (dphi/dm_bar)
therefore
extra_dpdm(i)= - beta*(clight**2) *rho_face(i) dphi/dm_bar

which i have asked about the quantities and the shape of this new term
in previous emails.

 My question is how does this extra acceleration change the pressure
scale height in the code. I can not find the relation. I know that
pressure scale height depends on p_face and rho_face and grav. The
extra_dpdm has impact on density and pressure but I can not find the
relation between grav and extra_dpdm in the code.

Thank you very much
Sheyda Najafi

On 5/7/19, sheyda najafi <sh.najafi8590 at gmail.com> wrote:
> Dear Mesa users
>
> I have another question regarding my previous emails.
>
> As i  said I wanted to study the evolution of 1 solar mass star in the
> presence of
> new extra acceleration term :(-beta*(clight**2)*(dphi/dr)).
>  i have changed the
> extra_dPdm subroutine of MESA code to add
> the extra term (-beta*rho(clight**2)*(dphi/dm)) ==>( in lagrangian
> framework) .
>  where phi depends on local density of each layer(cell density: s% rho).
>  which is named (xn) in the code. and beta is a constant.
> i have added the extra term as follows
>
>            do k = 1,2
>              s% extra_dPdm(k) = -( ((   xn(k)-    xn(k+1))/(s%
> dm_bar(k+1)))*((beta*(clight**(2.0d0)))*(s% rho_face(k))) )
>
>          end do
>
>          do k = 2, s% nz-1
>           s% extra_dPdm(k) = -( (((   xn(k)-    xn(k+1))/(s% dm_bar(k+1)))+
> &
>          ((   xn(k-1) -   xn(k))/(s% dm_bar(k))))
> *((beta*(clight**(2.0d0)))*(s% rho_face(k))*(1.d0/2.d0)) )
>
>          end do
>
>          do k = s% nz-1, s% nz
>           s% extra_dPdm(k) = -( ((   xn(k-1) -    xn(k))/(s%
> dm_bar(k)))*((beta*(clight**(2.0d0)))*(s% rho_face(k))) )
>
>         end do
>
> I noticed that  s% dm_bar(s%nz) = s%dm(s% nz) + s%dm(s%nz-1)/2
> and not (s% dm(s% nz) + s% dm(s%nz-1))/2.
>
> i am not sure that dm_bar must be used in extra_dpdm or  (s% dm(s% nz)
> + s% dm(s%nz-1))/2.
>
> Thank you very much
> Sheyda Najafi
>
>
> On 10/21/18, sheyda najafi <sh.najafi8590 at gmail.com> wrote:
>> Dear Josiah and MESA users
>>
>>
>> I have another question regarding my previous emails. i have asked it
>> before but unfortunately no one answered.
>>
>> As i  said I wanted to study the evolution of 1 solar mass star in the
>> presence of
>> new extra acceleration term :(-beta*(clight**2)*(dphi/dr)).
>>  i have changed the
>> extra_dPdm subroutine of MESA code to add
>> the extra term (-beta*rho(clight**2)*(dphi/dm)) ==>( in lagrangian
>> framework) .
>>  where phi depends on local density of each layer(cell density: s% rho).
>>  which is named (xn) in the code. and beta is a constant.
>>
>>            do k = 1,2
>>              s% extra_dPdm(k) = -( ((   xn(k)-    xn(k+1))/(s%
>> dm_bar(k+1)))*((beta*(clight**(2.0d0)))*(s% rho_face(k))) )
>>
>>          end do
>>
>>          do k = 2, s% nz-1
>>           s% extra_dPdm(k) = -( (((   xn(k)-    xn(k+1))/(s%
>> dm_bar(k+1)))+
>> &
>>          ((   xn(k-1) -   xn(k))/(s% dm_bar(k))))
>> *((beta*(clight**(2.0d0)))*(s% rho_face(k))*(1.d0/2.d0)) )
>>
>>          end do
>>
>>          do k = s% nz-1, s% nz
>>           s% extra_dPdm(k) = -( ((   xn(k-1) -    xn(k))/(s%
>> dm_bar(k)))*((beta*(clight**(2.0d0)))*(s% rho_face(k))) )
>>
>>         end do
>>
>> Since i did not know the value of s%rho(0) i did'nt have xn(0) to
>> write extra_dpdm it in 1 part therefore i have split it in three
>> parts. how it is possible to write this term in 1 part. i mean if i
>> want to write the extra term as :
>>
>>          do k = 1, s% nz
>> s% extra_dpdm(k)= -((   xn(k-1)-    xn(k)) / s% dm_bar(k)) * beta*
>> (clight**2)* ( s%rho_face(k)).
>>         end do
>>
>> I will have problem in finding the value of xn(0). because xn is
>> related to cell density s%rho(0) via
>>
>>    xn(i)= (( n1*(lambda**(4.d0+n1)) )/(((mpl)**(n1))*clight*beta*s%
>> rho(i)))**(1.d0/(n1+1.d0))
>>
>> lambda, mpl, n1 and beta are constant.
>>
>>  which density should i use as rho(0) ??
>>
>> (for your kindly information my codes are attached in the email)
>> Thank you in advance
>>
>>
>> Sheyda
>>
>