danieljprice / phantom

Phantom Smoothed Particle Hydrodynamics and Magnetohydrodynamics code
https://phantomsph.github.io
Other
103 stars 223 forks source link

Running difference between using setup_grtde.f90 and setup_star.f90 #399

Closed Shigaharuki3012 closed 1 year ago

Shigaharuki3012 commented 1 year ago

When I simulated TDE using _setup_grtde.f90_,an error "FATAL ERROR! eos: utherm < 0 : u = -1.148E-10" always appeared when t reaches 800 (tmax=8.0E+04 and dtmax=160.0 were setted ),while the error would not appear when using "star" setup and the command "make moddump MODFILE=moddump_tidal.f90" .In the latter case the program evoluted normally.I don't know whether it's a bug or not.

To note ,the input parameters were the same when using the two settings.In "grtde" setup they are as follows(): mhole = 1.000E+06 ! mass of black hole (solar mass) mstar = 1.000 ! mass of star (solar mass) rstar = 1.000 ! radius of star (solar radii) stardensprofile = 1 ! star density profile (1=adiabatic, 2=kepler) beta = 5.000 ! penetration factor ecc = 0.9500 ! eccentricity (1 for parabolic) norbits = 5.000 ! number of orbits dumpsperorbit = 100 ! number of dumps per orbit nr = 50 ! particles per star radius (i.e. resolution) theta = 60.000 ! inclination of orbit (degrees) ~

while in "star" setup case,the parameters in tde.tdeparams are(all code units are setted as solar units): beta = 5.000 ! penetration factor mh = 1.000E+06 ! mass of black hole (code units) ms = 1.000 ! mass of star (code units) rs = 1.000 ! radius of star (code units) theta = 60.000 ! stellar rotation with respect to y-axis (in degrees) r0 = 780.000 ! starting distance (code units) ecc = 0.9500 ! eccentricity of stellar orbit (1 for parabolic) use binary = F ! true/false

So besides the same mass and radius of the BH and star,theta,e,beta are also the same. Note that the r0 above is calculated from the apocentre relation "Ra=Rp(1+e)/(1-e)" and transferred into solar radius.

danieljprice commented 1 year ago

A big difference is that you can relax the star with SETUP=star but with grtde the star is not relaxed at all

Shigaharuki3012 commented 1 year ago

So we aren't able to use a relaxed star when using "grtde" setting. But I don't know whether the fatal error is related to just whether the star is relaxed and the default parameters are overwritten. In fact I guessed it'd be no problem if I just changed some of the basic parameters:)Maybe something wrong with my other settings? Anyway thank you for your reply.

danieljprice commented 1 year ago

if you compile SETUP=star without also specifying GR=yes then you are also using Newtonian gravity instead of GR. Relaxing a star in the GR code is a bit more complicated as you have to relax in the Minkowski metric then compile the moddump with METRIC=Kerr

Shigaharuki3012 commented 1 year ago

if you compile SETUP=star without also specifying GR=yes then you are also using Newtonian gravity instead of GR. Relaxing a star in the GR code is a bit more complicated as you have to relax in the Minkowski metric then compile the moddump with METRIC=Kerr

Thank you for your reply.Please forgive my lack of knowledge in the area of Makefile and my carelessness :(...Now I know how to deal with it. But I noticed that there is a "grstar" setup,where METRIC has been setted as Minkowski initially when wanting to relax the star.So maybe I can directly use "grstar" setup,after which I can then compile the moddump with METRIC=Kerr?

danieljprice commented 1 year ago

yes, this is what grstar is for. The plan is to make all of this automated in the near future to enable a single "relax-and-place-into-orbit" like grtde, but this requires a bit of other work first (namely #358 so we can switch the metric from Minkowski to Kerr at runtime).

We also need some better documentation on the whole procedure, it would be great if you might consider contributing this -- just add an rst file to the docs directory :)

Shigaharuki3012 commented 1 year ago

I will have a try after I've supplemented my knowledge.Thank you for your patient help and reply:)