Magnetosphere of an orbiting neutron star

We conduct force-free simulations of a single neutron star undergoing orbital motion in flat spacetime, mimicking the trajectory of the star about the center of mass on a compact binary system. Our attention is focused on the kinetic energy being extracted from the orbit by the acceleration of the magnetic dipole moment of the neutron star and, particularly, on how this energy gets distributed within its surrounding magnetosphere. A detailed study of the resulting magnetospheric configurations in our setting is presented, incorporating as well the effects due to neutron star spin and the misalignment of the magnetic and orbital axes. We find many features resembling those of pulsar magnetospheres for the orbiting neutron star—even in the absence of spin-of particular interest being the development of a spiral current sheet that extends beyond the light cylinder. Then, we use recent advances in pulsar theory to estimate electromagnetic emissions produced at the reconnection regions of such current sheets.


I. INTRODUCTION
A new era of multimessenger astronomy has started with the detection of gravitational waves (GW) from a binary neutron star merger (GW170817) by advanced LIGO/advanced Virgo [1,2], followed by broadband electromagnetic (EM) observations [3]. Binary systems involving a neutron star (NS) are the most likely sources for such simultaneous detections of GW and EM signals. In this context, EM emissions from the relatively cleaner environment preceding the merger could provide crucial information about the merger process, sky localization of the source, and the physical parameters of the system, which cannot be accurately obtained only by the gravitational wave observation.
The NS in the binary is generally expected to possess an approximately dipolar magnetic field, and they are likely to be surrounded by a force-free magnetosphere. The source for the precursor EM counterparts comes fundamentally from the orbital and rotational energy of the binary and its individual constituents. This kinetic energy is first electromagnetically extracted from each compact object, by means of the surrounding plasma, and later reprocessed within the magnetosphere to produce the observable emissions. Starting from the pioneering work of Goldreich and Julian [4] and Blandford and Znajek [5], the mechanisms to tap rotational energy from a compact object (immersed in a tenuous plasma) were unveiled: while a NS admits a classical electrodynamics interpretation as a Faraday disk, for a spinning black hole (BH), the energy is instead extracted in a form of generalized Penrose process, widely known as the Blandford-Znajek mechanism. For the orbital motion, on the other hand, there are two such extraction mechanisms identified. The first one is associated with a purely classical effect, produced by the acceleration of a magnetic (e.g., dipole) moment [6]. The EM energy loss rate from the motion of a single NS with a dipole magnetic field, as well as for two NSs in a binary, has been estimated (assuming vacuum) in Ref. [7]. The second mechanism is known as unipolar induction, which essentially captures the transfer of linear momentum into EM energy, due to the motion of a conductor (or compact object) through a magnetized medium. This effect was studied a long time ago for moving conductors such as satellites (see, e.g., Refs. [8,9]), later extended to NS binaries [10][11][12][13], and also generalized for BHs, relying on ideas from the membrane paradigm [14] to build a circuit model in black hole-neutron star (BHNS) binaries [11,12,15,16].
For a BHNS binary, one might associate the moving magnetic dipole (MD) effect with orbital energy being extracted from the NS, while the unipolar induction (UI) mechanism would operate to remove energy from the BH as it moves across the magnetic field of the NS. In a neutron star-neutron star (NSNS) system, on the other hand, both mechanisms are expected to operate together at each NS, with intensities depending mainly on their relative magnetizations. That is, each star could produce the MD effect due to the orbital motion of its own magnetic moment, and could operate as UI as it moves across the magnetic field of its companion. Of course, this is only schematic, and nontrivial superposition of these two mechanisms, along with other possible effects (like, e.g., magnetospheric flares [17]), would complicate the picture. The relative orientations of the magnetic moments (with respect to the orbital plane, and among them in the case of NSNS) and the role played by the spin of each compact object complicates things even further.
Simulations of compact NS binaries in full general relativity, paying attention to precursor EM signals, have been carried out (see, e.g., Ref. [18] for a review). Late-time inspiral phases of NSNS were considered in Refs. [19][20][21] using general relativistic force-free (GRFF) simulations, broadly matching ideal magnetohydrodynamics stellar interiors with an exterior force-free magnetosphere. Also, BHNS binary systems at fixed orbital separation-with and without BH spin-were studied in Ref. [22], relying on a similar GRFF numerical approach. Overall, the EM luminosity found was consistent with estimations from the UI model, 1 and the Poynting flux distributions were analyzed in both scenarios. Even though the MD mechanism has been mentioned in these articles as contributing to the total computed luminosity, typically less attention has been devoted to this effect in the literature as compared to the UI mechanism. For NSNS binaries, further magnetospheric properties such as the formation of current sheets (CSs) were reported [19][20][21], considering different relative magnetic strengths and orientations. However, in these previous studies, the details of the magnetosphere and the EM signals associated with its structure have not been considered in great depth.
One of the most challenging aspects of the problem is how-and how much of-the available energy transferred to the plasma produces the emissions on the different EM bands. Even for pulsars, for which the attempts to understand the principal magnetospheric properties has been done for many decades (e.g., Refs. [23-26]), an ultimate answer to this question still remains elusive, although, of course, a lot can be learned from the development of pulsar theory in this respect (e.g., Refs. [27,28]). In recent years, there has been significant progress on particle-in-cell (PIC) simulations, which self-consistently model the regions of plasma production and particle acceleration (see, e.g., Refs. [29][30][31]).
In this paper, we aim to further clarify the properties of the magnetosphere around a NS in a compact binary system. We consider a single NS with dipolar magnetic fields and surrounded by a force-free plasma, in orbital motion. Our NS follows a trajectory, in flat spacetime, which mimics that of a NS about the center-of-mass (CoM) on a particular binary system. This way, we pay our attention to the MD energy extraction process and how this EM energy gets distributed within the magnetosphere, in a sense, decoupling it from the UI mechanism (or other curvature effects) involved in the binary. Such a simplified setting allows us to conduct rather inexpensive, very accurate, numerical simulations for a detailed study of these systems. In particular, we investigate their magnetospheric features in close analogy to those of pulsars. And then, we use recent results of pulsar theory to infer possible EM signals from our numerical results. We see the expected γ and x rays luminosities estimated from the orbiting NS are rather weak, rendering their possible detections by current and near future facilities quite unlikely. And thus, as observed in Ref. [32], the best chance to detect EM precursor signals from the inspiral phase of compact binary systems is from magnetospheric pulsarlike configurations producing coherent radio emission.
The code used here to evolve the equations of force-free electrodynamics was first described in Ref. [33] and later extended in Ref. [34], where a careful treatment to handle the boundary conditions on the NS surface was presented, in contrast to the matching procedure used in previous GRFF simulations. 2 Since then, our code has been further tested and employed in other astrophysical scenarios [35,36], as well.
The paper is organized as follows. In Sec. II, we set up the problem and describe our numerical implementation. The results are presented in Sec. III, first focusing on the magnetospheric properties of a NS in circular orbits and then following inspiral trajectories associated with BHNS and NSNS binaries. Then, possible observational implications of our results are discussed. We summarize and conclude in Sec. IV. Throughout this paper, G and c denote the gravitational constant and the speed of light, respectively.

A. General setting
The purpose of this paper is to clarify the magnetosphere around a NS in a compact binary system. Here, we suppose that only the NS has a strong magnetic field and the magnetic field of its companion is much weaker. This is trivial for the BHNS case and would be a good approximation for the NSNS case because the first born NS in NSNS systems is likely to be weakly magnetized [37]. 1 This holds only prior to the last few orbits before the merger in Refs. [19][20][21], in which the dynamics becomes more violent and nonlinear. 2 It is not clear whether such matching employed in the GRFF simulations can accurately represent the (approximately) perfectly conducting NS surface, which is the key condition for describing the magnetosphere around the NS.
We shall also assume that the spacetime is flat and that the NS moves on a given trajectory centered around the origin of a Cartesian coordinate system x a ¼ ft; x; y; zg, which would represent the CoM of the binary system. The line element in this coordinates reads thus, with α ¼ 1, β i ¼ 0, and γ ij ¼ δ ij representing the lapse, shift, and spatial metric, respectively. The trajectory is defined here by the radial distance R o ðtÞ and the phase φ o ðtÞ [with Ω o ðtÞ ≡ _ φ o ðtÞ being the associated angular velocity]. However, our numerical domain will be centered on the NS instead, thus describing the dynamics from an adapted foliation with coordinatesx a ¼ ft;x;ŷ;ẑg (see Fig. 1). The coordinates transformation into this "comoving" foliation is defined by t ¼t; z ¼ẑ; x ¼x þ R o ðtÞ cosðφ o ðtÞÞ; and thus Therefore, the line element (1) in the new coordinates is where againα ¼ 1 andγ ij ¼ δ ij . However, there is now a nonvanishing shift vector accounting for the difference among the NS worldline tangentt a and the hypersurface normaln a in this foliation. That is, where we have dropped all time dependencies and defined v o ≔ R o Ω o to simplify the notation. We will consider essentially two types of trajectories in this work: (i) a purely circular motion at a constant radius R o and orbital frequency Ω o and (ii) quasicircular inspiral orbits which mimic those of a NS in a binary system close to merger. The boundary condition at the stellar surface is derived by assuming the perfectly conducting condition, which can be easily generalized to incorporate the NS spin at frequency Ω Ã by whereφ a ≡ ð∂φÞ a . Thus, the resulting condition on the electric field measured by a fiducial observer in this adapted foliation (i.e.,Ê a ≔ F abn a ), can be written The initial configuration is taken to be a magnetic dipolar field (of dipole-moment μ) and vanishing electric field. The NS is gradually set in motion (at fixed R o;0 ) until it reaches the desired initial orbital frequency Ω o;0 after some time t ¼ t 0 ; from then on, the system follows its prescribed trajectory. We focus primarily on the cases where the magnetic axis is aligned with the orbital angular momentum, but other scenarios in which these two axis are not aligned are considered as well. In our setup, the z axis is always perpendicular to the orbital plane, so the misalignment is attained by just tilting the magnetic moment by an angle χ along the x-z plane.

B. Numerical implementation
We evolve a particular version of force-free electrodynamics derived in Ref. [38], which has some improved properties in terms of well-posed-ness and involves the full force-free current density. More concretely, we shall consider the evolution system given by Eqs. (8)-(10) in Ref. [34]. Our numerical scheme to solve these equations is based on the multiblock approach [39][40][41][42], in which the numerical domain is built from several nonoverlapping grids where only grid points at their boundaries are sheared. The equations are discretized at each individual subdomain by using difference operators constructed to satisfy summation by parts. In particular, we employ difference operators which are eighth-order accurate on the interior and fourthorder accurate at the boundaries. Numerical dissipation is incorporated through the use of adapted Kreiss-Oliger operators. These compatible difference and dissipation operators were both taken from Ref. [43]. A fourth-order Runge-Kutta method is used for time integration.
We solve the force-free equations in a region between an interior sphere at radiusr ¼ R Ã that represents the NS surface (i.e., R Ã denotes the NS radius) and an exterior spherical surface located atr ∼ 75R Ã . The domain is covered by a total of 6 × 12 subdomains, with 6 patches to cover for the angular directions and 12 being the number of spherical shells expanding in radius. These spherical shells do not cover regions of identical radial extension, having more resolution near the inner boundary than in the asymptotic region; from layer to layer, the radial resolution is decreased by a factor 1.3. Typically, we adopt a resolution with total grid numbers of Nˆθ × Nφ × Nr with Nφ ¼ 2Nˆθ ¼ 240, while Nr is taken so as to satisfy Δr ≲ 0.7rΔθ everywhere in the domain. Here, Δr and Δθ denote the grid spacing forr andθ, respectively.
As already mentioned, the stellar surface is assumed to behave as an idealized perfect conductor. Thus, the normal component of the magnetic field is set to its dipole value assumed from the stellar interior, and the electric field is prescribed according to Eq. (6) in order to represent its orbital motion and spin. The electric field components are imposed, by means of the penalty method [40][41][42], fixing the incoming physical modes to a particular combination of outgoing modes. At the outer boundary, on the other hand, we set maximally dissipative (no-incoming) conditions to allow all perturbations to propagate away. The numerical implementation of such boundary conditions has been detailed in Ref. [34] (in particular, Sec. II C and the Appendix), so we recommend the interested readers to refer there for further details. In order to handle CSs, for which the force-free approximation breaks down, we use a standard approach in which the electric field is effectively dissipated to maintain the condition that the plasma is magnetically dominated (i.e., B 2 − E 2 > 0), as discussed in Ref. [33] (see also Ref. [44]).

C. Analysis quantities
We would like to monitor the EM energy and its associated fluxes. In force-free electrodynamics, the fourmomentum, p a ¼ −T ab EM t b , is conserved (i.e., ∇ a p a ¼ 0) 3 in the stationary spacetime. In the comoving coordinates fx a g, it reads Hence, we measure where the Poynting luminosity L is integrated on spherical surfaces of radiusr around the NS and withŜ i ≔ ϵ ijkÊ jBk being the spatial Poynting vector. We are also interested in monitoring the charge distribution and electric currents present during the dynamics. Thus, we shall look at the force-free current density along the magnetic field, as seen by a fiducial observern a , withρ ¼ D jÊ j being the charge seen by this observer. Actually, it is worth mentioning at this point that since in this case n a andn a represent exactly the same vector field, the splitting of the electromagnetic tensor in its magnetic and electric components is the same. Hence, the switch from one description to the other in terms of the electric and magnetic fields is quite direct. Since the vector transformation for the spatial index is trivial, one is only left with the appropriate displacement of the point where the field is evaluated. Therefore, although we employ the NS frame fx a g to evolve the fields, we use the "CoM" coordinates fx a g to plot all the relevant quantities and describe our results.
We also consider the Lorentz invariant quantity, where the sign is chosen plus for timelike and minus for spacelike currents, as in Ref. [27]. Counterstreaming of different signs of charge is required at regions having negative values of this quantity, which may lead to plasma instabilities and dissipation (see, e.g., Refs. [27,45,46]).

III. RESULTS
Our interest in this work is centered on the last few orbits of a NS on a compact binary system, until the orbit reaches an innermost stable circular orbit or the NS gets tidally disrupted. Thus, we start by a detailed study of the magnetospheric properties of the NS in a circular orbit, choosing parameters in a relevant range: i.e., We first focus on the case in which the magnetic moment is aligned to the orbital angular momentum (i.e., χ ¼ 0) and the NS is not spinning (i.e., Ω Ã ¼ 0). Our simulations always relax to stationary states presenting similar features to those of pulsar magnetospheres, with strong equatorial CSs. We analyze the electric charge/current distributions of the surrounding plasma, along with the resulting Poynting flux luminosity inside the orbital light cylinder. A comparison with vacuum magnetospheres, within the same setting, is also included here. Then, both the NS spin and the misalignment of the magnetic axis are incorporated into the picture. We vary these parameters (i.e., Ω Ã and χ) independently, and observe their impact on the solutions. Later, we consider representative quasicircular inspiral orbits, taking relevant values for the parameters in the contexts of BHNS and NSNS binaries, and connecting with our previous luminosity estimations. Finally, we elaborate on the implication of our results to EM observations.

A. Circular orbits with no spin
First, we consider the case in which the NS is in a circular orbit with no spin. In numerical computation, the NS is set into circular orbital motion by gradually bringing the angular frequency to its final value Ω o . The motion generates disturbances on the EM field that propagates into the surrounding plasma; when these waves return to the NS surface, they are reflected due to the perfectly conducting boundary condition. After an initial transient of about two orbits, these perturbations--continuously injected from (and reflected at) the stellar surface-equilibrate within the magnetosphere and lead to a quasistationary solution. Such configuration is illustrated by the three-dimensional plots in Fig. 2, where several representative magnetic field lines are shown. The injected Alfven waves twist the magnetic field lines, as the nonspinning NS follows the circular orbit, producing a pattern of alternate signs on their toroidal components, like the one depicted in the right panel of Fig. 2. Waves launched from magnetic footprints in opposite hemispheres of the star meet at the dipole equator, inducing sharp discontinuities on the toroidal field across the orbital plane. These discontinuities form an spiral CS, which is represented by the contour plots in the figures. Electric field is being effectively dissipated at these regions, in order to locally maintain a state in which B 2 − E 2 ≳ 0. This is a rather standard strategy employed in forcefree electrodynamics simulations to avoid violations of the magnetic dominance condition B 2 −E 2 >0 (see Ref. [44] for a physical justification). Thus, a small value of this Lorentz invariant quantity is a convenient indicator to illustrate CSs in this context (right panel). On the other hand, left panel of Fig. 2 shows that intense charge density (exceeding the Goldreich-Julian value) develops at the CS, possessing two main components of opposite sign over the representative magnetic field lines (with the colors indicating its magnitude in logarithmic scale), along with two contours of the electric charge density, normalized by the Goldreich-Julian value Ω o B=2πc. Right panel: normalized toroidal component of the magnetic field, B ϕ =B, at the y ¼ 0 plane; together with a contour of B 2 − E 2 ∼ 0.04B 2 (in gray), signaling the presence of a spiral CS over the orbital plane. Same representative magnetic field lines are also indicated in the figure (in solid black). spiral arms. As can be noticed, magnetic reconnection takes place near the transition from negative to positive charge, producing and ejecting closed magnetic loops outward.
The spiral structure of the magnetic field, clearly manifest on the equatorial CS, is also found for other relevant quantities such as electric charge/currents and Poynting flux density distributions; thus, reflecting the helical symmetry of the problem in this particular setting. The EM solution looks essentially static from a corotating frame.
As we shall see in more detail later, these configurations resemble in many aspects those of pulsar magnetospheres, with CS outside the wave zone ∼c=Ω o , and energy being extracted (in this case, kinetic energy from the orbital motion) and carried away by the surrounding plasma. Figure 3 shows the topology of the magnetic field lines after 2.5 orbits on the y ¼ 0 plane (i.e., the one containing both the NS and CoM). Vacuum and force-free solutions are included for a comparison. We remark that, even though the driving mechanism is the same (i.e., the effect of a moving MD moment) and both systems produce comparable radial Poynting flux intensities (shown in the color scale in the figure), the presence of the plasma changes the topology of the magnetic field and the form in which the EM energy is distributed and transported. In particular, the force-free plasma is an essential ingredient for CS formation.
We further analyze the Poynting flux distribution of the force-free solutions over spherical surfaces located in the outer wave zone. One such EM flux density, after 2.5 orbits at a sphere of radiusr ∼ 60R Ã , is shown in Fig. 4. As can be seen from the plot, the flux is concentrated along a broad beam of approximately 60°-75°in the azimuthal direction and within approximately 60°from the orbital plane. There is a visible feature in the distribution atθ ∼ 90°, reflecting the magnetic reconnections occurring at the equatorial CS. Once the magnetosphere has settled, these structures in the Poynting flux become stationary, just corotating with the orbit. And thus, the flow of EM energy produces a lighthouse effect at orbital frequency. Besides the specific details of the distribution presented here, our results are in good qualitative agreement with the ones obtained from GRFF simulations of binary systems involving a nonspinning BH companion [22], and also for a weakly magnetized NS companion in Refs. [19,20].
We integrate the luminosity for the late time solutions at different radius and normalize them with the EM luminosity by the MD radiation formula (see, e.g., Ref. [7]), The radial Poynting flux density after 2.5 orbits is plotted on a spherical surface at radiusr ∼ 60R Ã (i.e., in the wave zone). The flux is normalized with its peak value.
A typical radial distribution of the luminosity can be seen in Fig. 5, where we plot integrated Poynting flux through concentric spheres around the NS, as a function of radius. As in pulsars, the luminosity is constant up to the light cylinder R LC ≡ c=Ω o , where dissipation at the CS begins. Different numerical resolutions were considered in order to test its convergence. Throughout this work, we have employed the intermediate resolution N θ ¼ 120 for all the simulations, which is practically converging (it differs in less than 4% with respect to the higher resolution one).
We note that such a resolution required in the present orbital setting is larger than the one we typically needed for pulsars (see Fig. 2 in Ref. [34]), where N θ ¼ 80 was enough to resolve even for the misaligned configurations. Notice that dissipation taking place in the regionr ∼ ð1-4ÞR LC represents here 65% of the luminosity. These are the typical percentages that we get for most of the circular orbits explored, while the values obtained for pulsars were instead closer to 40% (within the same region). 4 We have also tested our outer boundary conditions, by considering different radial locations for the outer surface of the computational domain. The solutions found at the overlapping regions are essentially identical, both qualitatively and quantitatively. The scaling of the total luminosity as a function of the orbital velocity is analyzed and summarized in Fig. 6. For the vacuum magnetosphere, we find that L=L 0 ≈ γ 11 is the Lorentz factor of the orbital velocity) approximates the numerical results quite well (perhaps "accidentally"). On the other hand, the following expression is used as a guideline for the luminosity in the force-free simulations: The error bars in the force-free data indicate certain dispersion found in the values of L within the light cylinder. These errors become more significant at lower values of v o -as can be expected due to the enlargement of the wave zone-for which a higher resolution (than the one employed here for all the runs) would be needed. We note that the luminosity does not depend on the other dimensionless parameter of the problem, i.e., R o =R Ã . Within our numerical error, we could not find significant deviations on the ratio L=L 0 , when taking different values of this parameter at fixed v o .
The magnetospheric structure after 2.5 orbital periods is displayed in Fig. 7, in which relevant aspects of the solution are represented on the y ¼ 0 plane (the plane containing both the NS and the CoM). There, electric charge and current distributions are displayed (left and middle panels), as is the Lorentz invariant quantity ϱ from Eq. (13). Magnetic field lines projected onto the plane are also shown in all these plots for a reference.
First, we notice that, similarly to aligned pulsar solutions, there are strong current layers that form a Y-point with the equatorial CS (see middle panel). The magnetic field lines  beyond this region are not necessarily open, since the equatorial CSs do not extend continuously outside, but instead posses a spiral structure (as previously seen in Fig. 2). In contrast to a spinning NS, we find here that the orbital motion induces a pattern of alternating signs in the charges and parallel currents inside the light cylinder. This is qualitatively consistent with the current distributions shown in Fig. 12 of Ref. [20], especially for the case of a weakly magnetized NS companion "U=u" where we can interpret this alternate pattern as arising from the orbit (MD) effect and the strong currents connecting the two NSs as representing the UI mechanism. A similar pattern is also found for the quantity ϱ (right panel), meaning that the four-current alternates character from timelike to spacelike, the later being the regions where counterstreaming of charged particles would be expected. Finally, and arguably the most relevant difference of the orbital case with respect to pulsar magnetosphere solutions is the fact that there are almost no charges present nor currents flowing along the magnetic field lines in the polar region of the NS.

B. Spin effects
In order to understand the impact of the NS spin, we shall fix a set of representative parameters for the orbital motion and consider different values for Ω Ã . In particular, we pick up an orbit with frequency Ω o ¼ 0.25c=R o and separation R o ¼ 2.5R Ã , and set the spin to Ω Ã ¼ κΩ o , varying κ from −1 to 1. Note that choosing negative values of κ corresponds to antialigned orbital and spin angular momenta. We measure the luminosity and normalize it with L orb , which corresponds to its value for a purely circular orbit at angular velocity Ω o (i.e., χ ¼ 0, Ω Ã ¼ 0). Figure 8 shows the results, together with a curve (dashed black) that represents an estimation from the simple addition of the orbital and spin contributions (i.e., L orb þ L spin ), and a fitting (red solid line) of the form where w 0 ¼ 4 15 ðv o =cÞ 2 fðv o =cÞ ≈ 0.042 comes from the pure orbital part, while the last term represents the pure spin contribution. The other coefficients are fitted from the numerical data, giving w 1 ∼ 0.06, w 2 ∼ 0.35, and w 3 ∼ 0.23. They account for the nontrivial superposition of the two dynamical effects. We note that for the case in which the spin and orbital angular momenta are antialigned, the resulting luminosity is very close to the direct sum of each contribution, whereas for aligned cases, there is an extra enhancement. Also, it is worth mentioning that the orbital motion produces typically much weaker luminosity than the NS spin alone (i.e., pulsars), for this particular setting being approximately 5% of the pulsar spindown luminosity. The reason for this is that the luminosity by the orbital motion is ðv o =cÞ 2 smaller than that by the spin motion from the post-Newtonian viewpoint. On the other hand, we also notice that when the orbit and spin are synchronized, then the pulsar luminosity can be significantly enhanced (for this case, on about 70%).
In Fig. 9, we display the charge and current distributionslike the previous plots in Fig. 7-for the case in which the spin is synchronized with the orbit (Ω Ã ¼ Ω o ). As seen before, for the synchronized motion, the spin effect tends to "dominate" over the orbital one, thus resulting in a magnetospheric configuration very similar to that of pulsars. Both electric charge and parallel current distributions are consistent with typical aligned pulsar results (see, e.g., Refs. [27,47]): strong current layers manifest as a local enhancement of the parallel current components j k (middle panel), with further currents (of opposite signs for each hemisphere) flowing along the magnetic field lines at the polar regions. The intensity of these currents is approximately twice stronger than those in the pure spinning setting, in line with the luminosity being also about two times the pulsar spin-down luminosity. Another significant difference is again the location of the Y-point at each side of the NS, which we find to be approximately given byr ≈ c=Ω o AE R o . Spacelike currents are present (in red color, right panel) along the current layers, surrounding the equatorial CS and also within a small region inside the polar cap.

C. Misalignment effects
We focus now on the cases in which the magnetic axis is not aligned with the orbital angular momentum. To get the main aspects of the influence of this angle, without exploring the whole parameter space, we shall fix the same parameters Ω o and R o as in Sec. III B, and consider several angles χ ¼ f0°; 30°; 60°; 90°g for two NS spins, Ω Ã ¼f0;Ω o g. These two limiting cases display very different behaviors, as shown in Fig. 10: for the synchronized motion with Ω Ã ¼ Ω o  (bottom panel), the spin effect completely dominates--as expected from the previous section results-and thus, the resulting luminosity resembles that of an oblique pulsar, i.e., L ≈ Lðχ ¼ 0Þð1 þ sin 2 χÞ (see, e.g., Ref. [26]). On the other hand, in the pure orbital setting, Ω Ã ¼ 0 (top panel), the luminosity at a given radius exhibits a strong phase dependency. 5 Inside the region given byr ≲ c=Ω o , this can be modeled quite well with L ≈ Lðχ ¼ 0Þð1 þ 0.5 sin 2 χÞ × ð1 − 0.93 sin 2 ðΩ o tÞ sin 2 χÞ. Figure 11 illustrates the magnetosphere with Ω Ã ¼ 0, after 2.5 orbits, for three misalignment angles. Although the configurations are intrinsically three dimensional, some insight can be gained from considering them at the co-orbiting plane (defined by the orbital axis and the vector pointing to the NS from the CoM, in this case, the y ¼ 0 plane). In color, we represent the electric current component along the magnetic field, together with some representative magnetic field lines projected to the plane. Alternate patterns on the currents can be seen again inside the close zone, like in the aligned setting. Such currents become particularly intense (j k ∼ 4Ω o B=2π) at the CS, that oscillates here about the dipole equator. These shapes are reminiscent of those of oblique pulsars, although there the CSs oscillate about the rotational equator instead. Finally, we notice that for the orthogonal case (i.e., χ ¼ 90°), the CS looks steady and smooth, in contrast to the intermediate inclinations for which there are signs of magnetic reconnection activity and plasmoids. 6

D. Inspiral orbits
In this section, we analyze the evolution of the luminosity during the late stages of an inspiral orbit. Even though in this paper we are not including curvature effects, we shall mimic the inspiral phase of a binary system by taking relevant astrophysical parameters to set the NS trajectory around the CoM. To that end, we consider parameters similar to those used in Ref. [22]. Specifically, we take the mass ratio of the BHNS binary to be The total mass is thus given by M ¼ 5.6 M ⊙ . Instead, for the BNS trajectory, we take an equal mass system, i.e., q ¼ 1, and hence, M ¼ 2.8 M ⊙ . Assuming quasicircular inspiral orbits, the angular frequency evolves approximately according to Ref. [48], with Ω o;0 being the initial orbital frequency at t ¼ t 0 . Our orbital radius then reads ð18Þ with the orbital separation, a, through R o ≡ q 1þq a. We start from a quasistationary initial configuration with R o;0 ≃ 5.5R Ã (and its corresponding Ω o;0 ) at t ¼ t 0 , following the inspiral evolution from that moment on, until reaching a final radius R o ∼ 2R Ã . The expected disruption radius in our setting would be R dis o ≃ 2.4R Ã , as estimated [18] from a dis ≃ 3M BH ðC=0.2Þ −1 ðq=7Þ −2=3 .
The evolution of the luminosity during the trajectory for several extraction radii is shown in Fig. 12 (top panel), where time is presented in physical units and relative to t d , the time for which the estimated disruption radius R dis o is attained. The luminosity is given in physical units, as well, for a typical magnetic field strength at the stellar pole of 10 12 G. Since luminosity scales exactly as B 2 in force-free FIG. 11. Effects of misalignment on the magnetosphere with Ω Ã ¼ 0. Parallel electric currents, 2πj k =Ω o B, on the y ¼ 0 plane after 2.5 orbits are shown in color scale, for angles χ ¼ 30°(left panel), χ ¼ 60°(middle panel), and χ ¼ 90°(right panel). Black lines show some representative magnetic field lines along the plane. 5 Note that such phase modulations also appear in the analytic estimations for vacuum magnetospheres, arising from a term proportional to ð⃗ μ · _ ⃗ aÞ 2 (see the second term of Eq. (A6) in Ref. [7]). 6 This fact should be taken with some caution, since it may be indicating a lack of numerical resolution/dissipation. electrodynamics, their values can be rescaled by B 2 12 ≡ ð B 10 12 G Þ 2 , as indicated in the plot. From the initiation of the inspiral motion, the luminosity is no longer constant inside the light cylinder approximately c=Ω o ðtÞ, since the rising of the luminosity initiates at the star and it takes some time to propagate outward. This is clearly seen from different curves representing integration of the flux at several spheres enclosing the NS. Note that the luminosity measured atr ¼ 50R Ã departs from the other values from the very beginning, reflecting the fact that a significant fraction of the Poynting flux gets dissipated at the CS. Such a difference is further enhanced later, reaching approximately 2 orders of magnitude (with respect to the value at r ¼ 5R Ã ), by the time-delay effect already mentioned. The extra (black dashed) curve, included in the plot, depicts the predicted/estimated values from our previous results of stationary circular orbits [see Eq. (15)]. It can be noted that such estimated curves match with the dynamical results quite well (especially when the comparison is made taking the values at the integration radius close to approximately c=Ω o ðtÞ), suggesting that the timescale of the magnetospheric response is comparable to the inspiral one during the whole evolution.
In Fig. 13, we illustrate the three-dimensional (radial) Poynting flux density along the inspiral BHNS trajectory. The values have been normalized with L 0 ðtÞ=4πr 2 , for better visualization and to facilitate the comparison among the different snapshots. We find that-as in the stationary circular orbits-the flux distribution through a sphere of radius r ¼ c=Ω o ðtÞ remains broadly unchanged, only rotating at the orbital frequency Ω o ðtÞ. The outgoing fluxes outside the "instantaneous light cylinder," i.e., c=Ω o ðtÞ, form spiral structures over the x-y plane(s), producing the equatorial CS (with magnetic reconnections) as these Alfven fronts propagate outward. We now estimate how the inclusion of curvature (associated to both the NS and the BH) can impact on the luminosity. To that end, we shall pick up a state at an orbital separation of a ≈ 6.6R Ã (represented in the middle panel of Fig. 13) and compare its luminosity with the results reported in Ref. [22] for essentially the same parameters. Our results show L ≈ 2.4 × 10 39 B 2 12 erg=s, which should be associated in our case exclusively to the MD acceleration mechanism. On the other hand, they obtained L ≈ 6.2 × 10 40 B 2 12 erg=s, for the case where both compact objects were nonspinning. Hence, taken at face value, one might say that the luminosity is enhanced by a factor approximately 25 due to curvature effects. However, a few important considerations are necessary. First, in Ref. [22], the luminosity is measured at large distances, r ≳ 90R Ã , and not inside the light cylinder where we compute it. As we found dissipation at CSs outside the light cylinder, the comparison is not direct, and perhaps even slightly larger luminosity could be expected at the light cylinder when curvature is accounted for. On the other hand, it is not very clear how to disentangle from different curvature effects like unipolar induction, estimated in this context to give L UI ≈ 6.4 × 10 39 B 2 12 erg=s [22]. Similarly to the BHNS binary scenario, we also follow the evolution of the luminosity for an equal mass BNS system (see the bottom panel of Fig. 12). The trajectory is evolved from an orbital separation a i ∼ 6R Ã up to a final one a f ∼ 3R Ã (which sets the reference time t f ). Thus, we consider a complementary dynamical range with respect to the one studied in Refs. [19][20][21], where the simulations begin at about our final orbital separation and proceed though the late-inspiral phase until (and after) the merger. Although a comparison of the luminosity here is not as direct as before, we can estimate from their case "U=u" of a weakly magnetized companion, L ≈10 42 B 2 12 erg=s at orbital separation a ∼ 3R Ã , whereas we obtain L ≈ 10 41 B 2 12 erg=s. Therefore, the enhancement attributed to the effects of curvature is again approximately an order of magnitude.

E. Implication to observation
Here, we briefly discuss the implications of our numerical results. As shown in Sec. III A, CSs are always developed in the local wave zone, near the radius c=Ω o , for orbiting NSs. This feature is shared with the pulsar magnetosphere, for which CSs are also developed just outside the light cylinder at radius. This suggests that some EM signals similar to the ones in pulsars are likely to be emitted from the orbiting NSs, even in the absence of NS spin. In the following, we estimate the luminosity of the EM signals referring to the latest studies of magnetospheric emissions from pulsars.
Before going ahead, we show approximately how much EM energy can be emitted in total during the inspiral phase. For this estimate, we simply employ Eq. (14) with the orbital evolution determined by the Newtonian gravity plus gravitational radiation reaction via quadrupole formula as (e.g., Ref. [48]) where η is the symmetric mass ratio, i.e., the ratio of the reduced mass to the total mass, M, of a binary and written as q=ð1 þ qÞ 2 . Then, the maximum total energy dissipated in EM waves is calculated as Z where B is the magnetic field strength at the pole (i.e., B ¼ μ=2R 3 Ã ); a i and a f denote the initial orbital radius at t ¼ t i and the orbital radius at the onset of merger at t ¼ t f , respectively; and we assume that a i ≫ a f . In this equation, we suppose that M BH ≈ 4.2 M ⊙ , M NS ≈ 1.4 M ⊙ , and a f ≈ 5M. As we have shown in this paper, the luminosity is enhanced in the close binary orbits, and hence, the total energy emitted could be larger by a factor of several than Eq. (20). Nevertheless, Eq. (20) indicates that for the typical magnetic field strength at the NS pole of B ¼ 10 12 G, the available energy for the EM signals is at most approximately 10 40 erg, and thus, for the observation being possible, the presence of an efficient emission mechanism or unusually high magnetic field strength would be necessary.
The latest PIC simulations (e.g., Ref. [30]) show that near the so-called Y-point at the light cylinder of the pulsar magnetosphere, CSs are developed, and as a result, the reconnection of magnetic field lines is enhanced. In such a region, the magnetic field strength and number density of the electron-positron pair are significantly increased. Highenergy electrons and positrons accelerated in the strong magnetic field become the sources of the synchrotron and inverse Compton radiations. Because the particle energy is high, the energy can be extended to the MeV and GeV gamma-ray bands. Such an emission is consistent with the presence of the gamma-ray pulsars observed by the Fermi satellite for the isolated pulsar case [49]. The luminosity of the gamma-ray pulsars could be 1%-100% of the spindown luminosity of the pulsars. In the present context, the spindown luminosity should be replaced by the total luminosity of the orbiting NS, which is L ≲ 10 42 erg=s for typical magnetic field strength of order 10 12 G at the poles. Thus, it is reasonable to consider that the gamma-ray luminosity would be at most approximately 10 42 erg=s, for which the duration is ≲10 ms. It would be quite difficult to detect such a low-luminosity gamma-ray source with very short duration using the current and near-future gamma-ray telescopes, if we suppose that the typical distance to the source is ≳100 Mpc. This would be also the case for the observation of x rays. If the magnetic field strength is as high as that of magnetars, approximately 10 14 G at its pole, the luminosity could be enhanced, and the gamma-ray signals may be observable. However, such high magnetic fields may not be very likely for inspiraling NSs, as the observational results of NSNS in our Galaxy indicates [50].
As discussed in Sec. III D, Ref. [22] suggests a significant enhancement of the EM luminosity if curvature effects, i.e., interaction of the magnetic field with the BH, are included: at least an order of magnitude seems possible. Although it is not clear how this extra luminosity is distributed in terms of the magnetospheric configuration, that is, how to disentangle from the unipolar induction mechanism, in which the two compact objects are connected through a twisted bundle of magnetic field lines and thus a dc circuit established in the near zone, localized plasma winds in between the two compact objects are expected to produce curvature radiation dominating in the γ-rays band, and synchrotron emissions in the hard x-rays and soft γ-rays bands [15]. It was also suggested that a fraction of this flux within the bundle, carried in the form of plasma kinetic energy, will reach the primary NS surface and form a hot spot emitting thermal energy as x rays [15]. This flux-tube UI configuration was further proposed to model fast radio bursts arising from coherent curvature radiation in the late-inspiral phase (but, prior to the last few orbits) in NSNS binaries [51]. All these speculations encourage us to perform force-free simulations including general relativistic effects, and also to model this radiation along the lines of Ref. [27] to obtain light curves from our solutions.
The latest work also shows that in the reconnection region, radio waves may be emitted by coalescence of magnetic islands in the CS that produces magnetic perturbations propagating away [28,31,52]. For this process, the predicted spectral flux density of the radio waves with the frequency ν is where B L denotes the magnetic field strength for the reconnection region, Γ and l are the bulk Lorentz factor and characteristic scale of the magnetic islands, D is the distance to the source, and τ ∼ 10=ω p with ω p the plasma frequency. Here, l is determined by the force balance and energy conservation in the magnetic islands, and Γ is inferred as 10-100 [28,52]. In the present context, we obtain a quite small flux for orbits close to the innermost stable circular orbit, Here, we suppose that the total mass of the system is 5.6 M ⊙ , a ≈ 100 km, and the magnetic field strength at the NS pole is B ¼ 10 12 G. We also simply set τν ¼ 1. The expected value of l in the reconnection region of the pulsar magnetosphere is about ten times larger than the thickness of the CS, which is proportional to B −3=2 L as [28,52] ∼0.1 cm This is the reason why, for the strong magnetic fields in the reconnection region, the intensity of this type of the radio waves is low (i.e., in this scenario, the luminosity decreases with the decrease of the orbital separation, a). If the magnetic field strength for the late-inspiral stage of NS binaries is B ≈ 10 10 G, the luminosity would be 10 5 times lager than that of Eq. (22) because the width of the magnetic islands becomes larger. However, the flux is still approximately 1 mJy.
Since the size of each magnetic island is likely to be quite small for the strong magnetic field case, a large number of the magnetic islands may be simultaneously generated. As discussed in Ref. [52], fast waves from many merging magnetic islands may nonlinearly interact and transfer the energy into the plasma. If this happens and the emission occurs in an optimistically coherent way, the predicted luminosity is written as [52] which approximately leads to the following flux density: Thus, a much higher flux than Eq. (22) can be predicted, although the flux is still low for the detection by the current wide-field-of-view radio telescopes like CHIME [53] and OVRO-LWA [54]. We note that in a small fraction of pulsars [55], a substantial fraction of the dissipation energy of order 10 −4 is emitted in the radio band. If the orbiting NSs can have such a high efficiency for the radio emission, the radio waves may be observable.
As indicated above, it will not be very easy to detect an EM precursor of NS mergers. However, the following point should be kept in mind: in the absence of tidal disruption in the system of BHNS binaries, no EM counterparts are expected after the merger. Even in such cases, a precursor associated with the moving NS with magnetic fields can be emitted as an EM counterpart of the merger of BHNS binaries. Distinguishing a BHNS binary with a fairly large chirp mass from binary black holes is not an easy task only in the detection of gravitational waves, because the gravitational waveforms for two cases are quite similar. The precursor for BHNS binaries is likely to be emitted always if the magnetic field strength of the NS in binaries is as strong as that for the typical isolated NSs. The observation of the precursor for the BHNS binaries will play an important role for surely identifying the BHNS binaries.
It is also worth mentioning that if a precursor is detected for NS binaries, the magnetic field strength of the NS could be estimated. This will provide us important information for the evolution of the magnetic fields in an old NS that has not experienced significant mass accretion.

IV. CONCLUSIONS
In this paper, we have considered force-free magnetosphere induced by a NS orbiting in a binary system, aiming at capturing the EM effects produced by the orbital motion about the CoM of the binary (BHNS or NSNS) system, in a sense, isolating the role played by the acceleration of the MD moment of the NS, from the effects of curvature. The inclusion of curvature-more specifically a BH companion-has been deferred for a subsequent work. Our present approach, however, has allowed a detailed and systematic study of the magnetospheric properties of these systems in close analogy to pulsars. In particular, the existence of strong return current layers and CS enabled us to connect with known EM emission mechanisms from pulsar theory.
Before considering realistic inspiral orbits associated with the binaries, we first analyzed the properties of circular orbits, noticing--as later confirmed-that the stationary configurations attained here are reliable approximations of the system's states throughout the inspiral. The results can be summarized as follows. For the aligned and nonspinning scenario, a strong spiral CS develops along the orbital plane, beginning at approximately c=Ω o . An alternate pattern of charges/currents is found in the near zone of the NS, enclosed by thin return current layers that form a Y-point with the equatorial CS. In contrast to pulsars, there is almost no charge flowing over the polar region of the NS. The luminosity can be well described by L≈ 4 15c 5 μ 2 R 2 o Ω 6 o fðv o =cÞ, where f ∼ ð1-7Þ represents relativistic corrections in v o =c, to the analytic formula (14). The orbital effect on the luminosity is much weaker than the one produced by the NS spin (comparing at a given angular frequency), although it can produce a considerable enhancement (factor approximately 2) on the pulsar spindown luminosity if the motion is synchronized. Finally, for the case that the magnetic and orbital axis are not aligned (i.e., χ ≠ 0), the magnetosphere resembles that of oblique pulsars, with a strong CS fluctuating about the dipole equator (instead of the rotational equator, as in misaligned pulsars [26]). Interestingly, in the nonspinning case, the total luminosity acquires an orbital-phase dependency.
We consider inspiral binary trajectories using relevant parameters for both BHNS and NSNS scenarios, finding that the evolution of the luminosity follows quasiadiabatically from our previous estimations based on the circular orbits configurations. The Poynting fluxes emanating from the orbiting NS form spiral structures orthogonal to the orbital axis and produce magnetic reconnections inside the CS as these Alfven fronts propagate outward. Overall, the radial Poynting flux distribution in the wave zone establishes a lighthouse effect, rotating at orbital frequency Ω o ðtÞ and peaking around the orbital plane, which is in qualitative agreement with previous GRFF studies for the cases where the companion to the primary NS is either a weakly magnetized NS [19,20] or a BH [22]. This suggests that the EM flows beyond the light cylinder approximately c=Ω o ðtÞ are dominated by the MD effect, while UI occurs almost exclusively between the two compact objects in the near region. One might interpret these two mechanisms as constituting two approximately independent sources of EM energy for the plasma (albeit both mined from the kinetic energy of the orbital motion). However, by the interaction of the magnetic field with the curvature of the companion (if curvature is included), the field strength can be amplified by further twisting, resulting in enhanced luminosity of the MD effect-valued here to represent 1 order of magnitude-.
Furthermore, the contribution from each of these two mechanisms (i.e., MD/UI) will be channeled into separate processes within the magnetosphere, which then derive in the actual EM emissions. In particular, we have devoted our attention to the existence of strong CSs induced by the NS's orbital motion alone. These CSs were not reported before, in the BHNS scenario of Ref. [22], this may be due to the fact that reflection symmetry about the orbital plane was imposed in their simulations, whereas for the weakly magnetized NS companion case studied in Refs. [19,20], we see evidence of such CSs forming (last panel of Fig. 10 in Ref. [20]), although it was omitted in their discussions. 7 We argue these CSs are very important due to their role in explaining radio emissions in pulsars, where recent studies [28,31,52] have shown that the coalescence of magnetic islands (or plasmoids) in the reconnection regions can 7 Only a trailing dissipation region, making a CS tail behind the weakly magnetized NS companion, was described in these works. produce radio waves, as a coherent superposition of many such individual pulses. We have borrowed these ideas and applied them to get concrete estimations of spectral flux densities of radio waves in our context. The predicted fluxes are not large enough to be detected by the current wide-field-of-view radio telescopes like CHIME [53] and OVRO-LWA [54]. We conclude that only if a substantial fraction of the luminosity (e.g., ≳10 −4 ) is emitted in the radio band, as in some pulsars [55], then radio waves may be observable.