[Mesa-users] some problems in using Mesa to calculate wind transfer

Pablo Marchant pamarca at gmail.com
Thu Dec 12 10:48:07 EST 2019


Great,

but again, a reminder for next time, keep mesa-users on your replies. Even
if your issue is resolved, that way everyone else knows it is. Remember to
hit reply-all!

Cheers

On Thu, Dec 12, 2019 at 2:08 PM Lingrui <llrw at 163.com> wrote:

> Thanks a lot again!
> It is my good luck that you answer my question with such patience. And I
> finally understand the 2014work!
> Wish you good!
>
>
>
> Sent from my Mi phone
> On Pablo Marchant <pamarca at gmail.com>, Dec 12, 2019 19:59 wrote:
>
> Hi Lingrui,
>
> as Warrick said, please keep mesa-users on the loop. Kept the discussion
> going but missed mesa-users had been taken out a couple of posts above.
>
> You're trying to reproduce the results of Ablimit et al. 2014, but I don't
> think you're interpreting that work correctly. Or perhaps what you're
> trying to use in MESA is not what you expect. The wind mass transfer
> routines are meant to model when a fraction of the wind of a star is
> captured by its companion and accreted into it. The model you're trying to
> reproduce deals with an irradiated star, that begins to lose mass from
> that. Irradiation comes from accretion into a compact object in the system.
> This is the sort of thing that you need to implement as a wind on the star,
> not as wind accretion. So you should look for the other_wind routines for
> single stars instead.
>
> Cheers
>
> On Thu, Dec 12, 2019 at 11:43 AM Lingrui Lin <llrw at 163.com> wrote:
>
>>
>> Thanks a lot.
>> I am trying to use MESA to take into accout the two models in ApJ_780_80,
>> which is attached.
>>
>>
>> At 2019-12-12 18:20:42, "Pablo Marchant" <pamarca at gmail.com> wrote:
>>
>> Thanks Lingrui,
>>
>> I just tested this and noticed you have not specified a mass loss
>> prescription for the donor star in the system. I'm not sure what you mean
>> when you say that b%s_donor%mstar_dot is not zero, it is indeed zero at the
>> beginning of the simulation. That value will be non zero when the system
>> initiates mass transfer, but that mdot will not correspond to winds, but
>> rather to loss of mass through L1. So the routine you're using won't affect
>> that.
>>
>> It would be useful if you describe what exactly you're trying to achieve.
>> Are you just trying to use a different model for mass transfer efficiency
>> rather than a model for wind accretion?
>>
>> Cheers
>>
>> On Thu, Dec 12, 2019 at 11:02 AM Lingrui Lin <llrw at 163.com> wrote:
>>
>>>
>>> Oh, thanks a lot.
>>> attached is my inlists and run_star, run_binary_extras.
>>> Hope they can run properly.
>>> Thank you!
>>> Sincerely
>>>
>>>
>>>
>>>
>>> At 2019-12-12 17:55:35, "Pablo Marchant" <pamarca at gmail.com> wrote:
>>>
>>> Hi Lingrui,
>>>
>>> Please include all necessary files (inlists, run_star_extras and
>>> run_binary_extras) as attachments rather than pasting them into the email.
>>> That ensures us we are using exactly your setup.
>>>
>>> Cheers
>>>
>>> On Thu, Dec 12, 2019, 10:25 AM Lingrui Lin <llrw at 163.com> wrote:
>>>
>>>>
>>>> oh, I am so sorry. I am new to mesa and new to do research
>>>> here is my inlist_project:
>>>> &binary_job
>>>>
>>>>    inlist_names(1) = 'inlist1'
>>>>    inlist_names(2) = 'inlist2'
>>>>
>>>>    evolve_both_stars = .false.
>>>>    change_ignore_rlof_flag = .true.
>>>> / ! end of binary_job namelist
>>>>
>>>> &binary_controls
>>>>
>>>>    m1 = 1.0d0  ! donor mass in Msun
>>>>    m2 = 1.2d0 ! companion mass in Msun
>>>>    initial_period_in_days = 1d0
>>>>
>>>>    fr = 0.05
>>>>    fr_limit = 1d-2
>>>>
>>>>    !do_enhance_wind_1 = .true.
>>>>    !do_enhance_wind_2 = .true.
>>>>
>>>>    !accretion_prowered_irradiation = .true.
>>>>    do_wind_mass_transfer_1= .true.
>>>>    limit_retention_by_mdot_edd = .true.
>>>>    use_other_binary_wind_transfer = .true.
>>>> / ! end of binary_controls namelist
>>>>
>>>>
>>>> and here is my other_wind subroutine:
>>>> and of cource I have set  b% other_binary_wind_transfer =>
>>>> wind_transfer_routine
>>>> The print* in the subroutine work all properly
>>>>  subroutine wind_transfer_routine(binary_id,s_i,ierr)
>>>>               integer,intent(in)::binary_id,s_i
>>>>               integer,intent(out)::ierr
>>>>
>>>> real(dp)::etas,etaa,phi,etahydr,etahe,lim_mtransdot,C1,loghe
>>>>               type(binary_info),pointer::b
>>>>               ierr=0
>>>>               call binary_ptr(binary_id,b,ierr)
>>>>               if(ierr/=0)then
>>>>                       write(*,*)'failed in binary_ptr'
>>>>                       return
>>>>               end if
>>>>               etahydr=0.15
>>>>               phi = 0.5
>>>>               etas=0.9
>>>>               etaa=1
>>>>               if (b%m(b%d_i).lt.b%m(b%a_i)) then
>>>>                    b%wind_xfer_fraction(b%d_i) = &
>>>>                    -b%s_donor%mstar_dot*secyer*1d7/3.5**2* &
>>>>                    b%m(b%d_i)**(-5/3)* &
>>>>                    (b%m(b%a_i)+b%m(b%d_i))**(2/3) &
>>>>                    /(etas*etaa)/phi**2
>>>>               else
>>>>                    b%wind_xfer_fraction(b%d_i) = &
>>>>                    -b%s_donor%mstar_dot*secyer*1d7/3.5**2* &
>>>>                    b%m(b%d_i)**(-1.9)* &
>>>>                    (b%m(b%a_i)+b%m(b%d_i))**(2/3) &
>>>>                    *b%m(b%a_i)**0.24/(etas*etaa)/phi**2
>>>>               end if
>>>>               print*,'b%m(b%d_i)/Msun',b%m(b%d_i)/Msun
>>>>               print*,'b%m(b%a_i)/Msun',b%m(b%a_i)/Msun
>>>>               print*,'lim_mtransdot000',lim_mtransdot
>>>>               lim_mtransdot=7.2d-6*(b%m(b%a_i)/Msun-0.6)
>>>>               print*,'lim_mtransdot',lim_mtransdot
>>>>               !C1=3
>>>>               if (-b%mtransfer_rate*secyer/Msun &
>>>>                    .gt.lim_mtransdot) then
>>>>                    print*,'critical accretion'
>>>>                    b%mdot_wind_transfer(b%d_i)=-lim_mtransdot* &
>>>>                    Msun/secyer
>>>>               end if
>>>>
>>>>                loghe=log10(-b%mtransfer_rate*secyer/Msun* &
>>>>                            b%s_donor%center_he4)
>>>>                    !print*,'centerhe4',b%s_donor%center_he4
>>>>                    print*,'loghe',loghe
>>>>                    print*,'b%mdot_wind_transfer(b%d_i)', &
>>>>                     b%mdot_wind_transfer(b%d_i)
>>>>                    print*, 'b%mtransfer_rate',b%mtransfer_rate
>>>>                    print*,'b%s_donor%mstar_dot', &
>>>>                     b%s_donor%mstar_dot
>>>>                    print*,'b%wind_xfer_fraction(b%d_i)', &
>>>>                     b%wind_xfer_fraction(b%d_i)
>>>>                if(loghe.gt.-7.06d0.and.loghe.lt.-5.95d0) then
>>>>                        etahe=0.54*loghe+4.16
>>>>                else if (loghe.gt.-5.95d0.and.loghe.lt.-5.76d0)then
>>>>                        etahe=0.54*(loghe+5.6)**2+1.01
>>>>                else
>>>>                        etahe=1
>>>>                end if
>>>>               end subroutine wind_transfer_routine
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>
>>>
>>>
>>>
>>
>>
>> --
>> Pablo Marchant Campos
>> M.Sc on Astrophysics, Universidad Católica de Chile
>> PhD on Astrophysics, Argelander-Institut für Astronomie, Universität Bonn
>>
>>
>>
>>
>>
>
>
> --
> Pablo Marchant Campos
> M.Sc on Astrophysics, Universidad Católica de Chile
> PhD on Astrophysics, Argelander-Institut für Astronomie, Universität Bonn
>
>

-- 
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/20191212/471ad341/attachment.htm>


More information about the Mesa-users mailing list