With high resolution hydrodynamics simulations, we show that the optimal values of domain radius and grid resolution for the software Sailfish when simulating time-based eccentricity evolution of equal mass, non-circular accreting binaries in a circumbinary disk to be $r_{\rm out} \leq 15a$ and $\delta x / a \le 0.01 $. These values provide a useful guideline for optimizing the performance of simulation runs while maintaining scientific accuracy. Each artificial parameter is probed with 15 runs of 2000 orbits each

Hu, Zhongtian

English

Clemson University: TigerPrints

Clemson University TigerPrints All Theses Theses 8-2023 Simulating the Eccentricity Evolution of Accreting Equal-Mass Binaries: Numerical Sensitivity to the Computational Domain Size and Grid Resolution Zhongtian Hu zhongth@clemson.edu Follow this and additional works at: https://tigerprints.clemson.edu/all_theses  Part of the Cosmology, Relativity, and Gravity Commons, and the Physical Processes Commons Recommended Citation Hu, Zhongtian, "Simulating the Eccentricity Evolution of Accreting Equal-Mass Binaries: Numerical Sensitivity to the Computational Domain Size and Grid Resolution" (2023). All Theses. 4094. https://tigerprints.clemson.edu/all_theses/4094 This Thesis is brought to you for free and open access by the Theses at TigerPrints. It has been accepted for inclusion in All Theses by an authorized administrator of TigerPrints. For more information, please contact kokeefe@clemson.edu. Simulating the eccentricity evolution of accretingequal-mass binaries: numerical sensitivity to thecomputational domain size and grid resolutionA Master’s ThesisPresented tothe Graduate School ofClemson UniversityIn Partial Fulfillmentof the Requirements for the DegreeMaster of PhilosophyPhysicsbyZhongtian HuAugust 2023Accepted by:Dr. Jonathan Zrake, Committee ChairDr. Marco AjelloDr. Dieter HartmannAbstractWith high resolution hydrodynamics simulations, we show that the optimal values of domainradius and grid resolution for the software Sailfish when simulating time-based eccentricity evolutionof equal mass, non-circular accreting binaries in a circumbinary disk to be rout ≤ 15a and δx/a ≤0.01. These values provide a useful guideline for optimizing the performance of simulation runs whilemaintaining scientific accuracy. Each artificial parameter is probed with 15 runs of 2000 orbits each.Keywords: hydrodynamics, accretion, binaries.iiAcknowledgmentsI would like to thank my advisor Dr. Jonathan Zrake for his committed support and guidancethroughout this journal. I would also like to thank Dr. Marco Ajello and Dr. Dieter Hartmann fortheir professional guidance.I would like to thank my parents for their never-fading love and support. I would also liketo thank my significant other, Nicole Maric, for her continuous support and care throughout. I woudlike to thank Madeline Clyburn and Antonio Circiello for their friendship and support.iiiTable of ContentsTitle Page . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiAcknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32.1 The orbital mechanics of binary eccentricity evolution . . . . . . . . . . . . . . . . . 33 Simulation Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73.1 Equations of motion and initial conditions . . . . . . . . . . . . . . . . . . . . . . . . 73.2 Simulation parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95 Conclusion and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15ivChapter 1IntroductionCircumbinary disks have been a source of scientific interest due to their potential influenceon the motions of objects embedded within or around them, such as the motion of binaries. Diskscan be found in young stellar and protoplanetary systems, such as in HD 142527 (Garg et al. (2020),Boehler et al. (2017), and Hunziker et al. (2021)). Binary black holes are expected to form in theaccretion disks of Active Galatic Nuclei, where they are driven to mergers through gas dynamicsand multi-body encounters (Bartos et al. (2017), Antonini and Rasio (2016), Tagawa et al. (2020)).AGN themselves could also be the result of binary mergers from supermassive black holes (SMBHs)during galactic mergers (Begelman et al., 1980). Since gravitational interaction is the primary way inwhich the binary influences the disk and vice versa, such as through momentum and mass transfervia binary accretion, there has been numerous studies conducted on this topic using numericalsimulations (e.g. Macfadyen and Milosavljevic (2006), Shi et al. (2012), and D’Orazio et al. (2012)).When looking at long term orbital evolution of eccentric binaries in circumbinary disks, Muñozet al. (2019) has found that the orbits of eccentric binaries settles around eb ≈ 0.5 if their initialeccentricity is eb ≥ 0.1 and eb ≈ 0 if their initial eccentricity is eb ≤ 0.1, respectively. Zrake et al.(2021) confirmed the existence of the eccentricity ”attraction” from Muñoz et al. (2019) and refinedthe values to eb ≈ 0.45 and the lower threshold e at eb ≈ 0.08 for equal mass eccentric binaries.D’Orazio and Duffell (2021) expanded on studies from Zrake et al. (2021) with results on binarysemi-major axis evolution, finding that eb ≤ 0.1 binaries circularize with their semi-major axesexpanding, while those that tend to eb ≈ 0.5 have seen their semi-major axes decay. Siwek et al.(2023) studied the eb evolution in unequal mass eccentric binaries and found that all large mass-ratio1binaries evolved to steady state with eccentricities attracted to eb ≈ 0.5, and binaries with q ≤ 0.3went to coalescence, while smaller mass-ratio binaries expanded. For a comprehensive review oncircumbinary disk accretion, please see the review from Lai and Muñoz (2022).The purpose of this study is to probe the simulation code Sailfish’s sensitivity to the artificialparameters domain radius and grid effective resolution. These parameters affect the performance ofsimulations since they determine the number of cells computers have to process. Finding optimalvalues that let the simulation give converging results while minimizing performance workload canhelp further improve the code as well as optimization.2Chapter 2Background2.1 The orbital mechanics of binary eccentricity evolutionThe definition of the orbital eccentricity of a gravitational binary system is:eb =√1 +2ϵh2µ2, (2.1)where ϵ is the specific orbital energy; h2 is the magnitude squared of the specific orbital angularmomentum; µ is the standard gravitational parameter GM , with M being the reduced mass of thebinary. It is a dimensionless parameter that characterizes how circular the orbit is: e = 0 means theorbit is a perfect circle; 0 < e < 1 means the orbit is an ellipse, and so on. To turn this equationinto a more useful form, we first examine the total energy of the binary without the disk when thebinary components are at pericenters or apocenters:E = −µr+h22Mr2, (2.2)where it can be rewritten as:3E +µr− h22Mr2= 0orr2 +µEr − h22mE= 0.The negative sum of the solutions to this quadratic equation is the coefficient of the µE r term, hencewe can define the semi-major axis of the orbit as:a = − h2E. (2.3)The eccentricity of the orbit thus can be rewritten with the definition of semi-major axis as:e =√1− h2µa. (2.4)To study the eccentricity evolution of the binary in a disk, the eccentricity vector associatedwith the binary orbit must be used, and we take into account the gravitational interaction betweenthe binary and all the gas parcels present. To derive the eccentricity vector, we first consider agravitational central force as the time derivative of the total momentum of the system:ṗ = − µr2r̂. (2.5)Knowing that the eccentricity of an orbit can be seen as how oblate the orbit is, we want to find thevector that is parallel to the semi-major axis:ṗ× L = −µ(r̂ × (r× ṙ)). (2.6)With a bit of vector manipulation, the equation becomes:ṗ× L = − µr3(rrṙ − r2ṙ)= −µ(rṙr2− ṙr).(2.7)Noting the conservation of angular momentum, we can rewrite the LHS as taking the time derivative4of p× L:ddt(p× L) = ddtµrr, (2.8)we immediately arrive at the defintion of the Laplace-Runge-Lenz vector by integrating both sideswith respect to time:A = p× L− µrr. (2.9)Dividing both sides by µ, the eccentricity vector of the binary orbit is thus:e =v × hµ− rr. (2.10)where h = r × v, and v = ṙ2 − ṙ1, where r1 and r2 are the respective positions of the binarycomponents with respect to the origin. In a conserved orbit, the eccentricity vector has a fixeddirection pointing along the semi-major axis, with its magnitude being the eccentricity of the orbit.In this study, however, the circumbinary disk takes angular momentum away from the binary, whichresults in a rate of change on the eccentricity vector’s direction and magnitude. For this project,the change in magnitude is explored.The rate of change of the binaries orbital eccentricity can be determined by taking the timederivative of the eccentricity vector we just found:ė =ddt(v × hµ− rr)=v̇ × h+ v × ḣµ− ( ṙr− rr2ṙ)=v̇ × h+ v × ḣµ− v⊥r(2.11)where v⊥ is the radial velocity. By rewriting the force terms as v̇ = − µr2 r̂ + fgas and µ = GM , theradial velocity cancels out with the gravitational force term, and we are left with:GM ė = fgas × (r× v) + v × (r× fgas)−·MM[v × (r× v)] (2.12)where M = M1 +M2 is the total mass of the binary; fgas is the gravitational force of the gas in the5circumbinary disk on the binary; Ṁ is the accretion of mass from the disk by the binary; r × v isthe specific angular momentum, where r = r2 − r1 and v = ṙ2 − ṙ1. We use this form to calculatethe ė evolution.6Chapter 3Simulation Setup3.1 Equations of motion and initial conditionsTo probe the sensitivity of the binary eccentricities to artificial parameters, we use our hy-drodynamics numerical simulation application named Sailfish, a high order GPU-powered Godunovcode and the spiritual successor of Mara3 (Zrake et al., 2021), to simulate a coupled, equal-massbinary that sits in a locally isothermal circumbinary disk in a two dimensional, Cartesian staticallyrefined mesh. Sailfish numerically solves the 2D vertically averaged Navier-Stokes equations withconstant kinematic viscosity ν:∂Σ∂t+∇ · (Σv) =·Σsink (3.1)∂Σv∂t+∇ · (Σvv + P I−Tvis) =·Σsinkv + F g (3.2)which are described in full in Tiede et al. (2020). The binaries accrete material from the disk byacting as sinks of mass and momentum, which are denoted by Σ̇ in the equations; v is the gasvelocity, and P is the vertically integrated gas pressure Tiede et al. (2020). The prescription for thesink is acceleration free, as seen in the Σ̇sinkv term in equation (2.1). We note from Dittmann andRyan (2021) and Dempsey et al. (2020) that implementing a torque-free sink model can significantlydecrease the dependence of density profiles around the sinks on the sink rate, but since only con-verging runs are conducted in this study, the results don’t share major differences. The thickness7of the disk is described by the orbital Mach number parameter M = (hr )−1, which comes from thevertically averaged thin disk approximation Tiede et al. (2020). For this study, the thickness of thedisks in all simulations are kept at M = 10.3.2 Simulation parametersThe two parameters studied are the size rd (in units of semi-major axis a) and the gridspacing ∆x of the computational domain. We measure the binary eccentricity evolution as a functionof the domain size rd, holding fixed the grid spacing at ∆x = 0.01a. When studying the gridresolution parameter, the domain radius is held constant at rd = 15a while the grid spacing is variedby selecting values of ∆x which correspond to 50, 100, 150, 200, and 250 zones per semi-major axis a.The derivative ė of the binary orbital eccentricity is measured as a function of e, the domain radiusrd, and the grid spacing ∆x, in order to probe how sensitive ė(e) is to the numerical parameters rdand ∆x. Each batch of runs with fixed ∆x or fixed rd consists of ė measurements taken at threebinary orbital eccentricities: e = 0.2, 0.4, 0.6. A total of 30 runs were conducted.8Chapter 4ResultsFigure 4.1 shows what the density profile of the inner region of the circumbinary disk lookslike while the binary interacts with it. With non-circular binaries, material from the inner edge ofthe disk gets flung to the other side of the central cavity, creating this ”seashell” density profile asshown in the figure. Figure 4.2 shows the results of the 30 runs conducted. Both domain radius(rout) varying runs and effective resolution (δx/a) varying runs are shown. The rout runs seems toindicate that a domain radius of 15a is needed to resolve the binary’s eccentricity evolution. Theδx/a runs show that a 0.01 ∆xa grid resolution is required, but the general trend for each eccentricityis inconsistent, meaning eccentricity is sensitive to this parameter. The time evolution data selectedfor these runs begin at t = 1500 orbits and ending at t = 2000 orbits, which is when the system hasevolved through 3 viscous relaxation times (Zrake et al., 2021). The following table contains thedata from the two blocks of runs:eccentricitydomain radius a5 10 15 20 250.2 7.70 8.06 7.92 8.61 13.440.4 4.12 -2.81 0.68 1.12 1.740.6 3.07 -2.39 -2.45 -2.68 -2.519eccentricityeffective resolution ∆xa0.02 0.01 0.0067 0.005 0.0040.2 11.67 10.22 12.10 8.92 9.380.4 5.01 -0.07 -0.72 -2.63 2.430.6 -3.20 -2.46 -2.25 -2.14 -2.061015 10 5 0 5 10 15151050510150.20.40.60.81.01.21.41.6/Users/jackhu/Documents/Astro_Research/sailfish_data/data/sensitivity_tests/grid-resolution-varying/e0.6_gr250.pkFigure 4.1: Example of simulation. X and Y axes are the domain axes, while the color bar indicatesthe density of material present at each cell.115 10 15 20 25dr5051015ėṀdomain radius varying, ∆xa= 0.01e = 0.2e = 0.4e = 0.60.01004 × 10 3 6 × 10 3 2 × 10 2∆xa5051015ėṀgrid effective resolution varying, dr = 15Binary components’ eccentricity evolutionFigure 4.2: Results of e-dot over m-dot as a function of domain size and effective grid resolution.12Chapter 5Conclusion and DiscussionWe have probed the sensitivity of the hydrodynamics code Sailfish to the artificial pa-rameters domain radius and effective resolution by examining the time-evolution signaturesof eccentricity of equal mass, accreting binaries in different initial eccentricity configurations andcomparing with known converging results from Zrake et al. (2021). The range of eccentricities weree = [0.2, 0.4, 0.6], and, at each eccentricity, the domain radii and the grid resolution values probedwere rout = [5, 10, 15, 20, 25] and δx/a = [0.02, 0.01, 0.0067, 0.005, 0.004], respectively. Each run wasran until t = 2000 orbits. The domain size runs appeared to have converged with a sufficiently largevalue. However, we found eccentricity evolution of the binary is sensitive to the effective grid resolu-tion, especially for the case of e = 0.4 and e = 0.6. The sensitivity seems to introduce random errorrather than a consistent trend, so we conclude that even the lowest resolution we performed, withδx/a = 0.02 yields a measurement in the right ballpark, however we caution that such measurementsmay carry an uncertainty on the other of 10 − 20% associated with the random error arising fromgrid resolution effects.The large domain size of 15a, that seems to be required for a consistent measurement ofthe eccentricity time derivative, indicates that an imposed outer boundary condition can influencemeasurements of the orbital evolution of the binary even when the domain is large. We understandthis is because the surface density in the disk controls the gravitational force exerted on the binary,and that the outer boundary condition can introduce spurious features in the mass distributionof the circumbinary disk when it is placed too close to the binary. This is because the boundarycondition imposes an axial symmetry at a finite radius, so unless the disk extends far enough from13the binary to erase the effects of the time-varying tidal field, the boundary condition can lead tosuppression of real non-axisymmetric disk features that influence the orbital evolution. For thesereasons, we recommend that the domain size of future simulations of binaries with disks be set to15a or larger.14BibliographyAntonini, F. and Rasio, F. A. (2016). MERGING BLACK HOLE BINARIES IN GALACTICNUCLEI: IMPLICATIONS FOR ADVANCED-LIGO DETECTIONS. The Astrophysical Journal,831(2):187.Bartos, I., Kocsis, B., Haiman, Z., and Márka, S. (2017). Rapid and bright stellar-mass binary blackhole mergers in active galactic nuclei. The Astrophysical Journal, 835(2):165.Begelman, M. C., Blandford, R. D., and Rees, M. J. (1980). Massive black hole binaries in activegalactic nuclei. Nature, 287(5780):307–309.Boehler, Y., Weaver, E., Isella, A., Ricci, L., Grady, C. A., Carpenter, J. M., and Pérez, L. M.(2017). A close-up view of the young circumbinary disk hd 142527. The Astrophysical Journal,840.Dempsey, A. M., Muñoz, D. J., and Lithwick, Y. (2020). Inner boundary condition in quasi-lagrangian simulations of accretion disks. The Astrophysical Journal Letters, 892.Dittmann, A. J. and Ryan, G. (2021). Preventing anomalous torques in circumbinary accretionsimulations. The Astrophysical Journal, 921.D’Orazio, D. and Duffell, P. (2021). Orbital evolution of equal-mass eccentric binaries due to a gasdisk: Eccentric inspirals and circular outspirals. The Astrophysical Journal Letters, 914:null.D’Orazio, D. J., Haiman, Z., and Macfadyen, A. (2012). Accretion into the central cavity of acircumbinary disc. Monthly Notices of the Royal Astronomical Society, 436:2997–3020.Garg, H., Pinte, C., Christiaens, V., Price, D. J., Lazendic, J. S., Boehler, Y., Casassus, S., Marino,S., Pérez, S., and Zuleta, A. (2020). Non-keplerian spirals, a gas-pressure dust trap, and aneccentric gas cavity in the circumbinary disc around hd 142527. Monthly Notices of the RoyalAstronomical Society, 504:782–791.Hunziker, S., Schmid, H. M., Ma, J., Ménard, F., Avenhaus, H., Boccaletti, A., Beuzit, J.-L.,Chauvin, G., Dohlen, K., Dominik, C., Engler, N., Ginski, C., Gratton, R., Henning, T., Langlois,M., Milli, J., Mouillet, D., Tschudi, C., van Holstein, R. G., and Vigan, A. (2021). Hd 142527:quantitative disk polarimetry with sphere. Astronomy & Astrophysics, 648:A110.Lai, D. and Muñoz, D. J. (2022). Circumbinary accretion: From binary stars to massive binaryblack holes.Macfadyen, A. and Milosavljevic, M. (2006). An eccentric circumbinary accretion disk and thedetection of binary massive black holes. The Astrophysical Journal, 672:83 – 93.Muñoz, D., Lai, D., Kratter, K., and Miranda, R. (2019). Circumbinary accretion from finite andinfinite disks. The Astrophysical Journal, 889:null.15Shi, J.-M., Krolik, J. H., Lubow, S. H., and Hawley, J. F. (2012). Three-dimensional magnetohy-drodynamic simulations of circumbinary accretion disks: Disk structures and angular momentumtransport. The Astrophysical Journal, 749(2):118.Siwek, M., Weinberger, R., and Hernquist, L. (2023). Orbital evolution of binaries in circumbinarydisks.Tagawa, H., Haiman, Z., and Kocsis, B. (2020). Formation and evolution of compact-object binariesin agn disks. The Astrophysical Journal, 898(1):25.Tiede, C., Zrake, J., MacFadyen, A., and Haiman, Z. (2020). Gas-driven inspiral of binaries in thinaccretion disks. The Astrophysical Journal, 900(1):43.Zrake, J., Tiede, C., MacFadyen, A., and Haiman, Z. (2021). Equilibrium eccentricity of accretingbinaries. The Astrophysical Journal Letters, 909(1):L13.16

Simulating the Eccentricity Evolution of Accreting Equal-Mass Binaries: Numerical Sensitivity to the Computational Domain Size and Grid Resolution

https://tigerprints.clemson.edu/cgi/viewcontent.cgi?article=5109&context=all_theses

Simulating the Eccentricity Evolution of Accreting Equal-Mass Binaries: Numerical Sensitivity to the Computational Domain Size and Grid Resolution

Abstract

Similar works

Full text

Available Versions

Clemson University: TigerPrints