Exploring the Jet Formation in Binary Systems Applying 3D MHD Simulations

We investigate the formation of an ideal magnetized jet that originates from a disk acting as a boundary by conducting axisymmetric MHD simulations. Our simulations demonstrate that the magnetized jet is consistently launched and reaches a stable state. To further advance our study, we extended the model setup to three dimensions and performed 3D MHD simulations of the jet launched from a disk surface, achieving a stable and appropriate model setup. Additionally, we expanded our study by incorporating the companion star and examining the influence of the Roche potential on the jet material. Specifically, we investigate whether including the companion star in the model significantly affects the dynamical evolution of the jet. Our findings reveal the formation of an arc-like structure in the density map of the jet cross section, which is attributed to the direct tidal effects. This implies that while the primary physical effects and characteristics of the outflow on a larger scale are attributed to the host accretion disk, the direct tidal effects on the jet dynamics have a substantial impact, particularly in the vicinity of the Roche lobe and toward the secondary star.


INTRODUCTION
Astrophysical jets are fast-moving streams of material that have been observed in various astrophysical systems, including young stellar objects (YSOs), microquasars, and active galactic nuclei (AGN).They transfer energy, angular momentum, and mass to larger scales, such as molecular clouds, and significantly impact environmental evolution and processes like star formation.Jets are a well-known phenomenon, and their propagation and interaction with the interstellar medium have been the subject of numerous studies for many years (Banerjee et al. 2007;Mignone et al. 2010;Perucho & Bosch-Ramon 2012;Cielo et al. 2014;Sheikhnezami & Abgharian-Afoushteh 2019).Furthermore, jets and outflows are observable phenomena that can be used to trace the physical properties of their host sources, such as the central object and underlying accretion disk.It is now well known that the large-scale magnetic field and its interaction with the host accretion disk play a major role in jet formation and acceleration (Blandford & Payne 1982;Pudritz & Norman 1983;Shibata & Uchida 1985;Ferreira 1997).As a result, many studies have applied numerical simulations to investigate jet launching and acceleration while considering the disk as a boundary (Ouyed & Pudritz 1997;Anderson et al. 2006;Vaidya et al. 2011;Porth et al. 2011;Porth 2013;Fendt & Yardimci 2022) or including the disk evolution (Casse & Keppens 2002;Zanni et al. 2007;Murphy et al. 2010;Sheikhnezami et Fendt & Sheikhnezami 2013;Stepanovs et al. 2014;Mattia & Fendt 2020).There are plenty of studies of astrophysical jets that have focused on jet formation in single stellar systems.However, there are significant observations of jet precession or curved ballistic motion that strongly suggest that the jet source is part of a binary or even a multiple system.Numerous studies investigating binary systems including jets have been performed utilizing observational data or simulation techniques (Fendt & Zinnecker 1998;Shepherd et al. 2000;Crocker et al. 2002;Mioduszewski et al. 2004;Agudo et al. 2007;Cunningham et al. 2009;Monceau-Baroux et al. 2014;Keppens et al. 2014;Paron et al. 2016;Beltrán et al. 2016;Millas et al. 2019;El Mellah et al. 2019;Mellah et al. 2019;McLoughlin et al. 2021;Britzen et al. 2021;Dominik et al. 2021;Paron et al. 2022;López-Miralles et al. 2022).
Recently, numerical simulations of the disk forming the jet have been developed in three dimensions, making it feasible to study the non-axisymmetric properties of the jet as part of a binary system (Sheikhnezami & Fendt 2015).Moreover, in our recent papers, we have studied the effects of including the companion star in the model setup (Sheikhnezami & Fendt 2018, 2022).One important result from our 3D numerical MHD simulations is the observation of spiral structures in both the disk and the jet.Sheikhnezami & Fendt (2022) have discussed and demonstrated that the spiral structure observed in the jet is injected from the disk into the outflow and then propagates along the jet.This suggests that the main structures seen in the jet are delivered from the underlying disk.However, at low altitudes where the disk wind is close to the disk, the tidal forces on the jet materials are expected to be significant.As a result, the flow of the jet across the Roche lobe can be directly influenced by these tidal forces, leaving an imprint on the jet dynamics.Consequently, the resulting structure inside the jet is not solely the outcome of the injection process.
To investigate the direct impacts of tidal forces on the evolution of the jet within a binary star system, our specific focus is on the dynamical evolution of the jet.In this context, we treat the disk as a boundary.Initially, we utilize an axisymmetric simulation to explore and determine an appropriate model setup that enables the effective launching of a fast and stable jet from the disk boundary.Subsequently, we advance our setup to "three dimensions", enabling us to conduct a comprehensive study of jet formation following the full 3D evolution of a jet launched in single star and binary star systems.This allows us to examine the effects of the secondary star on the dynamics of the jet within binary systems.In our 3D simulations, the disk acts as a surface, and we do not consider its evolution.
Our paper is organized as follows: In Section 2, we discuss our model setup.In section 3 , we present and discuss our 3D simulation of a jet launched in a single star system.In section 4, we present the 3D simulations of the jet formation in a binary system, and we compare the obtained results.In Section 5, we summarize the results.

MODEL SETUP
We investigate the launching of an MHD outflow from a sub-Keplerian disk, which serves as a boundary.Our primary focus is to examine the direct influence of a companion star on the structure and evolution of the outflow.This is achieved by incorporating the Roche potential and conducting 3D MHD simulation runs.
For all simulations presented in this study, we utilize version 4.4 of the MHD PLUTO code (Mignone et al. 2007(Mignone et al. , 2012)).This code allows us to solve the timedependent Ideal MHD equations, which include the conservation of mass, momentum, and energy, (2) (3) Here, ρ is the mass density, ⃗ v is the velocity, P is the thermal gas pressure, ⃗ B represents the magnetic field and Φ indicates the gravitational potential.To describe the evolution of the magnetic field the induction equation is applied, i.e., The gas pressure follows an equation of state P = (γ − 1)u with the polytropic index γ and the internal energy density u.The total energy density is, , the binary separation D, the orbital period of the binary T and the mass ratio of the secondary to the primary star is q = Ms Mp .a The first row is the axisymmetric run and the rests are 3D runs.

Run
In this study, the heating and cooling processes are not considered to avoid the complexity.

Units and normalization
In this section, we explain the normalization we use throughout the paper.The inner radius of the disk boundary is shown by r i , and we adopt r i = 1.All distances are expressed in units of inner disk radius r i .The gas pressure and density at this radius are P i and ρ i .Velocities are expressed in units of the Keplerian velocity at the inner radius of the disk v k,i which we adopt v k,i=1 .Time is estimated in units of t i = ri v k,i .The magnetic field is measured in units of the magnetic field at the inner disk radius B i which can be defined as, B i = 2Pi βp where β p denotes the plasma beta parameter of the materials 3 .

Initial Conditions
We define the initial conditions following the setup applied by Vaidya et al. (2011) to explore a magnetized jet emanating from a disk boundary prescribed in sub-Keplerian rotation.For our axisymmetric reference run, we define the potential due to the gravity of the central star, denoted by Φ, using the following equation: We implement the companion star by applying the Roche potential (see section 3.3).Here, (r g , z g ) are softening parameters included to avoid singularity at the origin and are taken as r g = z g = 0.21 in our setup.Initially, we prescribe a hydrostatic equilibrium with a density distribution, The thermal pressure follows the polytropic equation of states p = Kp γ where γ = 5/3 is the polytropic index.We set the initial radial and vertical velocities to be proportional to the sound speed at the radius of the inner disk boundary, denoted by C s = γp i /ρ i .These values are arbitrary.
The initial magnetic field is initially set to be purely poloidal and follows the prescription used in Vaidya et al. (2011).The vector potential, denoted as A Φ , is defined as follows: where the parameter z d represents the thickness of the disk.
Additionally, we have also investigated an alternative initial magnetic field distribution defined by Zanni et al. (2007).In this case, the magnetic flux function ψ is used and given by: where B z,i represents the vertical field strength at (r = r i , z = 0).The vector potential, denoted as A ϕ , is then given by: The latter definition is more consistent with our 3D computational grids which excludes the innermost part of the radial axis (r < 0.5).As a result, we do not encounter any cuts in our magnetic field distribution.

Numerical setup
In our axisymmetric model setup, we use a numerical grid of (600 × 600) cells applying cylindrical coordinates with domain size of (r × z) = (120 × 120)r i .The computational domain covers the rectangular grids (400 × 400) cells with uniform spacing in the radial and vertical direction for the domain size of (r × z) = (20 × 20)r i and applying extra stretched grids (200 × 200) for the rest of the domain.
The piece-wise parabolic method (Mignone 2014) is used for spatial integration with a 3 rd -order Runge-Kutta scheme for time evolution and HLL Riemann solver.For the magnetic field evolution, we apply the DIV-CLEANING a method that enforces the solenoidal Condition of the magnetic field ∇.B = 0 2.4.Boundary Conditions Figure 1 shows the numerical setup for our simulations.The left panel shows our axisymmetric computational domain and the corresponding boundary conditions.The right panel, on the other hand, provides a schematic representation of the 3D setup, which includes the boundary conditions and the orbital plane of the binary star.In this section, we consider the left panel which includes the four boundaries in cylindrical coordinates (r, z).At the upper vertical and outer radial directions, we use the standard PLUTO outflow boundary conditions (zero gradients), while on the rotation axis, we apply axisymmetry with ∂/∂ϕ = 0 (This condition is applied only for our axisymmetric simulation and not for our 3D setup (see section 4.1).

Equatorial Plane boundary
At the lower part of the vertical direction, we define the "disk boundary" that emulates the accretion disk and has an essential role in our model setup.At the disk boundary, the disk parameters are given as the inflow values (Ouyed & Pudritz 1997;Vaidya et al. 2011).
The disk boundary contains the inner gap of the area (r < r i ) and the disk region of (r > r i ).At the "disk boundary" the inflow parameters of the sub-keplerian disk are given which should be in equilibrium with the hydrostatic corona initially.The density component is prescribed by, and following the polytrophic law the pressure of the disk material is, Here χ is the Fermi function (Vaidya et al. 2011) and is defined as follows, Applying the Fermi function imposes a slightly subkeplerian velocity to fulfill the radial equilibrium at the disk boundary.The toroidal velocity at the disk boundary is given by, Here v Kep is the toroidal velocity at the inner disk boundary It should be noticed that we implement timedependent boundary conditions for the disk boundary, which means that Equation 15 is used in our setup to determine the poloidal velocity.Consequently, the velocity of the inflowing material at the disk boundary aligns with the poloidal magnetic field (v p ∥B p ) and depends on the evolution of the poloidal magnetic field.This approach differs from previous studies (Ouyed & Pudritz 1997;Vaidya et al. 2011) where fixed-in-time conditions are applied for the outflow inlet.
The Stationary axisymmetric MHD flows conserve some quantities (Pudritz et al. 2007) along the magnetic flux function which is defined by ψ = rA ϕ .For instance, the conservation of mass and the magnetic flux along a field line gives rise to the relation between the poloidal mass flux and the poloidal magnetic field, i.e., where η is called the mass load along the field lines and is one of the conserved quantities of the jet.
In addition, the toroidal field in the rotating flows like the jet follows the induction equation which is required to satisfy the following equation, Here Ω F is called the iso-rotation parameter.This parameter is interpreted as the angular velocity of the field lines and is another conserved value along the field lines in the jet.In our set up Ω F follows the Keplerian profile and we have Ω F (r) ∝ r −1.5 .To prescribe the consistent steady MHD jet, we impose the condition 16 to the disk boundary.
In the Equatorial boundary, we impose the ideal MHD approximation i.e., E ϕ = (⃗ v × ⃗ B) ϕ = 0 which implies that the poloidal velocity and poloidal magnetic field are parallel, v p ||B p .We impose this condition for the vertical magnetic field B z while we use the zero gradient for the radial magnetic field.
We have presented and discussed our (2.5D)axisymmetric simulation in the appendix to examine the 3D simulation of the jet launched in a single star system which can then be incorporated to non-axisymmetric 3D simulations of the jet in a binary star system.

3D MHD SIMULATION OF JET LAUNCHING
In principle, axisymmetry imposes a symmetry constraint that may result in the loss of certain threedimensional effects in our approach.Consequently, to investigate the formation of a magnetized jet in a binary system, it is necessary to enhance our axisymmetric model setup to incorporate three dimensions and accurately capture the full three-dimensional evolution of the magneto hydrodynamic (MHD) jet launching from a disk boundary.Therefore, in this section, we present and discuss our 3D reference run, which serves as a basis for our analysis.Subsequently, we implement the companion star into the model and examine the launching of the jet in a binary star system.

3D model setup
To enhance our study and extend the model setup to 3D, we utilize polar coordinates (r, ϕ, z), as depicted in Figure 1 (right panel).The computational domain in our 3D simulation spans from (0.5, 80) radially and (0, 200) vertically, measured in units of the inner disk radius r i .The azimuthal direction covers an angle range of (0, 2π).Our computational grid consists of 500 × 100 × 600 cells, with a uniform grid in the ϕ direction, as well as in the inner radial and vertical directions.Specifically, we incorporate 400 uniform cells for 0 < r < 20 and 400 uniform cells for 0 < z < 20, while applying a stretched grid for the remaining part of the grid.
It is important to note that, unlike in our axisymmetric simulation, the z-axis no longer serves as a symmetry axis.Instead, we implement an "inner boundary" (see Figure 1).At the inner boundary which is located at r = 0.5r i , we apply a zero gradient for all variables.However, we restrict the velocities to allow for"only" the inward motion of material into the jet domain.Specifically, we "only" allow for negative radial velocity (v r <= 0) and positive vertical velocity (i.e., in a situation where v r is indeed positive, it is set to 0).
Applying the zero gradient is a proper choice since the simulation is stable.We allow for the freely copying the variables into the ghost cells.In addition, the inner opening angle is really small and is about 1 degree in our setup.If the inner cone is sufficiently small, We would not expect a physical effect, as the flow can always flow besides the axis.Only material that actually hits the inner 1 degree or less in the cone, is lost from the grid.However, it does not affect the flow, since it is not injected again.Thus, the physical evolution of the jet is not affected.For other variables, we utilize a zero gradient condition, meaning that ∂/∂z = 0 and ∂/∂r = 0 at the z-axis is applied.At the outer radial and vertical boundaries, we apply a modified outflow boundary condition that only allows outward radial and vertical velocities.The lower boundary in the z direction represents the disk boundary, as explained in section 2.4.1.In the ϕ direction, we apply a periodic boundary condition.We utilize the same initial conditions and normalization as our axisymmetric model setup.
3.2.3D Simulations of Jet Formation in a Single Star System In this section, we first discuss our 3D reference run, which presents the formation of a magnetized jet from a disk boundary.The simulation run time is about 5000 dynamical times, which is long enough for the jet to reach steady states.Figure 2 displays the 3D views of the mass density distribution for run "3Dc0" at dynamical times t=0, 3400, 5000 t i .The solid lines represent the magnetic field lines.Figure 3 also demonstrates the 2D slices of the mass density maps at various planes at dynamical time of 5000t i .
Considering Figure 2 and 3, it can be observed that a smooth jet is continuously formed from the disk boundary and expands into a larger area.
Other physical variables, such as the velocity components, can also be studied to assess the quality of the 3D model setup.Figure 4 shows the different velocity components of run "3Dc0", including the radial velocity component (v r ), toroidal (v ϕ ), and vertical velocity components (v z ) of the formed jet, respectively.We observe that the jet reaches high velocities of about 20 times the local sound speed or about 2 times the inner radius Keplerian velocity, i.e., about 200km/s.The early evolution of the jet is quite fast and the initial shock wave reaches the upper boundary quickly.This can be estimated by considering the kinematic timescale of the jet, which is given by the ratio of the distance traveled (∆L) to the jet velocity (v jet ≃ 2v i ).In this case, we obtain t kin ≈ ∆L/v jet ≈ 200/2 = 100.
Considering all the physical variables of the jet in run "3Dc0" with the single star, the evolution of the jet material exhibits that the jet keeps the axial symmetry properly, although we do not use the axisymmetric boundary and this shows the quality of the 3D model setup.

3D simulations of jet formation in a binary system
In this section, we present additional simulation runs that are detailed in Table 1.These runs incorporate the Roche potential generated by the companion star.Table 1 provides the input parameters for our simulations.The key parameters of the binary system include the mass ratio, represented as q = M p /M s , which signifies the mass of the secondary star relative to that of the primary star.Additionally, since the origin of the coordinate system remains at the primary star, the binary separation is determined by the position of the secondary star which stays outside the computational domain.
Figure 1 (right panel) displays the 3D computational domain and the position of the companion star relative to the primary star.In our model, the primary star is situated at the origin of the coordinate system, while the companion star remains outside the domain at a distance of r s .To account for the gravitational effects of the binary system, we utilize an effective potential given by: where M p and M s denote the masses of the central and companion stars, respectively.Here D represents the position of the companion star and is defined by: ⃗ D = xD cos ωt + ŷD sin ωt cos δ + ẑD sin ωt sin δ. (18) The first two terms in the equation 17 account for the gravitational effects of the primary and companion stars, respectively.The third term is known as the indirect term and arises due to the non-inertial nature of the coordinate system (The origin of our coordinates system stays in the primary star).This term accounts for the acceleration of the origin of the coordinate system and is essential for accurately modeling the motion of the binary system (Sheikhnezami & Fendt 2018).The parameter δ shows the inclination angle of the binary orbit concerning the disk boundary.In all the simulations described in this paper, we have set the inclination angle to zero.This means that the binary orbit and the disk boundary are co-planar.It is important to note that the companion star is located outside the computational domain and is not directly included in the simulation.
We have performed the additional runs in three dimensions, referred to as '3Dc1', '3Dc2', and '3Dc3', with mass ratios of the secondary star to the primary being 0.1, 1, and 2, respectively.To reduce the orbital period of the secondary object and to enhance the tidal effects during the run time, we choose a separation of approximately 100 times the inner disk radius (r i ) between the secondary and primary objects.It is important to note that in all runs, the orbit of the secondary object remains coplanar with the primary object.Subsequently, we proceed to compare our other simulation runs with the 3D reference run.By doing so, we can examine the impact of the Roche potential on the physical properties of the jet and gain further insight into the dynamics of the jet present in a binary system.

JET FORMATION IN A BINARY STAR SYSTEM WITH A HIGH MASS RATIO
Here, we consider the evolution of the jet in run "3Dc3" including a binary star with a mass ratio of q = 2 for the secondary star compared to the primary star.In Figure 5, we can observe the evolution of the mass density of the jet material for this run.The figure displays 3D rendering and multiple snapshots of the mass density at different planes, namely the x − z plane and the crosssection of the jet at z=50, at different dynamical times.The "+" symbol determines the position of the inner Lagrange point L1 of the binary systems.
Regarding Figure 5, we find the direct influence of the secondary star on the mass distribution of the jet.It is evident that the mass density distribution does not maintain axial symmetry, and there is an increase in the amount of jet material.The increase of density is directed towards the inner Lagrange point and due to the orbital motion forms an "arc-like" structure. .Jet evolution in a binary system.Shown are 3D rendering and the snapshots of the mass density for run "3Dc3" with a mass ratio of the secondary to primary star q = 2 and the separation of D = 100r i .The snapshots are shown at different planes and at times 0, 1580, 3100, 5000 t i .Here the "+" symbol shows the projected position of the inner Lagrange point L1.
The arc-like structure is visible in the cross-section of the jet at z=50.This structure exhibits rotation which is synchronized with the orbital motion of the secondary star.To quantitatively observe the orbital motion of the feature caused by the secondary star, we have presented the angular profiles of the radial velocity v r and mass density ρ for the jet in run "3Dc3" including a binary star with mass ratio of 2. These profiles were obtained at a radius of r = 60 and z = 50, and at various time steps, as shown in Figure 6.
Analyzing the angular profile of the radial velocity, we observe different peaks occurring at different times, coinciding with the position of the inner Lagrange point ϕ L1 (the "+"symbol).This indicates that the jet material with higher radial velocity leaves the grid more rapidly resulting in the lower density at that region ("+" symbol in angular profile of density).In addition, in the crosssection of the jet and at larger angular distances from ϕ L1 (indicated by the "+" symbol), an increase in mass density is observed.The enhanced mass density, combined with the orbital motion of the companion, results in the formation of an arc-like structure.The angular profile of the radial velocity provides further confirmation that the affected regions within the jet by direct tidal effects, are rotating and synchronized with the orbital motion of the secondary.
It appears that the arc-like feature is predominantly observed at lower altitudes.As we move further out in the jet, the deviation from axial symmetry becomes less pronounced.This observation is supported by considering the jet at different cross-sections.To provide a clearer understanding and visual representation, we present the mass density of the jet at various cross-sections,i.e., z=5, 50, 100, 180, in Figure 7.The figure clearly shows that the enhanced density is formed towards the secondary and is visible at all altitudes.However, as we move far- ther away from the companion star and reach greater heights, the impact of this effect diminishes, although it still leaves a discernible imprint in all of the displayed cross-sections.Moreover, the arc-like structure exhibits a consistent azimuthal angle across various cross-sections of the jet, indicating that there is no noticeable delay in the formation of the arc-like structure along the jet.

Comparison of all the 3D simulations
For the sake of the comparison, we have performed different binary runs implementing a secondary with the mass ratios of q= 0.1, 1, 2. Figure 8 illustrates the snapshots of the mass density at a cross-section of z=50 for the jets formed in a single star system and a binary system with mass ratios of q = 0.1, 1, and 2 (from top to bottom) respectively at various dynamical times.The "+" symbol determines the position of the inner Lagrange point L1 of the binary systems.One can observe that the stronger tidal effects appear due to the larger mass ratio which can produce a denser region in the density map.It is seen that run "3Dc1" with a mass ratio of 0.1 represents a very small deviation from the axial symmetric distribution of the mass density (second row).In addition, as we increase the mass ratio the deviation becomes larger and a prominent arc-like feature is observed in the jet cross-section of z=50 in run "3Dc3".Additionally, we consider other physical variables, like the radial and toroidal velocity component of the jet, and check if there is any signature of the direct tidal effects is seen.For instance, Figure 9 illustrates the radial velocity of the jet material at the cross-section of z=50 for all the runs listed in Table 1.It is observed that the radial velocity disrupts the symmetry and increases towards the inner Lagrange point and the companion star.In Figure 6, we have demonstrated that the angular position of the peak in radial velocity aligns with the location of the inner Lagrange point, for jet in run '3Dc3".Now, Figure 9 provides further confirmation that the affected regions within the jet by direct tidal effects, are rotating and synchronized with the orbital motion of the secondary.
If we specifically examine the radial velocity in the x−z plane for the simulation that incorporates the binary star with a mass ratio of 2, identified as run "3Dc3" in Table 1, we can observe a prominent radial deviation of the jet (see Figure B1 in the appendix).This deviation indicates the bending of the jet towards the companion star, which changes over time due to the position of the inner Lagrange point and the motion of the secondary star.
The tidal effects on the jet material are primarily influenced by the extension of the Roche lobe, which, in turn, depends on the binary separation and the mass ratio between the secondary and primary star.In the vicinity of the Roche lobe, the direct tidal effects are significant, but their effectiveness may diminish on a larger scale.Our findings indicate that the companion star exerts significant direct influences, particularly in the region near the Roche lobe of the binary system.These influences are the bending and radial deviation of the jet, as observed in our binary simulations.

Outflow mass flux
The outflow mass flux is an important observable and it is essential to be measured properly.In principle, what comes in from the disk boundary (inflow boundary condition) must go out from the grids.So if we compare the inflow at the disk boundary with the total outflow mass flux (radial plus vertical mass flux), it must be the same.In some previous studies (Ouyed & Pudritz 1997;Fendt & Čemeljić 2002;Porth & Fendt 2010;Vaidya et al. 2011), the disk was taken as a time-independent boundary condition.However, Krasnopolsky et al. (1999) applied a more flexible and time-dependent conditions applied for the disk boundary.They have discussed that it is impossible to assign arbitrary fixed values to all the quantities at the disk boundary and some quantities needs to be left free to adjust according to the inflow inlet.As we mentioned before, in our approach, the disk boundary is defined as a time-dependent and is not fixed in time.In particular, we define the poloidal velocity according to the equation 15.Thus in our approach the v p ∥B p . and the poloidal velocity is dependent to the magnetic field evolution as well.Consequently, the total outflow mass flux evolves during the time and also is dependent on the field strength and distribution.
In order to check the conservation of the total mass flux (radial plus vertical mass flux) in our simulation, we apply different control volumes shown in Figure 10 where z m is the height of the control volume and r m is the outer radius of the jet.
We incorporate the control volume and estimate the integrated vertical mass flux of the outflow applying the following integral: where Ṁejec is the vertical mass flux, ρ is the mass density, v z is the vertical velocity and r is the cylindrical radius.The integral is taken over the entire cross-sectional area of the jet at z = z m and along the radial extension of r = 0 to r = 80r i .Subsequently, we measure the radial mass flux using the following integral: Here, r m indicates the outer radius of the jet and z m is the height of the control volume.
Figure 11 demonstrates the time evolution of the vertical Ṁz , radial Ṁr and the total mass flux Ṁtot of the jet in run "3Dc0" including the single star (top) and the jet in run "3Dc3" including the binary star of mass ratio q = 2.The profiles compare the measured mass fluxes for different control volumes of radius r m = 80 and z m = 100, 150, 200, respectively.
It is evident from the figure that the total mass flux measured for different control volume is similar for the jet in single star.In addition, for the jet in binary star system the total mass flux is obtained for different control volumes is similar.Thus the total mass flux is conserved in each runs.
In addition, it is seen that the total outflow mass flux is slightly larger for the jet in the binary star system compared to the one in single star system.This difference arises due to the evolution of the magnetic field strength and distribution, which deviate from the initial conditions in both runs over time.In each run the total mass flux is conserved and the time-dependent characteristic of the disk boundary determines the total outflow mass flux.
For example, Figure B2 illustrates the vertical magnetic field B z in the x − z plane for the jet in the "3Dc0" run with a single star (top) and the "3Dc3" run with a binary star system having a mass ratio of 2 (bottom) at different dynamical times.The figure demonstrates that the field strength is not symmetric and is larger in certain regions for the jet in run "3Dc3" with a mass ratio of 2. The stronger magnetic field in the binary system leads to a more efficient magnetic torque and thus a higher outflow mass flux for the jet in binary star.This accounts for the slightly larger total outflow mass flux observed in the binary run compared to the single star run.
We can also compare the radial and vertical outflow mas flux in all runs presented in Table 1.
For a better comparison, we estimate the vertical mass flux of the jet at two different heights, z = 100r i and z = 150r i .
Figure 12 (the first two panels) illustrates the time evolution of the integrated vertical mass flux of the jet at heights of z = 100r i and z = 150r i , encompassing all 3D runs listed in Table 1.It is evident from Figure 12 that the vertical mass flux in the reference run (including the single star) is lower than in the binary runs, confirming that the presence of a secondary star in the binary system affects the mass distribution within the outflow.As we include the binary star and increase the binary mass ratio, the vertical mass flux exhibits an increasing trend, particularly at later dynamical times.This indicates a deviation from axial symmetry in the jet structure, leading to the formation of a dense region inside the jet.Moreover, if we measure the outflow mass flux further out, we observe the same trend, although the overall mass flux is lower than that at z=50 in all runs.This is because the jet material flows out in a poloidal direction (along ⃗ v p = ⃗ v r + ⃗ v z ) and it is required to consider the radial mass flux of the outflow, as well.
Figure 12 (last panel) presents the time evolution of the integrated radial mass flux of the jet at a radius of r = 80r i across the cylindrical surface for all 3D runs listed in Table 1.We find that the radial mass flux in the reference run is once again lower than in the binary runs, confirming that the mass density and radial velocity change towards the companion star and increase.This radial mass flux is integrated over the cylindrical surface of the jet, indicating the radial deviation of the mass density towards the secondary star.This is another indication on the jet bending caused by the tidal effect on the jet structure.
The obtained trend for the outflow mass fluxes confirms, once again, that while the evolution of the accretion disk plays a significant role in forming the observed features within the jet, the Roche potential also exerts a direct impact on the jet and alters its mass flux.4.3.Comparison to our previous study In our previous study (Sheikhnezami & Fendt 2022), we conducted a comprehensive analysis by simultaneously considering the evolution of both the disk and the jet.It has been shown that the evolution of the disk in binary systems results in the formation of spiral arms within the disk, which is subsequently propagated into the jet.These spiral arms exhibit a small delay in forming along the jet due to the combined influence of the orbital motion of the binary star and the propagation of the jet.In addition it has been discussed that on a larger scale, the prominent features observed in the jet structure are attributed to the residual effects from the evolution of the host accretion disk.
In comparison, in the present study, we did not observe the formation of spiral structures within the jet.Furthermore, the region affected by the orbital motion of the binary system displayed a consistent azimuthal angle across different cross-sections of the jet.These findings indicate that the arc-like structure observed in the density or the enhancement of radial velocity is directly induced by the secondary star.
Our current results are entirely consistent with our earlier findings in Sheikhnezami & Fendt (2022), which demonstrated that the primary physical effects and features of the outflow at larger scales, such as the formation of spiral arms, originate from the host accretion disk evolution.Additionally, our current study reveals that the direct tidal effects on the dynamics of the jet make a substantial contribution, particularly in the region near the Roche lobe and towards the secondary star.

CONCLUSIONS
We have studied the formation of a magnetized jet within a single star and binary star systems.Our specific focus relies on the jet evolution and thus, we treated the disk as a surface and disk evolution is not considered.As a first step, we utilized an axisymmetric run to explore and determine an appropriate model setup that establishes the effective launching of a fast and stable jet from the disk boundary.Subsequently, we have advanced our setup to "three dimensions", enabling us to conduct a comprehensive study of jet formation following the full 3D evolution of a jet launched in single star and binary star systems.Thus, we have implemented the companion star (outside the computational domain) and performed various 3D simulations with the mass ratios of q = 0.1, 1, 2 to investigate the direct tidal effects on the jet evolution.
Our study yielded the following results; • We presented the axisymmetric MHD simulations of jet formation from a disk boundary (shown in the appendix).Through our analysis, we find that the global evolution of the jet reaches a steady state.Subsequently, we have enhanced our study and developed our axisymmetric model setup to three dimensions.Considering our 3D reference run "3Dc0", We find that the jet from the disk surface is continuously launched and reaches a high speed of about 20 times the local sound speed which is about the observed values of a typical jet from young stellar objects.
• Furthermore, we find that in our 3D reference simulation, although we do not employ an axisymmetric boundary for the z-axis, the formed jet exhibits a reasonable degree of axial symmetry, particularly after the early evolution (t ≃ 300).This observation indicates the quality of the 3D model setup which can be incorporated for studying the full 3D evolution of the jets launched in a binary system.
• By implementing the secondary star in the model setup, we find distinct variations in the physical parameters of the jet.For instance, we see the significant difference in the mass density distribution or the radial velocity v r components of the jet in binary star compared to the single star system, Specifically, inside the jet and toward the secondary star an increasing in the mass density is formed and represent an "arc-like" structure.This feature is rotating and is synchronized with the orbital motion of the secondary.Our findings suggest that the companion star has a substantial and direct influ-ence, especially in the region near the Roche lobe of the binary system.These influences manifest as the formation of a denser feature within the jet, as well as the bending and radial deviation of the jet, as observed in our binary simulations.
• Based on our analysis of the total mass flux for the jets in the single star and binary star systems, we can conclude that the total mass flux is conserved within each run when considering different control volumes.However, the time-dependent nature of the disk boundary plays a crucial role in determining the total outflow mass flux.Furthermore, we have observed that the total outflow mass flux is slightly greater in the jet of the binary star system compared to the single star system.This difference can be attributed to the evolution of the magnetic field strength and distribution, which deviate from their initial conditions as the simulations progress over time.
• Moreover, based on evaluating the vertical mass flux in all runs, we find that the integrated value in the reference run (single star) is lower than in the binary runs, and this trend increases as the mass ratio increases.It indicates that the presence of a secondary star, influences the mass distribution and the velocity field, especially at the outer part of the outflow, leading to a larger deviation from the vertical motion of the materials and a redistribution of mass within the jet.
• Our results are entirely consistent with our earlier findings in Sheikhnezami & Fendt (2022), which demonstrated that the primary physical effects and features of the outflow at larger scales, such as the formation of spiral arms, originate from the host accretion disk evolution.Additionally, our current study reveals that the direct tidal effects on the dynamics of the jet make a substantial contribution, particularly in the region near the Roche lobe and towards the secondary star.
We thank Andrea Mignone and the PLUTO team for the possibility of using their code.We would also like to acknowledge the helpful and insightful comments from an anonymous referee, which have greatly contributed to productive discussions and an improved presentation of our findings.We would also thank Christian Fendt for the helpful comments.Our simulations were performed on the TURIN cluster of Institute for Research in Fundamental Sciences (IPM) and the COBRA cluster of the Max Planck Society.

GENERAL EVOLUTION OF AXISYMMETRIC JET LAUNCHING
Here, we present and discuss the results of our axisymmetric reference run "2Dc0" which investigates axisymmetric jet launching.Figure A1 shows the time evolution of the mass density, with solid lines indicating the poloidal magnetic field lines.The continuous and powerful jet is seen to form from the lower part of the grid (at the disk boundary) and develop outward into a larger area.As the jet spreads out, it forms a bow shock that affects the corona and sweeps up materials.The resulting jet exhibits typical collimation, which is evident and caused by the toroidal component of the magnetic field.Our model setup and the resulting jet confirm the validity of the "disk-wind" theory originally proposed by Blandford & Payne (1982) and later supported by Pudritz & Norman (1983).
Figure A2  that the jet materials initiate with low velocities and achieve high speeds as time passes.The bow shock is visible in the velocity map and eventually leaves the grid.In addition, Figure A3 shows the snapshots of the ratio of the poloidal velocity v p to the Alfven v A and fast magneto-sonic v fm speed for the reference run in logarithmic scale and at dynamical time 1000, respectively.Referring to Figure A3, it is apparent that the outflow materials exhibit highly super-Alfvenic behavior and reach super-fast magneto-sonic speeds.
The global evolution of the axisymmetric reference run "2Dc0" demonstrates that the jet achieves a steady state, indicating that it can be employed as a basis for extending the model to 3D, in the next step of our study.

THE PHYSICAL VARIABLES OF A JET IN A BINARY SYSTEM WITH A MASS RATIO OF 2
Here, we present some other variables of the jet in run "3Dc3" with the mass ratio of 2.
Figure 1.Computational domain.The left panel shows our axisymmetric computational domain and the corresponding boundary conditions.The right panel, on the other hand, provides a schematic representation of the 3D setup, which includes the boundary conditions and the orbital plane of the binary star.

Figure 2 .
Figure 2. 3D reference run evolution.Shown are the 3D view of the mass density maps of run "3Dc0" listed in Table1 at dynamical times t = 0, 3400, 5000t i .The overlaid lines represent the magnetic field lines.

Figure 3 .Figure 4 .
Figure3.3D reference run evolution 2. The displayed images represent 2D slices of the mass density maps at various planes for the "3Dc0" run, taken at dynamical times of 5000t i .
Figure5.Jet evolution in a binary system.Shown are 3D rendering and the snapshots of the mass density for run "3Dc3" with a mass ratio of the secondary to primary star q = 2 and the separation of D = 100r i .The snapshots are shown at different planes and at times 0, 1580, 3100, 5000 t i .Here the "+" symbol shows the projected position of the inner Lagrange point L1.

Figure 6 .Figure 7 .
Figure 6.Angular profile.Shown are the angular profiles of the radial velocity vr and mass density ρ for the jet in run '3Dc3' with the mass ratio of 2, at radius r = 60 and z = 50 and at times t = 316, 1580, 3193, 5000 t 0 .Here the "+" symbol shows the angular position of the inner Lagrange point (ϕ L 1 ).

Figure 8 .
Figure8.Comparison of jets in binary and single star systems.The snapshots depict the mass density at a cross-section of z=50 for the jets formed in a single star system and a binary system with mass ratios of 0.1, 1, and 2 (from top to bottom) respectively at various dynamical times.The "+" symbol determines the position of the inner Lagrange point L1 of the binary systems.

Figure 9 .
Figure9.Radial velocity comparison.The provided snapshots display the radial velocity at the cross-section of z=50 for the jets in the runs labeled as '3Dc0','3Dc1, '3Dc2', and '3Dc3', respectively, at various dynamical times.Here the "+" symbol shows the position of the inner Lagrange point.

Figure 10 .
Figure 10.Control volume.Shown is the control volume incorporated to measure the radial and vertical mass flux of the jet.Here zm is chosen to be 100,150,200 and r is 80 (the outer radius of the jet).

Figure 11 .
Figure 11.Mass flux conservation.Shown are the time evolution of the vertical Ṁz, radial Ṁr and the total mass flux Ṁtot for control volume at rm=80 and zm=100,150, 200.The top panels correspond to the jet in single star system and the bottom to the jet in binary star system, with mass ratio of 2.

Figure 12 .
Figure 12.Comparison of mass fluxes.Shown are the time evolution of the integrated vertical mass flux of the jet at a height of z = 100r i (first) and z = 150r i (middle) across the entire cross-section of the jet.The last panel demonstrates the time evolution of the integrated radial mass flux of the jet at a radius of r = 80r i across the cylindrical surface, for all 3D runs presented in the paper.

Figure A1 .Figure A2 .
Figure A1.Reference run.Shown are the maps of the mass density in the Logarithmic scale for the dynamical times t=0, 50, 100,1000.The solid lines indicate the poloidal magnetic field lines.

Figure A3 .
Figure A3.Reference run.Shown are the snapshots of the ratio of the poloidal velocity to the Alfven and fast magneto-sonic speed for the reference run in logarithmic scale and at dynamical time 1000, respectively.

Figure B1 .
Figure B1.Radial velocity evolution.The provided snapshots display the radial velocity at plane x − z for jet in run "3Dc3" with the mass ratio of 2, at various dynamical times.

Figure B2 .
Figure B2.Vertical magnetic field.The presented snapshots depict the vertical magnetic field Bz in the x − z plane for the jet in the "3Dc0" run with a single star (top) and the "3Dc3" run with a binary star system having a mass ratio of 2 (bottom) at different dynamical times.

Table 1
Characteristic simulation parameters: The input parameters for our simulation runs are defined in this table.The parameters shown in the table are plasma beta βp = 8πP B 2 , the mass load coefficient η = ρV P