format: text

text:
# Force Free Waves and Black Hole Magnetospheric Causality

## 1 Introduction

Applications of force free electrodynamics are common in the astrophysical literature. The assumption that the plasma inertia is negligible compared to the magnetic energy density is an accurate depiction of the plasma state in many contexts. This condition justifies the force free approximation based on a local physical frame analysis. The most notable successes of such an approach are descriptions of equilibria in the solar magnetosphere. Motivated by the solar work, it has also been extended to studies of the physics of black hole magnetospheres (see Blandford and Znajek (1977) and Thorne et al (1986)) and more recently gamma ray bursts (see Blandford (2002)).

However, the implementation of the force free condition to black hole magnetospheres is not so trivial because general relativity is not a local theory. In particular, the concepts of an event horizon and ergosphere are globally defined (Hawking, 1973). Black hole magnetospheres are the main focus of this treatment, yet all of the results that are derived in this article apply equally well to any force free calculation. A salient example of the global nature of general relativity is the fact the physical boundary condition at the event horizon is that all particles that are near the black hole must corotate with the event horizon and go into the black hole on trajectories that approach the ingoing light-like principal null congruence, as seen by any external observer (see Chapter 3 of Punsly (2001)). This applies to any plasma, no matter how tenuous and is independent of the magnetic field strength. There are no stable equilibria near the horizon, gravity always overwhelms the other forces. This property is precisely what distinguishes a black hole from another type type of compact object such as a neutron star or white dwarf. All plasma near the horizon is inertially dominated: a global condition. In a sense, the plasma inertia is infinite near the horizon as viewed globally. By contrast, a local analysis would lead one to believe that the force free approximation is valid (i.e., the plasma inertia is negligible). In spite of this defining property of the black hole, force free treatments of black hole magnetospheres that extend all the way to the event horizon have been proposed by many authors (most notably Blandford and Znajek (1977) and Thorne et al (1986)).

In this paper, we consider the plasma waves in a force free magnetosphere. This is important in order to clarify certain anecdotal defenses of the force free approximation near the event horizon which infer that force free plasma waves can transport different information than their counterparts in the full perfect magnetohydrodynamic (MHD) treatment of plasma waves in magnetically dominated magnetospheres (see Chapter 2 of Punsly (2001) for a detailed MHD treatment and Blandford (2001) for the force free discussion). In this article, it is shown by means of a straightforward explicit calculation that this claim is erroneous. This article elucidates the current structure within the individual wave modes for the first time (even the lengthy, complicated and very detailed analysis of certain topics of Uchida (1997) in the context of black hole magnetospheres does not address these most relevant issues in a force free magnetosphere).

## 2 The Force Free Approximation and Magnetically Dominated MHD

This section compares and contrasts the force free approximation with a magnetically dominated perfect MHD system. Generally, the force free approximation is considered as a zero mass limit of the full MHD theory. It is usually considered an expedience that reduces the dimension of the full MHD system of hyperbolic differential equations to five. It is very rarely argued to be a more accurate physical description than the full MHD theory (but, see Blandford (2002) for such an argument). Generally, the force free approximation is invoked is to simplify the mathematical treatment not to improve on the physical depiction. If the force free approximation is not pathological then one would expect it to occur as the natural limit of zero inertia in the full MHD theory. This section begins to lay the groundwork to show that this is indeed the case. The set of force free equations are the current constraint written covariantly in terms of the Maxwell tensor and in a local inertial coordinate system,

$`F^{\mu \nu }J_\nu =0,`$ (2-1a)
$`\rho _e𝐄+{\displaystyle \frac{𝐉\times 𝐁}{c}}=0,`$ (2-1b)
and Maxwell’s equations (written in local coordinates)

$`\times 𝐁={\displaystyle \frac{1}{c}}{\displaystyle \frac{𝐄}{t}}+{\displaystyle \frac{4\pi }{c}}𝐉,`$ (2-2)
$`\times 𝐄={\displaystyle \frac{1}{c}}{\displaystyle \frac{𝐁}{t}},`$ (2-3)
$`𝐁=0,`$ (2-4)
$`𝐄=4\pi \rho _e.`$ (2-5)
Equation (2.1) also implies the conditions:

$`𝐄𝐁=0,`$ (2-6)
$`𝐄𝐉=0.`$ (2-7)
The force free condition in equation (2.1) arises from the zero mass density limit of perfect MHD, the plasma has no inertia to stress the magnetic field lines.

## 3 Force Free Plasma Waves

The causal structure of a force free system is governed by the waves supported in the force free plasma. In this section we derive the nature of the information that can be transported along each of the wave characteristics. Recall that in perfect MHD the Alfven mode can transport charge and current that was either field aligned (force free) or cross-field (inertial). The fast wave cannot transport charge or field aligned currents, only inertial currents (see Punsly (2001), Chapter 2). It turns out that the effects of the force free constraint are equivalent to the zero mass density limit of perfect MHD plasma waves, as expected. Since the fast wave still cannot transport force free current in the zero inertia limit, when the force free constraint on the current is imposed, the only value of current consistent with these two concepts (effectively two mathematical equations) is zero current. The fast mode in the force free limit carries no current or charge.

As a consequence of equation (2.6), one can find a frame in which the electric field vanishes if the Lorentz invariant $`F^{\mu \nu }F_{\mu \nu }>0`$ (i.e., the field is magnetic). There exists a coordinate transformation to this frame that resembles a Lorentz boost in the azimuthal direction in an axisymmetric magnetosphere except that the boost parameter can exceed the speed of light. Thus, this is not a Lorentz transformation in general. However, the angular velocity as observed from asymptotic infinity in an asymptotically spherical coordinate system is referred to as the field line angular velocity, $`\mathrm{\Omega }__F`$, even when the azimuthal velocity is faster than the speed of light (Blandford, 2002; Punsly, 2001). Thus, this frame is often called the corotating frame of the magnetic field. The light cylinder is defined as the surface at which the azimuthal boost velocity equals the speed of light (in all physical frames). The implementation of this frame facilitates the calculations and one can perform a general coordinate transformation in order to evaluate plasma quantities in other physical frames. Since, this is not a Lorentz transformation in general, it is best to find coordinate invariant expressions. Another frame in which the electric field vanishes is the plasma rest frame. However, this frame is poorly defined in the force free limit since there is no accounting of the plasma state. Frames in which the electric field vanishes are sometimes called proper frames.

The calculations are performed in the short wavelength approximation. This is equivalent to the analysis found in Blandford (2002). No matter what the radius of curvature of the magnetic field, one can find plasma mode wavelengths that are short compared to this value. This is not very restrictive, even close to the event horizon of a black hole. One should be able to find proper frames in which the radius of curvature of the magnetic field is on the order of the geometrized mass of a black hole ( $`10^{14}`$ cm for a supermassive black hole). The unperturbed value of the four current, $`J_0^\mu `$, effectively vanishes in the short wavelength approximation that the waves vary on much smaller scales than the background magnetic fields (the currents are effectively created by the gradients of the magnetic field). As a result of Ampere’s law, this approximation is equivalent to treating the magnetic field of the unperturbed state as approximately a constant in the proper frame,

$`{\displaystyle \frac{}{x^\alpha }}𝐁0.`$ (3-1)
We analyze short wavelength linear (small amplitude) perturbations of equations (2.1)-(2.7). The waves are chosen in a plane wave representation to have an oscillatory behavior, $`\mathrm{exp}i(𝐤𝐫\omega t)`$. We define the $`x`$ direction to be parallel to $`𝐤`$, so all perturbed quantities are Fourier analyzed in $`(x,t)`$ space. For example, the four current is

$`J^\mu =J__0^\mu +e^{i(kx\omega t)}\delta J^\mu `$ (3-2a)
$`e^{i(kx\omega t)}\delta J^\mu .`$ (3-2b)
The interested reader can find a discussion of the effects of field line curvature in a black hole magnetosphere in Chapter 6 of Punsly (2001). The nature of the resultant hybrid modes is beyond the scope of this treatment that has been chosen to parallel the analysis of Blandford (2002). In summary, field line curvature mixes the properties of the plasma modes as it does in the earth’s magnetosphere (Lanzerotti and Southwood, 1979). Long wavelength modes are neither fast nor Alfven modes, but have hybrid properties (Appendix A is a derivation of the first order effects of magnetic field line curvature in the force free limit). The faster mode propagates slower as the wavelength increases and is no longer related to the critical point structure of the winds that reflect short wavelength mode speeds. The implication for black hole magnetospheres is that long wavelength force free fast modes travel slower than the speed of light and therefore can not propagate information from near the event hroizon to asymptotic infinity.

The analysis of long wavelength plasma modes is very complicated and one is referred to the in depth analysis in Chapter 6 of Punsly (2001) to see the details of outgoing mode propagation resulting from field line curvature and spacetime curvature near the event horizon of a black hole. The results of those calculations are in accord with the causality arguments in this article. In particular, one finds that long wavelength plasma modes are extremely inefficient in their ability to transport magnetic stresses, field aligned current and charge outward from the spacetime near an event horizon.

In a frame with a vanishing electric field, the perturbed equations (2.1) - (2.7), with the condition above become:

$`𝐁\times \delta 𝐉=0,`$ (3-3)
$`\omega \delta 𝐄=c𝐤\times \delta 𝐁i4\pi \delta 𝐉,`$ (3-4)
$`\omega \delta 𝐁=c𝐤\times \delta 𝐄,`$ (3-5)
$`𝐤\delta 𝐁=0,`$ (3-6)
$`i𝐤\delta 𝐄=4\pi \delta \rho _e,`$ (3-7)
$`\delta 𝐄𝐁=0.`$ (3-8)
It is useful to combine (2.2) and (2.3) to get the second order Maxwell’s equation,

$`^2𝐄(𝐄)={\displaystyle \frac{1}{c^2}}{\displaystyle \frac{^2𝐄}{t^2}}+{\displaystyle \frac{4\pi }{c^2}}{\displaystyle \frac{𝐉}{t}}.`$ (3-9)
This equation gives the following relation amongst perturbed electromagnetic quantities,

$`k^2\delta 𝐄𝐤\left(𝐤\delta 𝐄\right)={\displaystyle \frac{\omega ^2}{c^2}}\delta 𝐄+{\displaystyle \frac{i4\pi \omega }{c^2}}\delta 𝐉.`$ (3-10)
At this point we explicitly write out the components of the perturbed equation above in a local Cartesian basis. The x-axis was previously defined to be parallel to $`𝐤`$ and we define the z-axis such that $`𝐁`$ always lies in the x-z plane. The components of the second order Maxwell’s equation become:

$`\delta E_y={\displaystyle \frac{4\pi i\omega }{\left(c^2k^2\omega ^2\right)}}\delta J_y,`$ (3-11a)
$`\delta E_z={\displaystyle \frac{4\pi i\omega }{\left(c^2k^2\omega ^2\right)}}\delta J_z,`$ (3-11b)
$`\delta E_x={\displaystyle \frac{4\pi i}{\omega }}\delta J_x.`$ (3-11c)
Similarly an expansion of the components of the perturbed force free condition in (3.2) yields the following two conditions on the current:

$`\delta J_y=0,`$ (3-12)
$`\delta J_z=\mathrm{tan}\theta \delta J_x,`$ (3-13)
where $`\theta `$ is the angle between the propagation vector and the magnetic field direction. There is also a third constraint on the current that is obtained by substituting the perturbed second order Maxwell’s equation (3.11) into the perturbed degeneracy condition (3.8):

$`\delta J_x\left(c^2k^2\omega ^2\right)=\omega ^2{\displaystyle \frac{k_{}}{k_{}}}\delta J_z,`$ (3-14)
where $`k_{}`$ and $`k_{}`$ are the components of the propagation vector parallel and perpendicular to the magnetic field direction, respectively.

Equation (3.14) can be substituted into the relation (3.13) in order to eliminate $`\delta J_x`$. This yields the dispersion relation for plasma waves in a force free medium,

$`\left[\omega ^2c^2k^2\right]\left[\omega ^2c^2k_{}^2\right]=0.`$ (3-15)
### 3.1 The Fast Mode

The dispersion relation (3.15) has two solutions for the phase velocity. The fast mode is characterized by the phase velocity,

$`{\displaystyle \frac{\omega }{k}}=\pm c.`$ (3-16)
The group velocity is

$`{\displaystyle \frac{\omega }{𝐤}}=\pm c{\displaystyle \frac{𝐤}{k}}.`$ (3-17)
Inserting the value of the phase velocity into the perturbed induction equation, (3.5), and the perturbed degeneracy condition, (3.8), along with the vanishing divergence of the magnetic field, (3.6), and the constraint on the perturbed current, (3.14), we can solve for the electromagnetic components of the fast wave:

$`\delta B_x=0,`$ (3-18a)
$`\delta B_y=0,`$ (3-18b)
$`\delta B_z=\delta E_y,`$ (3-18c)
$`\delta E_x=0,`$ (3-18d)
$`\delta E_y=\delta B_z,`$ (3-18e)
$`\delta E_z=0.`$ (3-18f)
Substituting the values of the electric field from (3.18d-f) into the second order Maxwell’s equation, (3.11), and the perturbed Gauss’ law, (3.7), the light-like phase velocity of the fast wave requires that all of the components of the four-current density that are propagated along the characteristics of the force free fast wave vanish identically,

$`\delta J^\mu =0.`$ (3-19)
This is a coordinate invariant expression and is therefore true in all coordinate systems. The fast mode is linear in the force free limit, since there is a unique wave propagation speed given by (3.17). Thus, fast wave packets are linear superpositions of the oscillatory solutions. Therefore, the four-current density vanishes for fast wave packets as well.

From (3.18), one also has the following coordinate independent statements:

$`\delta F^{\mu \nu }\delta F_{\mu \nu }=0,`$ (3-20a)
$`\delta F^{\mu \nu }\delta F_{\mu \nu }=0,`$ (3-20b)
$`\delta F^{\mu \nu }F_{\mu \nu }=0,`$ (3-20c)
$`\delta F_{\mu \nu }U__G^\nu =0,`$ (3-20d)
$`{\displaystyle \frac{\omega }{c}}U__G^\nu =k^\nu ,`$ (3-20e)
where $`U__G^\nu `$ is the group four velocity of the wave, which is expressed in the coordinates of the proper frame as:

$`U__G^\nu (c,{\displaystyle \frac{\omega }{𝐤}}).`$ (3-21)
Not all of the relations in (3.20) are independent, for example (3.20d) implies (3.20a).

The primary result of this analysis, (3.19), that demonstrates the vanishing of the current in a fast mode is not a result of a biased interpretation by this author, a miscalculation or an unfounded assumption. The same result actually follows from the results of Blandford (2002). First of all, the wave structure that is described in a coordinate independent manner in (3.20) can be rewritten in the “laboratory frame coordinates” used by Blandford (2002). One then finds that, $`𝐤`$, $`\delta 𝐁`$ and $`\delta 𝐄`$ form an orthogonal triad as well as equations (7) and (8) of Blandford (2002). These are the results that are used in that treatment to describe the causal relevance of a fast wave to a force free magnetosphere. Inserting (7) and (8) into Ampere’s law, either (2.2) or (3.4), yields the fact that $`\delta J_y=0`$ and $`\delta J_z=0`$. By the orthogonality of the triad, $`\delta E_x=0`$, as a consequence of the definition of the propagation direction. Then Ampere’s law can be used again to show that $`\delta J_x=0`$. Finally, the vanishing of the charge density follows from either Gauss’ law or current conservation. Thus, (3.19) can be found from the analysis of Blandford (2002) quite readily.

### 3.2 The Alfven Wave

The other root of the dispersion relation (3.15) gives the phase velocity of the of the Alfven wave:

$`{\displaystyle \frac{\omega }{k}}=\pm c\mathrm{cos}\theta .`$ (3-22)
One can also compute the group velocity,

$`{\displaystyle \frac{\omega }{𝐤}}=\pm c{\displaystyle \frac{𝐁}{B}}.`$ (3-23)
Thus, the Alfven mode also propagates information at the speed of light, but only along the direction of the magnetic field.

Combining, the vanishing divergence of the perturbed magnetic field, (3.6), along with the constraint on the perturbed current, (3.14), and the induction equation, (3.5), inserted into the second order Maxwell’s equation, (3.11), one can solve for the components of the perturbed electromagnetic field transported along the characteristics of the Alfven mode:

$`\delta B_x=0,`$ (3-24a)
$`\delta B_y=i{\displaystyle \frac{4\pi k}{ck_{}^2}}\delta J_z=i{\displaystyle \frac{4\pi }{ck_{}}}\delta J_{},`$ (3-24b)
$`\delta B_z=0,`$ (3-24c)
$`\delta E_x=i{\displaystyle \frac{4\pi }{\omega }}\delta J_x=i{\displaystyle \frac{4\pi }{ck}}\delta J_{},`$ (3-24d)
$`\delta E_y=0,`$ (3-24e)
$`\delta E_z=\mathrm{cos}\theta \delta B_y=i{\displaystyle \frac{4\pi k_{}}{ckk_{}}}\delta J_{}.`$ (3-24f)
The component of the current along the magnetic field in the proper frame is defined as

$`J_{}=𝐉{\displaystyle \frac{𝐁}{B}}.`$ (3-25)
The perturbed current can be found from (3.22), (3.24) and (3.14),

$`\delta J_x=i{\displaystyle \frac{ck_{}k_{}}{4\pi k}}\delta B_y={\displaystyle \frac{k_{}}{k}}J_{},`$ (3-26a)
$`\delta J_y=0,`$ (3-26b)
$`\delta J_z=i{\displaystyle \frac{ck_{}^2}{4\pi k}}\delta B_y={\displaystyle \frac{k_{}}{k}}J_{},`$ (3-26c)
$`\delta 𝐉_{}=\delta 𝐉=i{\displaystyle \frac{ck_{}}{4\pi }}\delta B_y{\displaystyle \frac{𝐁}{B}}.`$ (3-26d)
From Gauss’ law, (3.7), and (3.24) the perturbed charge density is

$`\delta \rho _e=i{\displaystyle \frac{k_{}}{4\pi }}\delta B_y={\displaystyle \frac{J_{}}{c}}.`$ (3-27)
The second equality in (3.27) utilized the perturbed law of current conservation which is obtained from the combination of (3.7) and (3.4). The equations (3.24)-(3.26) can be combined to give the following useful relations amongst dynamical quantities that are propagated along the Alfven characteristics:

$`\delta 𝐉\times 𝐁=0,`$ (3-28a)
$`\delta 𝐄𝐁=0,`$ (3-28b)
$`\delta 𝐄\delta 𝐁=0,`$ (3-28c)
$`\delta 𝐁𝐁=0.`$ (3-28d)
Covariantly speaking, (3.28) implies the coordinate independent relations:

$`\delta J^\mu F_{\mu \nu }=0,`$ (3-29a)
$`\delta F^{\mu \nu }F_{\mu \nu }=0,`$ (3-29b)
$`\delta F^{\mu \nu }\delta F_{\mu \nu }=0,`$ (3-29c)
$`\delta F^{\mu \nu }F_{\mu \nu }=0,`$ (3-29d)
$`\delta F^{\mu \nu }\delta F_{\mu \nu }=0,`$ (3-29e)
$`\delta F_{\mu \nu }U__G^\nu =0,`$ (3-29f)
$`F_{\mu \nu }U__G^\nu =0.`$ (3-29g)
Equations (3.29a-d) are the coordinate independent analogs of the local expressions (3.28a-d). From (3.29a), the current flows along the magnetic field lines in a proper frame (this expression incorporates the assumption, (3.2), that field line curvature radii are large compared to the plasma mode wavelength). One can also get the orientation of the propagation vector relative to the perturbed and unperturbed fields in all frames by expressing (3.29a) in terms of the perturbed covariant version of Ampere’s Law,

$`\delta F^{\mu \alpha }k_\alpha F_{\mu \nu }=0.`$ (3-30)
The Alfven mode transports a charge and a field aligned current (covariantly speaking, a force free current).

## 4 Causality in Force Free Magnetospheres

The calculation in the last section allows one to compute the characteristics of short wavelength waves in a force free magnetosphere. There are four modes: both ingoing and outgoing fast and Alfven waves. The characteristics by themselves are not enough to time evolve a force free system, one must also know what information can be propagated along these curves. The force free waves carry only two independent pieces of information. In order to see this, note that the system has ten parameters: three magnetic field components, three electric field components, three current components and the charge. There are eight independent relations connecting these quantities - the force free equation (2.1) has three constraints, the second order Maxwell’s equation, (3.10), has three constraints as well and there are the two divergence relations involving the field components in (2.4) and (2.5). This leaves two independent quantities that are readily identified from the results of section 3 as the parallel (force free) current, $`𝐉_{}`$ and the displacement current orthogonal to the $`𝐤𝐁`$ plane, $`\delta E_y`$.

### 4.1 Fast Wave Characteristics

The fast wave characteristics can be used to propagate finite changes in the displacement current orthogonal to the $`𝐤𝐁`$ plane. Maxwell’s equations can be used to find the corresponding induced magnetic field. By contrast, no information on the force-free current can be propagated along the fast characteristics. If one were to cover a region of spacetime with a mesh of short wavelength fast wave characteristics (i.e., the boundary of the spacetime region is a essentially a pure fast wave emitter and absorber), and time evolve the system using the method of characteristics then one would find that the current and charge in the spacetime would remain unchanged for all time, in all frames.

There is a claim in Blandford (2001, 2002) that since the fast wave can transport a toroidal magnetic perturbation, “it carries information about the poloidal current.” This conjecture is inconsistent with the nature of the fast wave. The fast wave, like the light wave carries only displacement current. The time dependent version of Ampere’s law does not require a coupling between physical current and the magnetic field.

### 4.2 Alfven Wave Characteristics

By contrast, changes in the parallel current can be propagated along the Alfven wave characteristics. However, no changes in in the displacement current orthogonal to the $`𝐤𝐁`$ plane propagate along these characteristics. From equation (3.27) this current is essentially a pure charge moving along the magnetic field lines at the speed of light in the proper frame. The Alfven wave can thought of as a charge discontinuity or perturbation propagating at the speed of light along the magnetic field (and its associated parallel current) and the resulting changes in the electromagnetic field can be found from Maxwell’s equations.

Note that in the degenerate case of wave propagation parallel to the magnetic field, there is no $`𝐤𝐁`$ plane. In the degenerate case, (3.24) implies that the Alfven wave carries no physical current, only a displacement current, $`\delta E_z`$.

#### 4.2.1 The Force Free Plasma-Filled Waveguide

A very important application of the nature of the Alfven wave characteristics is a semi-infinite cylindrical waveguide threaded by a uniform axial magnetic field and filled with a force free plasma. A rotating conductive disk is suddenly attached to one end of the plasma-filled waveguide (see Chapter 2 of Punsly (2001) for pictures and complete details). The conductive disk behaves like a unipolar inductor or Faraday wheel. When the Faraday wheel is connected, a Poynting flux begins to be radiated down the plasma-filled waveguide. The Poynting flux is supported by a force free current. In the language of force free waves, the Faraday wheel has launched an axisymmetric force free Alfven step wave (see the Appendix B for details).

Note that by (3.26d), the Alfven wave needs to be oblique if it supports a field aligned current (i.e., a component of the propagation vector that is orthogonal to the field, $`k_{}`$, is required). In the plasma-filled waveguide, the continuity of the tangential electric field at the surface of the rotating disk and the degeneracy condition in the force free plasma imply that there is a poloidal electric field, $`E_{}`$, that is launched by the Faraday wheel orthogonal to the poloidal magnetic field, $`B^P`$,

$`E_{}={\displaystyle \frac{\mathrm{\Omega }__DB^P\rho }{c}}.`$ (4-1)
The coordinate, $`\rho `$, is the radial coordinate in a cylindrical coordinate system and the field line angular velocity, $`\mathrm{\Omega }__F`$ is identical to the disk angular velocity, $`\mathrm{\Omega }__D`$. The radial variation of the electric field in (4.1) and in the oscillatory modes of (B1) is tantamount to a standing wave component of the propagation vector that is orthogonal to the magnetic field (note that the same wave polarization phenomenon occurs in vacuum electromagnetic waveguides). The step Alfven wave in the cylinder (constructed in Appendix B) is a superposition (Fourier integral) of the oblique, axisymmetric Alfven normal modes of the waveguide. The standing wave properties (obliquity) of the axisymmetric normal modes in the cylindrical waveguide are a consequence of the boundary conditions on the cylinder walls.

#### 4.2.2 Pulsar Magnetospheres

One would expect a strong analogy between the plasma-filled waveguide and an idealized, force free, field aligned pulsar magnetosphere. The neutron star is the unipolar inductor, the magnetosphere is the plasma-filled waveguide and the wind is the Poynting flux in this analogy. Thus, one would expect that the pulsar wind is essentially an Alfven wave radiated by the neutron star in this idealized limit.

Unfortunately, the causal nature of MHD in a relativistic magnetosphere was obfuscated by anecdotal comments in Blandford (2001) through an analysis in the force free approximation. If the force free analysis is not pathological then it must agree with the magnetically dominated limit of MHD. There is a claim that for axisymmetric modes, the Alfven wave can only support variations in the poloidal magnetic field, not toroidal variations. However, this cannot be correct as the example of the plasma-filled waveguide shows by the explicit construction of the axisymmetric Alfven modes that transport toroidal magnetic fields. The concept that a standing wave component of the propagation vector that is orthogonal to the magnetic field can exist in an axisymmetric mode (as in the plasma-filled waveguide) was not understood in Blandford (2001). The existence of such a wave vector polarization implies, with the aid of (3.30) that axisymmetric Alfven modes that are capable of transporting toroidal field perturbations along the direction of poloidal magnetic field most certainly exist in a pulsar magnetosphere. This oversight allows for the claim that the force free analysis leads to an “opposite conclusion” to that of the MHD analysis. In reality, if one is prudent, it can be shown as it has be done here that the two analyses agree, as they must.

## 5 Applications to Black Hole Magnetospheres

The force free approximation has been used in models of black hole magnetospheres in which spin energy is extracted electromagnetically from the black hole (Blandford and Znajek, 1977; Thorne et al, 1986). The force free wind systems in these treatments are characterized by two constants in each magnetic flux tube: the poloidal current (equivalently the toroidal magnetic field or angular momentum flux) and the field line angular velocity (this is quantity is equivalent to specifying the ratio of energy flux to angular momentum flux, the Goldreich-Julian charge density and the cross field potential).

The mathematical method of solution implements the force free condition at the horizon as a boundary condition that is used to determine the field aligned current density and the field line angular velocity in each flux tube. Since, by (3.17), the fast wave can propagate in any direction, at the speed of light, outgoing fast waves can propagate from anywhere outside of the event horizon (the fast critical surface for outgoing waves is the event horizon). By contrast, the ”light cylinder” located within the ergosphere is the Alfven critical surface, which no outgoing Alfven wave can traverse (Blandford, 2002; Punsly, 2001). For rapidly rotating black holes, this surface is sufficiently far from the event horizon that the constraints imposed by the spacetime metric near the horizon are irrelevant. Thus, the use of the horizon boundary condition can only be justified from a causal perspective if the fast mode is involved in determining the wind constants. However, by (3.19), changes in the poloidal current or charge can not be transmitted along the fast wave characteristics in the black hole magnetosphere. Thus, the spacetime near event horizon can not affect changes in the wind parameters. The use of the event horizon as a causal wind boundary is not permitted in the force free approximation. Consequently, the mathematical method of solution is flawed. The same causal structure is known to exist in magnetically dominated perfect MHD (see Chapter 6 of Punsly (2001)). Going to the force free limit was conjectured in (Blandford, 2001, 2002) to obviate any need for concern with this method of solution that were revealed in the full MHD analysis. This claim is clearly not supported by the results of this article.

The causal structure of the black hole magnetosphere in the force free limit has been shown, in this article to behave formally as the zero mass limit of MHD. As in perfect MHD, the wind parameters must be created outside of (farther out from the black hole than) the ergospheric Alfven critical surface.

A further critique of the MHD analysis of black hole causality was that the Alfven wave really does not transport changes in current and charge, only oscillations (Blandford, 2001). However, this statement is not an accurate depiction of the nature of propagating disturbances in MHD. It is well known that small amplitude MHD step waves obey the same dispersion relations as oscillatory waves. Alfven step waves and Alfven wave packets in general obey the same constitutive relations amongst field and dynamical quantities and propagate along the same characteristics as the oscillatory waves (Kantrowitz and Petschek, 1966). This also follows by taking the appropriate Fourier integral to create a step wave from sinusoidal waves. The Alfven wave packets can transport finite changes in the charge density and field aligned currents in a force free magnetosphere. This is shown by explicit construction in Appendix B.

## 6 Discussion

There are three main reasons for not considering black hole magnetospheres that are force free everywhere to be physically realizable:

1. As mentioned in the Introduction, the physical plasma state near the event horizon is always inertially dominated. This physical state is exactly the opposite of the zero inertia limit that is inherent to the force free assumption.
2. From equation (2.7) there is no region within the force free magnetosphere where energy flux (Poynting flux) and angular momentum flux can be created. Energy and angular momentum flux (essentially, the poloidal current and Goldreich-Julian charge density) must be injected into the wind system at the boundaries. If the magnetosphere is force free everywhere then the energy flux (current and charge) must emerge from the horizon. However, the charge and current necessary to support the energy and angular momentum flux can not be propagated along the plasma wave characteristics emanating from the spacetime near the horizon as discussed in the last section. Thus, such a wind system is inherently acausal.
3. There are three MHD computer simulations of magnetically dominated black hole magnetospheres (the magnetic energy density exceeds the plasma energy density in the initial state) that have been published. All of them evolve to a state with strong cross-field (inertial, not force free) current flow in the ergosphere (Camenzind and Khanna, 2000; Koide et al, 2002; Semenov et al, 2001). The inertial current flow initiates in regions far enough from the event horizon that it is not interpretable as a horizon surface current. The strong cross-field current is a necessary component for extracting black hole energy in the theory of black hole gravitohydromagnetics (GHM) as well in Punsly (2001).

It should be noted that a wind system beyond (farther from the black hole than) the ergospheric dynamo region for the toroidal magnetic field (poloidal current) could in principle be well described by the force free approximation. From the nature of the force free characteristics, this dynamo region can exist only in the spacetime farther from the black hole than the ergospheric light cylinder (Alfven critical surface). The Poynting flux can propagate outward from the dynamo region to asymptotic infinity as an essentially force free wind. Such a wind system is described in detail in Chapter 9 of Punsly (2001).

I would like to thank Ferd Coroniti for his valuable and intelligent insights. I am also indebted to Serguei Komissarov for his useful comments regarding force free plasmas and propagating discontiuities.

## Appendix A First Order Curvature Effects - Hybrid Modes

In this appendix, we calculate the first order effects of field line curvature on the force free small amplitude plasma modes in a proper frame. One finds new hybrid modes as in MHD. The first order corrections to the short wavelength approximation can be found in the eikonal approximation. The sinusoidal solutions described in (3.2) are replaced by:

$$J^\mu =J__0^\mu +\delta J^\mu e^{ı[{\scriptscriptstyle 𝐤𝑑𝐫}\omega t]},$$
(A1)
where $`𝐤`$ is a slowly varying vector valued function. In the present circumstance, the unperturbed current no longer vanishes by assumption. However, by the force free condition, it flows parallel to the magnetic field in the proper frame. This can be expressed in terms of a slowly varying “radius of curvature”, $`R`$, as

$$𝐉__0=cR^1𝐁__0.$$
(A2)
The validity of the eikonal approximation in (A1) requires that:

$`𝐤R1,`$ (A3)
$`\frac{k}{k^2}1.`$ (A4)
The perturbed Maxwell’s equations and the perturbed degeneracy condition carryover from the sinusoidal case, (3.4) - (3.11), unchanged. The only perturbation equation that changes is the force free condition, (3.2), that now becomes:

$$𝐁__0\times \delta 𝐉+\delta 𝐁\times 𝐉__0=0.$$
(A5)
The perturbed force free relation found from expanding the ”x” and ”z” components of (A5) is

$$\frac{4\pi ik^2c^2}{c^2k^2\omega ^2}(kR)^1\delta J_z=\delta J_y.$$
(A6)
Note that by (A3), the large radius of curvature approximation, equation (A6) is approximately (3.12). The ”y” component of the perturbed force free relation (A5) combined with (A6) yields the relation:

$$\frac{16\pi ^2k^4c^4(kR)^2\left[c^2k^2\omega ^2\right]^2}{\left[c^2k^2\omega ^2\right]^2}\frac{k_{}}{k}B__0\delta J_z+\frac{k_{}}{k}B__0\delta J_x=0.$$
(A7)
Note that in the limit of a large magnetic field radius of curvature that (A7) becomes (3.13). Thus, the complete set of perturbed equations reduce to the equations found in section 3 as the curvature of the field gets small.

The dispersion relation is found by combining (A7) with the perturbed degeneracy condition (3.14),

$$\left[\omega ^2c^2k^2\right]\left[\omega ^2c^2k_{}^2\right]=16\pi ^2k^2k_{}^2c^4(kR)^2.$$
(A8)
There are two roots of the dispersion relation.:

$`\omega ^2=c^2k^2+{\displaystyle \frac{1}{2}}c^2\left[\sqrt{k_{}^4+64\pi ^2k_{}^2k^2(kR)^2}k_{}^2\right],`$ (A9)
$`\omega ^2=c^2k_{}^2{\displaystyle \frac{1}{2}}c^2\left[\sqrt{k_{}^4+64\pi ^2k_{}^2k^2(kR)^2}k_{}^2\right].`$ (A10)
Differentiating (A9) and (A10) yields the group velocities of the two modes. For parallel propagation, the group velocity is the speed of light up to terms that are $`O\left[(kR)^2\right]`$. The parallel propagation assumption is invalidated at second order since the magnetic field is curved. In the limit that $`(k_{}R)^21`$, one can obtain the following Taylor series expansions of the group velocity for oblique propagation:

$`\frac{d\omega }{dk}\pm \frac{c}{\sqrt{1+16\pi ^2\mathrm{cot}^2\theta \left(kR\right)^2}}\pm c\left[18\pi ^2(kR)^2\mathrm{cot}^2\theta \right],`$ (A11)
$`\frac{d\omega }{dk}\pm \frac{c\mathrm{cos}\theta }{\sqrt{116\pi ^2\mathrm{cot}^2\theta \left(kR\right)^2}}\pm c\mathrm{cos}\theta \left[1+8\pi ^2(kR)^2\mathrm{csc}^2\theta \right].`$ (A12)
Equations (A9) and (A11) approximate to (3.16) and (3.17) as the radius of curvature becomes infinite. Thus, these equations represent the fast hybrid mode. Similarly, (A10) and (A12) become (3.22) and (3.23) in the short wavelength limit. These equations represent the hybrid Alfven wave. Notice that the fast hybrid mode propagates slower than the speed of light for oblique waves propagating across a curved background magnetic field. From (A12), the Alfven mode does not propagate perpendicular to the field. Thus, there is not Alfvenic properties for the fast mode to couple into this case. Thus, the hybrid fast mode still propagates at the speed of light for perpendicular propagation as indicated in (A11).

There are implications of this result to the primary focus of this article. Since the fast wave no longer propagates isotropically, it can now carry a physical current (not just a displacement current). Denote the displacement current as, $`J__D^y`$. Then (3.11a), (A6) and (A9) imply that the ratio of displacement current to physical current for oblique fast hybrid waves is:

$$\frac{J__D^y}{J_{}}\frac{(kR)^3\mathrm{tan}^4\theta \mathrm{sin}\theta }{64\pi }.$$
(A13)
Consequently, the ratio of displacement current to the physical current diverges in the short wavelength limit as in (3.19). Similarly, from the perturbed Gauss’ law (3.7a), (3.11c) and (A9), one can find the ratio of displacement current to the charge density:

$$\frac{J__D^y}{\rho _e}c\frac{(kR)^3\mathrm{tan}^4\theta \mathrm{sin}\theta }{64\pi }.$$
(A14)
Comparing (A9) to (A13) and (A14), note that as the fast wave propagation speed slows down due to field line curvature, it carries more physical current and charge. It becomes increasingly more like the Alfven wave as the field line curvature increases. However, (A13) and (A14) show that in the eikonal approximation, the fast wave is still almost purely electrodynamic: a propagating displacement current with a minute physical current that is only $`O\left[(kR)^3\right]`$ of the displacement current.

## Appendix B The Alfven Wave Packet in the Cylindrical Waveguide

In this appendix it is shown by explicit construction that a superposition of oscillatory Alfven waves (an Alfven wave packet) does indeed transport field aligned currents and charge, not just oscillatory perturbations. We rely on the magnetically dominated MHD treatment of the plasma-filled cylindrical waveguide found in Chapter 2 of Punsly (2001) and pass to the force-free limit. In order to save space, we do not reproduce the results of Chapter 2 here, but refer to them in detail and the reader is urged to consult that reference. In summary, a cylindrical waveguide filled with a tenuous plasma is placed inside of an infinitely long solenoid so that it is threaded by a strong axial magnetic field. The waveguide is semi-infinite in the sense that one end is terminated by a rotating conductive wheel. This rotating conductor or Faraday wheel (also known as a homopolar generator or unipolar inductor) is an active element in the waveguide circuit that radiates plasma waves. In the first subsection, we construct the relevant axisymmetric, m=0, perfect MHD Alfven modes in the cylinder using the detailed results of Punsly (2001). In the next subsection, we compare and contrast a wave packet (a step wave) created from a linear superposition of the previously calculated oscillatory Alfven modes to the radiation emanating from the Faraday wheel (the nature of this radiation was derived previously in terms of a formal MHD wind solution in the waveguide in Punsly (2001)). Finally, it is noted that the depiction of the wind as a pure Alfven step wave becomes increasingly more accurate as one passes to the force-free limit.

### B.1 Components of the Oscillatory Waves

This section constructs an oscillatory solution from Chapter 2 of Punsly (2001) for an axisymmetric Alfven wave, $`m=0`$, in the plasma-filled waveguide that is consistent with the attachment of the Faraday wheel (in the realm of magnetically dominated perfect MHD).

#### B.1.1 The Electric Field

From equation (2.93b) of Punsly (2001), the radial electric field in an axisymmetric Alfven wave in a cylindrical waveguide is arbitrary and is determined by the boundaries (the field lines can shear relative to each other, hence these modes are often called “shear Alfven waves”). Considering the frozen-in condition within the Faraday wheel and the continuity of the tangential electric field at the Faraday wheel/plasma interface let us choose the radial electric field in the Alfven wave to be (in cylindrical coordinates):

$$E^\rho =\frac{\mathrm{\Omega }__D\rho }{c}B__0e^{i(k_zz\omega t)}.$$
(B1)
The other electric field components are given in Punsly (2001) by (2.86c) and (2.93a) as

$`E^\varphi =0,`$ (B2)
$`E^z=0.`$ (B3)
#### B.1.2 The Currents

From equations (B1)-(B3) and Gauss’ law, the charge density is;

$$\rho _e=\frac{\mathrm{\Omega }__DB__0}{2\pi c}e^{i(k_zz\omega t)}.$$
(B4)
From the Alfven wave dispersion relation, (2.82) of Punsly (2001), and the frozen-in form of the momentum equation, (2.86a) of Punsly (2001), along with equation(A1), the cross-field current density is

$$J^\rho =i\frac{v__Ik_z\mathrm{\Omega }__D\rho }{4\pi U__A^2c}B__0e^{i(k_zz\omega t)},$$
(B5)
where $`v__I`$ is the intermediate three speed and $`U__A`$ is the pure Alfven speed. In the force free limit, the pure Alfven speed diverges as the intermediate speed approaches the speed of light. From the law of current conservation,

$$\frac{\rho _e}{t}+𝐉=0,$$
(B6)
combined with (B4) and (B5) yields the axial current density:

$$J^z=\frac{\rho _e}{v__I}c^2=\frac{\mathrm{\Omega }__DB__0c}{2\pi v__I}e^{i(k_zz\omega t)}.$$
(B7)
#### B.1.3 The Magnetic Field

Using the value of the axial current in (B7) in Ampere’s law yields the toroidal magnetic field strength:

$$B^\varphi =\frac{B__0\mathrm{\Omega }__D\rho }{v__I}e^{i(k_zz\omega t)}.$$
(B8)
The other components are given by (2.92a) and (2.96) of Punsly (2001):

$`B^z=0,`$ (B9)
$`B^\rho =0.`$ (B10)
### B.2 Making Wave Packets

We can construct wave packets of the oscillatory solutions above by taking Fourier integrals. We are interested in the step wave. For the Alfven wave, one can use the dispersion relation (2.82) of Punsly (2001) in the cylinder to write,

$$e^{i(k_zz\omega t)}=e^{ik_z(v__Itz)}.$$
(B11)
The step function is given by,

$$\mathrm{\Theta }(v__Itz)=\frac{1}{2\pi i}_{\mathrm{}}^+\mathrm{}\frac{e^{ik_z(v__Itz)}dk_z}{k_z+iϵ}.$$
(B12)
We construct a wave packet of Alfven waves, $`\mathrm{\Psi }`$, with a spectral amplitude, $`A(k_z)`$,

$`\mathrm{\Psi }={\displaystyle _{\mathrm{}}^+\mathrm{}}A(k_z)\psi (k_z)𝑑k_z,`$ (B13)
$`A(k_z)={\displaystyle \frac{1}{(2\pi i)(k_z+iϵ)}},`$ (B14)
where $`\psi (k_z)`$ is one of the oscillatory wave function components in equations (B1)-(B10) and $`ϵ`$ is an arbitrarily small positive number in the usual sense. The wave packet of oscillatory solutions with field and current components given by (B1)-(B10) and spectral amplitude given by (B14) has the following field and current distributions (for the sake of demonstrating equivalence, the corresponding equation number from Punsly (2001) that were derived by the MHD wind calculation in the cylinder is placed on the right hand side of the equals sign in these relations):

$`E^\rho ={\displaystyle \frac{\mathrm{\Omega }__D\rho }{c}}B__0\mathrm{\Theta }(v__Itz)=(2.98),`$ (B15)
$`E^\varphi =0,`$ (B16)
$`E^z=0,`$ (B17)
$`\rho _e={\displaystyle \frac{\mathrm{\Omega }__DB__0}{2\pi c}}\mathrm{\Theta }(v__Itz)=(2.99),`$ (B18)
$`J^\rho ={\displaystyle \frac{v__I\mathrm{\Omega }__D\rho }{4\pi U__A^2c}}B__0\delta (v__Itz),`$ (B19)
$`J^z={\displaystyle \frac{\rho _e}{v__I}}c^2={\displaystyle \frac{\mathrm{\Omega }__DB__0c}{2\pi v__I}}\mathrm{\Theta }(v__Itz)=(2.119),`$ (B20)
$`B^\varphi ={\displaystyle \frac{B__0\mathrm{\Omega }__D\rho }{v__I}}\mathrm{\Theta }(v__Itz)=(2.118),`$ (B21)
$`B^z=0,`$ (B22)
$`B^\rho =0.`$ (B23)
These quantities agree with the downstream state of the plasma that was found near the flow front in Punsly (2001) in the limit that $`v_z=v__I`$. Notice that $`J^\rho `$ vanishes downstream of the wavefront and only has a surface component on the wavefront as it does in equation (2.110) of Punsly (2001).

### B.3 Conclusion

From equations (B18) and (B20), the Alfven wave packet actually transports charge and field aligned current.

The solution is not exact. There are errors associated with the small inertial terms in the magnetically dominated limit. At the wavefront, the relativistic MHD shock equations do not solve exactly. These equations are the conservation of the axial components of the stress energy tensor across the discontinuity in the frame of the wavefront. Pressure balance (magnetic pressure from the toroidal magnetic field in the frame of the wavefront) can not be achieved because the Alfven wave, unlike the fast wave, has no compressive properties. The errors in the stress-energy balance can be found in the frame of the propagating discontinuity by using (B15) and (B21) to compute the proper toroidal magnetic field. The errors in the shock equations are on the order of $`(\mathrm{\Omega }__D\rho U__A^1)^2`$. This is biquadratic in two small quantities: $`\mathrm{\Omega }__D\rho /c1`$ by construction and by the magnetically dominated condition, $`cU__A^11`$. The violation of the shock relations are therefore a very small second order effect in the magnetically dominated limit. Notice that the errors vanish completely in the force-free limit as the Alfven three speed approaches the speed of light. In magnetically dominated perfect MHD, these small errors are accounted for by a fast switch-on shock that creates a $`B^\rho `$ and an $`E^\varphi `$ downstream of the shock front. This is an MHD precursor, an infinitesimal distance upstream, to the Alfven rotational discontinuity at the terminus of the Alfven step wave. As demonstrated by explicit construction above, it is the Alfven rotational discontinuity that imprints the charge and field aligned current on the waveguide plasma. The interpretation of the waveguide wind solution as an Alfven wave is exact in the force-free limit and is extremely accurate to first order in the magnetically dominated MHD limit. However, in relativistic MHD it is important to note that both fast and Alfven modes are needed to solve this type of “Riemann problem.”

It was an oversight in Punsly (2001) that nonlinear fast waves such as an abrupt discontinuity can in principle transport charge and current. However, in practice, due to the limited polarization properties of such waves (the fast discontinuity can only affect changes in the magnetic field in the plane spanned by the normal to the wavefront and the magnetic field upstream) unless the magnetic field is perfectly homogenous (not the case in a magnetosphere) they must be accompanied by an Alfven rotational discontinuity in a magnetosphere in order for the discontinuities to be solutions of the relativistic shock equations (Kantrowitz and Petschek, 1966). The important point is that in general both modes are needed. The discussions of Punsly (2001) were designed to show that the Alfven wave is always required in MHD, so that the fast wave cannot determine the wind constants alone - and this conclusion is unchanged and hence the MHD causality arguments are unchanged. In the force-free limit, the situation is more extreme. Since the fast wave is linear in the force free limit (the fast speed is always c), all wave packets constructed from the complete set of oscillatory fast modes will also carry no current or charge in any frame (see section 3.1). Thus, the Alfven wave is responsible for setting the wind constants in a force free black hole magnetosphere.

This discussion corroborates the physical accuracy of the assumptions used in the one dimensional string analysis of flux tubes in Semenov et al (2001). They actually ignore the fast mode completely (because as they say, the fast wave is well known to primarily provide only pressure changes in an MHD flux tube) and only consider the Alfven mode. This seems like a very reasonable assumption based on the analysis provided here.


text:
# Acknowledgments

## Acknowledgments

We would like to thank G. Bonelli and P. Salomonson for discussions. This work is partly supported by EU contract HPRN-CT-2000-00122.


text:
# Properties of the Spin-flip Amplitude of Hadron Elastic Scattering and Possible Polarization Effects at RHIC

## I Introduction

Several attempts to extract the spin-flip amplitude from the experimental data show that the ratio of spin-flip to spin-nonflip amplitudes can be non-negligible and may be only slightly dependent on energy akch ; sel-pl .

For the definition of new effects at small angles and especially in the region of the diffraction minimum one must know the effects of the Coulomb-hadron interference with sufficiently high accuracy. The Coulomb-hadron phase was calculated in the entire diffraction domain taking into account the form factors of the nucleons prd-sum . Some polarization effects connected with the Coulomb hadron interference, including some possible odderon contribution, were also calculated z00 .

The model-dependent analysis based on all the existing experimental data of the spin-correlation parameters above $`p_L6GeV`$ allows us to determine the structure of the hadron spin-flip amplitude at high energies and to predict its behavior at superhigh energies yaf-wak . This analysis shows that the ratios $`Re\varphi _5^h(s,t)/(\sqrt{|t|}Re\varphi _1^h(s,t))`$ and $`Im\varphi _5^h(s,t)/(\sqrt{|t|}Im\varphi _1^h(s,t))`$ depend on $`s`$ and $`t`$ (see Fig.1 a,b). At small momentum transfers, it was found that the slope of the “residual” spin-flip amplitudes is approximately twice the slope of the spin-non flip amplitude. The obtained spin-flip amplitude leads to the additional contribution to the pure CNI effect at small $`t`$ (Fig. 1 c).

The dependence of the hadron spin-flip amplitude on $`t`$ at small angles is closely related with the basic structure of hadrons at large distances. We show that the slope of the so-called “reduced” hadron spin-flip amplitude (the hadron spin-flip amplitude without the kinematic factor $`\sqrt{|t|}`$) can be larger than the slope of the hadron spin-non-flip amplitude, as was observed long ago predaz .

## II The slope of the hadron amplitudes

For an exponential form of the amplitudes this coincides with the usual slope of the differential cross sections divided by $`2`$. At small $`t`$ ($`0÷0.1GeV^2`$), practically all semiphenomenological analyses assume: $`B_1^+B_2^+B_1^{}B_2^{}`$. If the potentials $`V_{++}`$ and $`V_+`$ are assumed to have a Gaussian form in the first Born approximation $`\varphi _1^h`$ and $`\widehat{\varphi _h}^5`$ will have the same form $`\varphi _1^h(s,t)exp(B\mathrm{\Delta }^2)`$, $`\varphi _5^h(s,t)=qBexp(B\mathrm{\Delta }^2)`$. In this special case, therefore, the slopes of the spin-flip and “residual”spin-non-flip amplitudes are indeed the same. A Gaussian form of the potential is adequate to represent the central part of the hadronic interaction. The form cuts off the Bessel function and the contributions at large distances. If, however, the potential (or the corresponding eikonal) has a long tail (exponential or power) in the impact parameter, the Bessel functions can not be taken in the approximation form and the full integration leads to different results.

If we take $`\chi _i(b,s)He^{a\rho },`$ we obtain

$`F_{nf}(s,t)=a/[(a^2+q^2)^{3/2}]1/[a\sqrt{a^2+q^2}]exp(Bq^2)`$ (1)
with $`B=1/a^2`$. For the “residual” spin-flip amplitude, on the other hand, we obtain tur1

$`\sqrt{|t|}\stackrel{~}{F_{sf}}(s,t)=(3aq)/[(a^2+q^2)^{5/2}](3aqB^2)/(\sqrt{a^2+q^2})exp(2Bq^2).`$ (2)
In this case, therefore, the slope of the “residual” spin-flip amplitude exceeds the slope of the spin-non-flip amplitudes by a factor of two. A similar behaviour can be obtained with the standard dipole form factor tur1 .

## III The determination of the structure of the hadron spin-flip amplitude

Note that if the “reduced” spin-flip amplitude is not small, the impact of a large $`B^{}`$ will reflect in the behavior of the differential cross section at small angles sel-sl . The method gives only the absolute value of the coefficient of the spin-flip amplitude. The imaginary and real parts of the spin-flip amplitude can be found only from the measurements of the spin correlation coefficient.

Let us take the spin nonflip amplitude in the standard exponential form with definite parameters: slope $`B^+`$, $`\sigma _{tot}`$ and $`\rho ^+`$. For the “residual” spin-flip amplitude, on the other hand, we consider two possibilities: equal slopes $`B^{}=B^+`$ and $`B^{}=2B^+`$. The results of these two different calculations are shown in Fig.2. It is clear that around the maximum of the Coulomb-hadron interference, the difference between the two variants is very small. But when $`|t|>0.01GeV^2`$, this difference grows. So, if we try to find the contribution of the pomeron spin-flip, we should take into account this effect. As the value of $`A_N`$ depends on the determination of the beam polarization, let us calculate the derivative of $`A_N`$ with respect to $`t`$, for example, at $`\sqrt{s}=500GeV`$.

If we know the parameters of the hadron spin non-flip amplitude, the measurement of the analyzing power at small transfer momenta helps us to find the structure of the hadron spin-flip amplitude. There is a specific point of the differential cross sections and of $`A_N`$ on the axis of the momentum transfer, - $`t_{re}`$, where $`|ReF_c^{++}|=|ReF_h^{++}|`$. This point $`t_{re}`$ can be found from the measurement of the differential cross sections addapr . At high energies at the point $`t_{re}`$ tur1  we obtain for $`pp`$-scattering

$`ReF_{sf}^h(s,t)={\displaystyle \frac{1}{2(ImF_{nf}^h(s,t)+ImF_{nf}^c(t))}}A_N(s,t){\displaystyle \frac{d\sigma }{dt}}ReF_{sf}^c(t).`$ (3)
We can again take the hadron spin-nonflip and spin-flip amplitudes with definite parameters and calculate the magnitude of $`A_N`$ by the usual complete form while the real part of the hadron spin-flip amplitude is given by (3). Our calculation by this formula and the input real part of the spin-flip amplitude are shown in Fig. 2 c. At the point $`t_{re}`$ both curves coincide. So if we obtain from the accurate measurement of the differential cross sections the value of $`t_{re}`$, we can find from $`A_N`$ the value of the real part of the hadron spin-flip amplitude at the same point of momentum transfer.

## IV The model predictions

The model zpc  takes into account the contribution of the hadron interaction at large distances and leads to the high-energy spin-flip amplitude. The model gives the large spin effects in the $`hh`$-elastic scattering and predicts non-small effects for the $`PP2PP`$ experiment at RHIC especially in the diffraction dip domain yaf-str . The additional pure CNI effects can be calculated using the Coulomb-nuclear phase prd-sum . These polarization effects will be present at RHIC energy, even though $`F_h^+0`$ at high energy. Our model calculations show on Fig.3 for both cases.

The model gives the standard $`t`$-dependence of $`ReF_h++`$ and $`ImF_h++`$. Instead of it, in a convenient parameterization of both the modulus and the phase one can obtain the alternative case, in which $`ImF_h(s,t)`$ has the zero at small $`t`$ (for details, see kun3 ). Such an approach enables one to specify the elastic hadron scattering amplitude $`F_h(s,t)`$ directly from the elastic scattering data. The difference between the phases leads either to central or peripheral distributions of elastic hadron scattering in the impact parameter space. The obtained form of $`A_N^{CN}`$ at small momentum transfers differs for the two variants beginning at $`|t|>0.05`$ GeV2 (Fig.4 a). The difference reaches $`2\%`$ at $`t=0.15GeV^2`$ and, in principle, can be measured in an accurate experiment. Now let us calculate the Coulomb-hadron interference effect \- $`A_N^{CN}`$ in the two alternatives for higher $`|t|`$: (i) the diffraction dip is created by the ”zero” of the $`ImF_h(s,t)`$ part of the scattering amplitude and $`ReF_h(s,t)`$ fills it; (ii) the diffraction dip is created by the ”zero” of the $`ReF_h(s,t)`$ part of the scattering amplitude and $`ImF_h(s,t)`$ fills it. The results are shown in Fig. 4 (b) for $`\sqrt{s}=50`$ GeV and in Fig. 4 (c) at $`\sqrt{s}=500`$ GeV.


text:
# Chandra observations of the recurrent nova CI Aql after its April 2000 outburst

## 1 Introduction

Recurrent novae (RN) are a rare subclass of cataclysmic variables, showing repeated outbursts with typical recurrence times of 50–100 yrs. These outbursts are thought to be thermonuclear runaways on the surface of a white dwarf that has accreted enough material from the companion star to ignite the surface hydrogen. If the hot white dwarf could be observed directly, it would have the characteristics of a supersoft X-ray source (SSS); i.e., its luminosity would lie between $`10^{36}`$ ergs/s and $`10^{38}`$ ergs/s, with $`kT`$ on the order of tens of eV (Nomoto 1982, Fujimoto 1982, Iben 1982). The short recurrence times require either or both small envelope masses and high accretion rates of the order of $``$10-8 M/yr (Warner 1995). Since the ejecta of RN typically are not metal-enriched, the white dwarfs in RN are not significantly eroded but may grow in mass (Starrfield et al. 1985). The fitting of RN optical light curves suggests that the white dwarf masses are high, near the Chandrasekhar limit (e.g. Hachisu & Kato 1999, Hachisu et al. 2000), suggesting that RN may be progenitors of type Ia supernovae.

Though only 8 other recurrent novae are known, three subclasses are defined (Warner 1995) according to the nature of the companion and, consequently, the length of orbital period: the T Pyx subclass with dwarf companions ($`P23`$ hrs), the U Sco subclass with evolved main-sequence or subgiant companions ($`P1`$ day), and the RS Oph subclass with red giant companions ($`P`$ $`>`$ 100 days).

CI Aql was discovered when it showed a $`\mathrm{\Delta }`$m = 4.6 mag outburst in 1917 (Reinmuth 1925). It was classified as a dwarf nova with long cycle length (Duerbeck 1987), but this was cast in doubt (Szkody and Howell 1992, Greiner et al. 1996) until a second outburst, discovered in April 2000 (Takamizawa 2000), revealed CI Aql to be a RN. A subsequent search in the Harvard Plate collection led to the discovery of another outburst in 1941/1942 (Schaefer 2001b). CI Aql is an eclipsing binary system with an orbital period of 14.84 hrs (Mennickent & Honeycutt 1995) and therefore belongs to the U Sco subclass.

Recurrent novae are expected to emit soft X-rays during a short interval after the ejected shell has become optically thin and before the hydrogen shell-burning ceases and during an interval of cooling. During this time, the effective radius should be somewhat larger than the radius of the white dwarf, and the luminosity should therefore be larger than 1037 erg/s. Thus it should fit the profile of a SSS. The only previous recurrent novae observed with X-ray detectors were T Pyx and U Sco. T Pyx has had its last outburst in 1967, but was argued to be a quasi-persistent supersoft X-ray source because optical data suggested a luminosity of $``$1036 erg/s (Patterson et al. 1998). T Pyx was not detected, however, in a ROSAT PSPC observation carried out in December 1998 (Greiner, unpublished). A BeppoSAX observation carried out 20 days after the U Sco outburst revealed a luminous ($``$1036 erg/s) soft X-ray component with a best-fit temperature of $``$9$`\times `$105 K, consistent with steady nuclear burning continuing for at least $``$1 month after the nova outburst (Kahabka et al. 1999). Further observations of the evolution of the soft and hard X-ray emission in other sources are crucial to better test theoretical models and also to test the conjecture that RN are SN Ia progenitors.

Based on the observed optical light curve of CI Aql over the first year after the 2000 outburst, Hachisu and Kato (2001) predicted that CI Aql should have been active at soft X-ray wavelengths between December 2000 until August 2001. This motivated us to observe CI Aql using Chandra.

## 2 Observations and Results

We obtained two Chandra observations using ACIS-S. The first one was performed on June 1, 2001 for 2.15 ksec, and the second one on August 1, 2001 for a total of 19.88 ksec (see Tab. 1 for details).

During the first Chandra observation, we detected CI Aql with a count rate of (6.0$`\pm `$2.3)$`\times `$10-3 cts/s, collecting a total of 15 photons. In the August observation, the mean count rate is (8.6$`\pm `$0.7)$`\times `$10-3 cts/s, thus CI Aql exhibited the same X-ray intensity as on 1 June 2000,

No single-component model fits the X-ray spectrum of CI Aql (second observation; Figs. 1–3). As a first approximation, and since the spectrum looks thermal, we have used a blackbody spectrum plus a power law, and find a reasonably good fit with $`kT`$=50$`\pm `$20 eV and a power law photon index of 1.78$`\pm `$0.72. Note, though, that due to the small number of photons the fit parameters are correlated, and therefore the reduced $`\chi _{\mathrm{red}}^2`$ is smaller than one. Using a white dwarf atmosphere model (van Teeseling et al. 1994) gives a similar value of $`kT`$=38$`\pm `$15 eV. The best-fit absorbing column is 3.6$`\times `$1021 cm-2 (for the blackbody model) or 5.5$`\times `$1021 cm-2 (for the white dwarf atmosphere model), consistent with the foreground extinction of E(B-V) = 0.85 (Kiss et al. 2001), but considerably less than the total galactic column in this BII=–0$`.^{}`$8 direction of 1.5$`\times `$1022 cm-2 (Dickey & Lockman 1990). The unabsorbed luminosity of the blackbody (white dwarf atmosphere) component is 1.0$`\times `$1033 (D/1 kpc)2 erg/s (8.0$`\times `$1032 (D/1 kpc)2 erg/s).

Physically, such a high-temperature optically thick emission at such a low luminosity is difficult to explain. The corresponding radius of the emitting area would be 40 km, a minute fraction of the white dwarf radius. We therefore turned to other spectral models.

In some recent X-ray grating observations of novae, the spectra show a wealth of narrow emission lines in the 0.2–2 keV range. Though most of the identifications of these lines and their origin are not settled, we considered the possibility that pure line emission at soft energies could model the spectrum of CI Aql: To model the main emission around 0.5 keV we need only one Gaussian line with a width corresponding to the energy resolution of ACIS-S. The best-fit line energy is 0.358 keV. The spectrum at high energies (above $``$4) can similarly be fit by models other than a power law. A Raymond-Smith model gives a best-fit temperature of $`kT`$ = 2.7$`\pm `$2.4 keV a thermal bremsstrahlung model gives $`kT`$ = 3.1$`\pm `$3.0 keV. Of course, the spectrum could also be explained by a superposition of several lines, but we do not have enough photons and/or spectral resolution to resolve this issue.

The first-epoch spectrum can be fit by the model parameters of the fit to the second-epoch observation. This is not a great surprise given the few photons in the first-epoch spectrum. Yet, we can conclude that (1) the intensity was seemingly steady, and (2) no major spectral change (change in the emission mechanism) occurred between the two Chandra observations both in the soft as well as in the hard component – this also excludes any dramatic decrease of absorption which is expected once the expelled shell gets optically thin for X-rays, or cooling of the soft component (though this cannot be directly proven with the 15 detected X-ray photons).

The second Chandra observation covers about 40% of the orbital period of CI Aql. We therefore analyzed the temporal behaviour of the X-ray emission, but did not find orbital variability.

## 3 Discussion

After our Chandra observations occurred, more optical data became available (Schaefer 2001a), leading Hachisu & Kato (2002a) to revise their modelling of the optical decay light curve of CI Aql. The important new result of their revision for the interpretation of our Chandra data is that the supersoft X-ray phase should have lasted only until May 2001, thus ending even before the first of our Chandra observations. But the source would have remained hot and several orders of magnitude more luminous than observed had nuclear burning ended as recently as May 2001.

The faintness of the detected X-ray emission and the fit results of a blackbody (or white dwarf atmosphere) model indeed suggest that CI Aql was not in the nova plateau phase during our Chandra observations. Instead, the luminosity derived for an optically thick emission model explaining the soft part of the X-ray spectrum is only of the order of 1033 (D/1 kpc)2 erg/s.

However, even beyond the plateau phase the system should have be observable as a SSS, because the cooling of the hot white dwarf and the concordant luminosity decay is expected to last months to years. Model calculations (Hachisu & Kato 2002b) show that the temperature should have declined by a factor of $``$2 (from 55 eV to 30 eV) between our two Chandra observations, whereas the Chandra data do not allow a temperature change by more than 10 eV. These model calculations also suggest that in the June–August 2001 time frame the luminosity of the cooling white dwarf would still be in the 10$`{}_{}{}^{37}10^{36}`$ erg/s range, at least a factor of 1000 more than the observed X-ray luminosity. Thus, we find it very unlikely that the soft X-ray component around 0.5 keV is an optically thick component from the (cooling) white dwarf. This implies that nuclear burning ended earlier than indicated by even the revised prediction.

Alternatively, the soft component could be composed of (a) non-resolved emission line(s). After fitting a Gaussian line to the soft component, one realizes that the spectral resolution of $``$100 eV makes this Gaussian line broad enough to fit perfectly the soft component, and not even a second Gaussian line is necessary to explain the X-ray emission below 0.7 keV. The best-fit line energy for this Gaussian is 0.358$`\pm `$0.058 keV, and the unabsorbed line luminosity 2.7$`\times `$1032 (D/1 kpc)2 erg/s. We suggest this line to be the H-like C VI line (nominal energy at 0.3676 keV).

Rebinning the spectrum into only 1$`\sigma `$ confidence level bins (Fig. 3), the search for possible further excess emission above the broad band continuum emission at higher energies reveals marginal evidence for possible lines at 0.8 keV, 1.28 keV and 6.3 keV. These could be interpreted as O VII, Ne IX or Mg X and Fe XXVI lines. Given the signal-to-noise ratio of the spectrum, this is speculative.

The hard spectrum is not constrained: it could be fit by either a power law, a Raymond-Smith model or a bremsstrahlung model. Also the marginal Fe detection does not provide a clue for its origin (e.g. reflection off the disk, or thermal plasma). Indeed, the Fe emission could be due to matter from the white dwarf, consistent with the strong FeII lines seen in the optical spectra (Kiss et al. 2001). The small temperature of order of 3 keV, if interpreted as thermal emission, is consistent with CI Aql having had a very slow nova outburst, and a rather low expansion velocity, a correlation which has been found already in earlier novae (e.g. Mukai & Ishida 2001 for Nova V382 Vel).

We note, on the other hand, that the unabsorbed luminosity of this hard component of $``$1030 (D/1 kpc)2 erg/s (0.5–10 keV) is in the range of X-ray luminosities usually found in cataclysmic binaries during quiescence. However, CI Aql during quiescence is definitely fainter: ROSAT observations in 1992 and 1993 failed to detect CI Aql, implying an upper limit on the unabsorbed X-ray luminosity of 5$`\times `$1029 erg/s (D/1 kpc)2 in the 0.1–2.4 keV range (Greiner et al. 1996). Thus, 16 months after its outburst, CI Aql had not yet returned to its quiescent level of X-ray emission.

Thus, it seems clear that we observed the CI Aql system when it had ended both the phase of shell burning and cooling as well as the phase in which the wind was obscuring the soft emission. In the cases of Nova V382 Velorum (Burwitz et al. 2002) and V1494 Aql (Krautter et al. 2002) it has recently been possible to observe the dramatic changes of the X-ray emission of novae from a highly absorbed spectrum during the very early times, over the development and subsequent fading of a luminous supersoft X-ray component towards a completely line dominated spectrum with hardly any detectable continuum. While the first two stages of this evolution are physically well understood in terms of the expanding shell blocking less and less of the X-ray radiation of the white dwarf until burning on its surface ceases, the origin of the line dominated spectrum is not clear. One possibility is that during the late phase of a nova, when the shell has become transparent, the medium around the system is ionized by the still UV-EUV bright white dwarf. It would require far less than one part in a thousand conversion efficiency of photoionization of the C VI line to explain, e.g., the spectrum of CI Aql as observed with Chandra, consistent with photoionization models around supersoft ionizers (Rappaport et al. 1994). However, because of the steep and rapidly changing density profile due to the expanding shell/wind, the recombination emission may be rapidly suppressed. Alternatively, the line dominated spectrum could result from the interaction of the expanding shell with the circumstellar matter, shocks within the expanding shell, or from shocks due to collisions of a fast wind with the interstellar matter. Detailed investigations of the complex emission line spectra have just begun, and promise a better understanding in the near future.

We are indebted to H. Tananbaum for granting DDT time for the Chandra observations. We thank I. Hachizu and M. Kato for prompt information on their modelling of CI Aql. This material is based upon work supported by the National Science Foundation under Grant No. 9815655 and by a NASA LTSA award, NAG5-10705.


text:
# Universal scaling behavior of directed percolation and the pair contact process in an external field

## 1 Introduction

In this paper we consider the scaling behavior of direction percolation (DP) and of the pair contact process (PCP) numerically in $`1+1`$ dimensions and analytically within a mean-field approximation (see for a review ). Both models display a continuous phase transition from an active to an inactive state. In contrast to the unique absorbing state of DP the PCP is characterized by infinitely many absorbing states. Despite this difference it was observed in numerical simulations that the critical exponents describing the steady state scaling behavior of the PCP agree well with those of DP . This result is in agreement with a renormalization group analysis of a phenomenological field theory approach of the PCP which reduces in the steady state to the corresponding DP equation .

In this work we consider the steady state scaling behavior of DP and the PCP in an external field which is conjugated to the order parameter. In particular we examine the scaling behavior of the order parameter (equation of state) as well as of the order parameter fluctuations. In contrast to previous works on the PCP we do not restrict our attention to the critical exponents but consider the so-called universal scaling functions too. In this way we are able to show that both models are characterized by the same universal scaling behavior, i.e., the steady state scaling behavior of the PCP belongs to the universality class of directed percolation.

In the next section we briefly review the scaling behavior of DP and define the scaling functions. We consider then in section 3 for the first time the PCP in a conjugated field and show that the scaling behavior is characterized by the DP critical exponents. The universal scaling behavior of DP and the PCP is presented in section 4. The obtained universal equation of state is compared to the results of a two loop renormalization group approach of the corresponding Langevin equation . Finally we derive the equation of state of the PCP within a mean-field approximation. Therefore we consider the PCP with particle creation at randomly selected sites. This random neighbor interaction suppresses long range correlations and the model is analytically tractable.

## 2 Directed percolation

In order to examine the scaling behavior of the $`1+1`$-dimensional DP universality class we consider the directed site percolation process using the Domany-Kinzel automaton . It is defined on a diagonal square lattice and evolves by parallel update according to the following rules. A site at time $`t`$ is occupied with probability $`p_2`$ ($`p_1`$) if both (only one) backward sites (time $`t1`$) are occupied. If both backward sites are empty a spontaneous particle creation takes place with probability $`p_0`$. Directed site percolation corresponds to the choice $`p_1=p_2=p`$ and $`p_0=0`$. The propagation probability $`p`$ is the control parameter of the phase transition, i.e., below a critical value $`p_\mathrm{c}`$ the activity ceases and the system is trapped forever in the absorbing state (empty lattice). On the other hand a non-zero density of (active) particles $`\rho _\mathrm{a}`$ is found for $`p>p_\mathrm{c}`$. The best estimate of the critical value of directed site percolation is $`p_\mathrm{c}=0.705489(4)`$  and we use this value throughout this work.

The density of active sites $`\rho _\mathrm{a}`$ is interpreted as the order parameter of the absorbing phase transition and vanishes at the critical point according to

$$\rho _\mathrm{a}\delta p^\beta ,$$
(1)
with $`\delta p=(pp_\mathrm{c})/p_\mathrm{c}`$. Furthermore the order parameter fluctuations $`\mathrm{\Delta }\rho _\mathrm{a}=L(\rho _\mathrm{a}^2\rho _\mathrm{a}^2)`$ diverge as

$$\mathrm{\Delta }\rho _\mathrm{a}\delta p^\gamma ^{}.$$
(2)
The fluctuation exponent $`\gamma ^{}`$ obeys the scaling relation $`\gamma ^{}=\nu _{}2\beta `$, where $`\nu _{}`$ describes the divergence of the spatial correlation length at the critical point. The critical behavior of the order parameter is shown in figure 1. The data are obtained from numerical simulations of systems with periodic boundary conditions. Considering various system sizes $`L131072`$ we take care that our results are not affected by finite-size effects.

Similar to equilibrium phase transitions it is possible in DP to apply an external field $`h`$ that is conjugated to the order parameter. Being a conjugated field its has to destroy the absorbing phase and the corresponding linear response function has to diverge at the critical point i.e.,

$$\chi =\frac{\rho _\mathrm{a}}{h}\mathrm{}.$$
(3)
In DP the external field is implemented (see for instance ) as a spontaneous creation of particles ($`p_0=h>0`$). Clearly the spontaneous creation of particles destroys the absorbing state and therefore the absorbing phase transition at all. Figure 1 shows how the external field results in a smoothening of the zero-field order parameter curve. The inset displays that the fluctuations are peaked for finite fields. Approaching the transition point ($`h0`$) this peak becomes a divergence signalling the critical point.

Close to the transition point the order parameter and its fluctuations obey the following scaling ansatzes

$$\rho _\mathrm{a}(\delta p,h)\lambda \stackrel{~}{r}_{\mathrm{sDP}}(\delta p\lambda ^{1/\beta },h\lambda ^{\sigma /\beta }),$$
(4)
$$\mathrm{\Delta }\rho _\mathrm{a}(\delta p,h)\lambda ^\gamma ^{}\stackrel{~}{d}_{\mathrm{sDP}}(\delta p\lambda ,h\lambda ^\sigma ),$$
(5)
with the scaling functions $`\stackrel{~}{r}_{\mathrm{sDP}}`$ and $`\stackrel{~}{d}_{\mathrm{sDP}}`$ where the index $`\mathrm{sDP}`$ denotes site directed percolation. The scaling function of the order parameter corresponds to the equation of state and we recover equation (1) by setting $`\delta \rho \lambda ^{1/\beta }=1`$ at $`h=0`$. On the other hand setting $`h\lambda ^{\sigma /\beta }=1`$ leads to

$$\rho _\mathrm{a}(\delta p=0,h)h^{\beta /\sigma }$$
(6)
at $`\delta \rho =0`$. Analogous we get the scaling behavior of the fluctuations at $`p_\mathrm{c}`$

$$\mathrm{\Delta }\rho _\mathrm{a}(\delta p=0,h)h^{\gamma ^{}/\sigma }$$
(7)
as well as equation (2). The field dependence of the order parameter and of its fluctuations at $`p_\mathrm{c}`$ is plotted in figure 2.

The above scaling forms imply that curves corresponding to different values of the external field collapse to a single one if $`\rho _\mathrm{a}h^{\beta /\sigma }`$ and $`\mathrm{\Delta }\rho _\mathrm{a}h^{\gamma ^{}/\sigma }`$ is considered as a function of the rescaled control parameter $`\delta ph^{1/\sigma }`$. Using the best known estimates of the critical exponents $`\beta =0.276486`$, $`\sigma =2.554216`$, and $`\gamma ^{}=0.543882`$  we get beautiful data collapses of our numerical data (see figure 3).

It is worth mentioning that the validity of the scaling behavior of the equation of state (4) implies the required singularity of the linear response function (3). Using the fact that the susceptibility is defined as the derivative of the order parameter with respect to the conjugated field we get

$`\chi (\delta p,h)`$ $`=`$ $`{\displaystyle \frac{}{h}}\rho _\mathrm{a}(\delta p,h)`$ (8)
$``$ $`\lambda ^{1\sigma /\beta }\stackrel{~}{c}_{\mathrm{sDP}}(\delta p\lambda ^{1/\beta },h\lambda ^{\sigma /\beta }).`$ (9)
Thus we get for the critical behavior

$$\chi (\delta p,h)=\frac{}{h}\rho _\mathrm{a}(\delta p,h)|_{h0}\delta p^\gamma $$
(10)
and

$$\chi (\delta p,h)=\frac{}{h}\rho _\mathrm{a}(\delta p,h)|_{\delta p0}h^{\gamma /\sigma },$$
(11)
respectively. Here the susceptibility exponent $`\gamma `$ is given by the scaling relation

$$\gamma =\sigma \beta $$
(12)
which corresponds to the well known Widom equation of equilibrium phase transitions. Notice that in contrast to the scaling behavior of equilibrium phase transitions the non-equilibrium absorbing phase transition of DP is characterized by $`\gamma \gamma ^{}`$.

Again, the validity of the scaling ansatz of the order parameter (4) can be used to check that the implemented field satisfies the condition (3), i.e., to check whether the field is conjugated to the order parameter or not. We apply this procedure in the following section where we examine the scaling behavior of the PCP.

## 3 Pair contact process

The PCP as introduced by Jensen  is one of the simplest models with infinitely many absorbing states showing a continuous phase transition. At time $`t`$ sites on a lattice of length $`L`$ with periodic boundary conditions can either be occupied ($`n_i(t)=1`$) or empty ($`n_i(t)=0`$). Pairs of adjacent occupied sites $`i,i+1`$, linked by an active bond, annihilate each other with rate $`p`$ or create an offspring with rate $`1p`$ at either site $`i1`$ or $`i+2`$ provided the target site is empty. The density of active bonds $`\rho _\mathrm{a}`$ is the order parameter of a continuous phase transition from an active state for $`p<p_\mathrm{c}`$ to an inactive absorbing state without particle pairs. The behavior of the PCP order parameter and its fluctuations are plotted in figure 4. The data are obtained from simulations on various system sizes $`L131072`$ with periodic boundary conditions. Our analysis reveals that the critical value is $`p_\mathrm{c}=0.077093(3)`$ which agrees with the value $`p_\mathrm{c}=0.077090(5)`$  obtained from a finite-size scaling analysis of the lifetime distribution.

In contrast to DP there is no unique absorbing state (empty lattice) but infinitely many, as any configuration with only isolated inactive particles is absorbing. Thus in the thermodynamic limit the system will be trapped in one of an infinite number of absorbing configurations for $`p>p_\mathrm{c}`$. Despite the different structure of the absorbing states the steady state scaling behavior of the PCP is believed to be characterized by the DP critical exponents $`\beta `$, $`\gamma ^{}`$, $`\gamma `$ etc. On the other hand the dynamical scaling behavior, associated with activity spreading of a localized seed depends on the details of the system preparation .

Recently, Dickman et al.  considered the PCP with an external field that randomly creates isolated particles. Thus the external field couples to the particle density but not to the order parameter itself, i.e., the external field is not conjugated to the order parameter. The authors observe that the external field shifts the critical values $`p_\mathrm{c}`$ continuously and that the critical exponents are unaffected by the presence of the particle source. In order to investigate the PCP in a conjugated field the implementation of the external field of  has to be modified. Several modifications of the external field are possible. For instance in absorbing phase transitions with particle conservation  the conjugated field triggers movements of inactive particles which can be activated in this way .

As shown below spontaneous particle creation with rate $`h`$ acts as a conjugated field analogous to DP. Figure 4 shows how the spontaneous particle creation smoothens the critical zero field curves similar to the DP behavior (see figure 1). We simulated the PCP at the critical value $`p_\mathrm{c}`$ for various fields. The order parameter and its fluctuations as a function of the external field $`h`$ are presented in figure 5. Approaching the transition point, $`\rho _\mathrm{a}`$ and $`\mathrm{\Delta }\rho _\mathrm{a}`$ scale according to the equations (6,7) where the exponents $`\beta /\sigma `$ and $`\gamma ^{}/\sigma `$ agree with the DP values.

Furthermore we assume that the order parameter and the order parameter fluctuations obey analogous to DP the scaling forms

$$\rho _\mathrm{a}(\delta p,h)\lambda \stackrel{~}{r}_{\mathrm{PCP}}(\delta p\lambda ^{1/\beta },h\lambda ^{\sigma /\beta }),$$
(13)
$$\mathrm{\Delta }\rho _\mathrm{a}(\delta p,h)\lambda ^\gamma ^{}\stackrel{~}{d}_{\mathrm{PCP}}(\delta p\lambda ,h\lambda ^\sigma ),$$
(14)
where the distance to the critical point is now given by $`\delta p=(p_cp)/p_\mathrm{c}`$. Using the DP values of the critical exponents $`\beta `$, $`\sigma `$ and $`\gamma ^{}`$ we get convincing data collapses (see figure 6). As pointed out above, the validity of the scaling ansatz (13) implies the singular behavior of the linear response function (3), i.e., the spontaneous particle creation in the PCP can be interpreted as an external field conjugated to the order parameter. Furthermore the data collapse confirms again that the steady state scaling behavior of the PCP is characterized by the DP exponents.

## 4 Universal scaling behavior

One of the most striking features of critical phenomena is the concept of universality, i.e., close to the critical point the scaling behavior depends only on few fundamental parameters whereas the interaction details of the systems do not alter the scaling behavior. In the case of systems with short range interactions these parameters are the symmetry of the order parameter and the dimensionality of space $`D`$ . Classical examples of such universal behavior are for instance the coexistence curve of liquid-vapor systems  and the equation of state in ferromagnetic systems . In the case of absorbing phase transitions the corresponding universality hypothesis states that systems exhibiting a continuous phase transition to a unique absorbing state generally belong to the universality class of directed percolation . A different scaling behavior is observed only in systems with additional symmetries (such as parity conservation or particle conservation) or in systems with quenched randomness.

Following the concept of universality two models belong to the same universality class if the critical exponents and the universal scaling functions are identical. The universal scaling functions $`\stackrel{~}{R}_{\mathrm{DP}}`$ and $`\stackrel{~}{D}_{\mathrm{DP}}`$ of the DP universality class can be easily determined by introducing non-universal metric factors $`c_i`$ and $`d_i`$ (see e.g. ) for each scaling argument in the scaling functions (4, 5, 13, 14), i.e.,

$$\rho _\mathrm{a}(\delta p,h)\lambda \stackrel{~}{R}_{\mathrm{DP}}(c_1\delta p\lambda ^{1/\beta },c_2h\lambda ^{\sigma /\beta }),$$
(15)
$$\mathrm{\Delta }\rho _\mathrm{a}(\delta p,h)\lambda ^\gamma ^{}\stackrel{~}{D}_{\mathrm{DP}}(d_1\delta p\lambda ,d_2h\lambda ^\sigma ).$$
(16)
Using the normalizations $`\stackrel{~}{R}_{\mathrm{DP}}(1,0)=\stackrel{~}{R}_{\mathrm{DP}}(0,1)=1`$ and $`\stackrel{~}{D}_{\mathrm{DP}}(1,0)=\stackrel{~}{D}_{\mathrm{DP}}(0,1)=1`$ the non-universal metric factors can be determined from the amplitudes of

$$\rho _\mathrm{a}(\delta p,h=0)(c_1\delta p)^\beta ,\rho _\mathrm{a}(\delta p=0,h)(c_2h)^{\beta /\sigma }$$
(17)
and

$$\mathrm{\Delta }\rho _\mathrm{a}(\delta p,h=0)(d_1\delta p)^\gamma ^{},\mathrm{\Delta }\rho _\mathrm{a}(\delta p=0,h)(d_2h)^{\gamma ^{}/\sigma }.$$
(18)
Like the value of the critical point $`p_\mathrm{c}`$ the non-universal metric factors may depend on the details of the considered system, e.g. the lattice structure, the boundary conditions, the used update scheme, etc. In the case of directed site percolation in the Domany-Kinzel automaton we have obtained the values $`c_1^{\mathrm{sDP}}=2.45`$, $`c_2^{\mathrm{sDP}}=0.112`$, $`d_1^{\mathrm{sDP}}=71.6`$ and $`d_2^{\mathrm{sDP}}=3984`$. For the PCP on a square lattice we have determined the values $`c_1^{\mathrm{PCP}}=0.665`$, $`c_2^{\mathrm{PCP}}=0.181`$, $`d_1^{\mathrm{PCP}}=3.21`$ and $`d_2^{\mathrm{PCP}}=62.11`$.

Analogous to the previous scaling analysis we set $`c_2h\lambda ^{\sigma /\beta }=1`$ and consider for both models the rescaled order parameter $`\rho _\mathrm{a}(c_2h)^{\beta /\sigma }`$ as a function of the rescaled control parameter $`c_1\delta \rho (c_2h)^{1/\sigma }`$ as well as the rescaled order parameter fluctuations $`\mathrm{\Delta }\rho _\mathrm{a}(d_2h)^{\gamma ^{}/\sigma }`$ as a function of $`d_1\delta \rho (d_2h)^{1/\sigma }`$, respectively. The corresponding data are presented in figure 7 and figure 8. In both cases we get a perfect data collapse of the curves showing that the one-dimensional PCP steady state scaling behavior belongs to the universality class of directed percolation.

Additionally to the universal scaling function $`\stackrel{~}{R}_{\mathrm{DP}}(x,1)`$ we plot in figure 7 the corresponding curve of an $`ϵ`$-expansion obtained from a renormalization group analysis of a Langevin equation (see equation (25) next section). Using the parametric representation  of the absorbing phase transition Janssen et al. showed recently that the equation of state is given by the remarkable simple scaling function

$$H(\theta )=\theta (2\theta )+\mathrm{O}(ϵ^3),$$
(19)
where $`ϵ`$ denotes the distance to the upper critical dimension $`D_\mathrm{c}=4`$, i.e., $`ϵ=D_\mathrm{c}D`$. Here the scaling behavior of the quantities $`\rho _\mathrm{a}`$, $`\delta p`$, and $`h`$ is transformed to the variables $`R`$ and $`\theta `$ through the relations

$$b\delta p=R(1\theta ),\rho _\mathrm{a}=R^\beta \frac{\theta }{2}.$$
(20)
The equation of state is given by

$$ah=\left(\frac{R^\beta }{2}\right)^\delta H(\theta )$$
(21)
with the metric factors $`a`$ and $`b`$. The whole phase diagram is described by the parameter range $`R0`$ and $`\theta [0,2]`$. Using the equations (19,20,21) we calculated the corresponding universal function and compare it in figure 7 with the numerically obtained universal scaling function $`\stackrel{~}{R}_{\mathrm{DP}}(x,1)`$. As can be seen the significant difference indicates that the $`\mathrm{O}(ϵ^3)`$ corrections to (19) are relevant, i.e., higher orders than $`\mathrm{O}(ϵ^2)`$ are necessary to describe the scaling behavior of directed percolation. Of course we expect that this difference will decrease with increasing dimension, i.e., for $`ϵ0`$.

Before closing this section we briefly mention that the dynamical scaling behavior of the PCP belongs to the DP universality class too if the spreading of a localized seed is considered at the so-called natural particle density . Examining spreading activity one usually considers the survival probability $`P_\mathrm{a}`$ of the activity as well as how the number of active particles $`N_\mathrm{a}=L\rho _\mathrm{a}`$ increases in time. In the case of DP the simulations are started with a single seed on an empty lattice. For the PCP an absorbing state at $`p_\mathrm{c}`$ is prepared to which a particle is added in order to create one seed (one active pair). At the critical point the following power-law behaviors are expected

$$N_\mathrm{a}t^\theta ,P_\mathrm{a}t^\delta .$$
(22)
Finite systems sizes limit these power-law behaviors and $`P_\mathrm{a}`$ and $`N_\mathrm{a}`$ obey the finite-size scaling ansatzes

$$N_\mathrm{a}(\delta p=0,L)\lambda \stackrel{~}{n}(\lambda ^{1/\theta }t,\lambda ^{1/\theta z}L)$$
(23)
$$P_\mathrm{a}(\delta \rho =0,L)\lambda \stackrel{~}{p}(\lambda ^{1/\delta }t,\lambda ^{1/\delta z}L).$$
(24)
where $`z`$ denotes the dynamical exponent. Analogous to the above analysis the universal scaling curves $`\stackrel{~}{N}`$ and $`\stackrel{~}{P}`$ are obtained by introducing appropriate non-universal metric factors. Using the values $`z=1.580745`$, $`\theta =0.313686`$, and $`\delta =0.159464`$  we get convincing data collapses and the universal scaling functions are plotted in figure 9.

## 5 Mean field scaling behavior

The mean-field equation of state of DP can be easily derived from the corresponding Langevin equation

$$_t\rho _\mathrm{a}=\delta p\rho _\mathrm{a}\lambda \rho _\mathrm{a}^2+\kappa h+D^2\rho _\mathrm{a}+\eta $$
(25)
which describes the order parameter field $`\rho _\mathrm{a}(\underset{¯}{x},t)`$ on a mesoscopic scale (see  for a detailed discussion). Here $`D`$ denotes the diffusion constant, $`\eta `$ denotes a multiplicative noise term with the correlator

$$\eta (\underset{¯}{x},t)\eta (\underset{¯}{x}^{},t^{})=\mathrm{\Gamma }\rho _\mathrm{a}(\underset{¯}{x},t)\delta ^d(\underset{¯}{x}\underset{¯}{x}^{})\delta (tt^{})$$
(26)
and $`\lambda >0`$, $`\kappa >0`$, and $`\mathrm{\Gamma }>0`$ are certain coupling constants. Neglecting spatial correlations and fluctuations ($`D=0`$ and $`\eta =0`$) one gets for the steady state behavior ($`_t\rho _\mathrm{a}=0`$)

$$\delta p\rho _\mathrm{a}\lambda \rho _\mathrm{a}^2+\kappa h=\mathrm{\hspace{0.33em}0}$$
(27)
from which it is easy to derive the universal scaling function

$$\rho _\mathrm{a}=\stackrel{~}{R}_{\mathrm{DP}}(c_1^{\mathrm{DP}}\delta p,c_2^{\mathrm{DP}}h)=\frac{c_1^{\mathrm{DP}}\delta p}{2}+\sqrt{c_2^{\mathrm{DP}}h+\left(\frac{c_1^{\mathrm{DP}}\delta p}{2}\right)^2},$$
(28)
where the non-universal metric factors are given by $`c_1^{\mathrm{DP}}=1/\lambda `$ and $`c_2^{\mathrm{DP}}=\kappa /\lambda `$, respectively. For zero-field we get the solutions $`\rho _\mathrm{a}=0`$ (absorbing state) and $`\rho _\mathrm{a}=c_1^{\mathrm{DP}}\delta p`$, i.e., the mean-field value of the order parameter exponent is $`\beta =1`$. On the other hand $`\delta p=0`$ leads to $`\rho _\mathrm{a}=(c_2^{\mathrm{DP}}h)^{1/2}`$ implying $`\sigma =2`$. As recently shown the universal scaling function (28) describes not only the scaling behavior of the DP universality class. It also occurs in the different universality class of absorbing phase transitions with a conserved field .

Let us now consider the following modification of the PCP. An active bond produces an offspring with rate $`(1p)`$ at an empty site selected at random. The rules for annihilation and action of the external field remain unchanged. This random neighbor interaction suppresses long range correlations and the model is therefore expected to be characterized by the mean-field scaling behavior. We denote the density of inactive bonds between an occupied and an empty site as $`\rho _\mathrm{i}`$. Bonds between empty sites are denoted as $`\rho _\mathrm{e}`$. Normalization requires $`\rho _\mathrm{e}+\rho _\mathrm{i}+\rho _\mathrm{a}=1`$.

Depending on the sites adjacent to the target site the number of active bonds $`n_\mathrm{a}`$, inactive bonds $`n_\mathrm{i}`$ or empty bonds $`n_\mathrm{e}`$ is changed. For instance if the adjacent sites are empty, for which the probability in absence of correlations is $`\rho _\mathrm{e}^2`$, the number of empty bonds decreases by two ($`\mathrm{\Delta }n_\mathrm{e}=2`$). On the other hand there are two new inactive bonds ($`\mathrm{\Delta }n_\mathrm{i}=+2`$). The total probability for this event is $`(1p)\rho _\mathrm{a}\rho _\mathrm{e}^2`$. A list of all possible processes and their mean field probabilities is given in table 1. Thus we obtain rate equations for the expectation values $`E[\mathrm{\Delta }n_x]`$ of the changes in active, inactive and empty bond numbers. These expectation values are zero in the steady state, i.e.,

$$E[\mathrm{\Delta }n_x]=\underset{\mathrm{\Delta }n_x}{}\mathrm{\Delta }n_xp(\mathrm{\Delta }n_x)=\mathrm{\hspace{0.33em}0}$$
(29)
with $`x\{a,i,e\}`$. In the case of $`E[\mathrm{\Delta }n_a]`$ we get

$`E[\mathrm{\Delta }n_a]`$ $`=`$ $`3p\rho _\mathrm{a}^34p\rho _\mathrm{a}^2\rho _\mathrm{i}p\rho _\mathrm{a}\rho _\mathrm{i}^2+2(1p)\rho _\mathrm{a}\rho _\mathrm{i}^2`$
$`+2(1p)\rho _\mathrm{a}\rho _\mathrm{i}\rho _\mathrm{e}+2h\rho _\mathrm{i}^2+2h\rho _\mathrm{i}\rho _\mathrm{e}=\mathrm{\hspace{0.33em}0},`$
whereas we get for the for inactive bonds

$`E[\mathrm{\Delta }n_i]`$ $`=`$ $`2(1p)\rho _\mathrm{a}\rho _\mathrm{i}^2+2(1p)\rho _\mathrm{a}\rho _\mathrm{e}^2+2p\rho _\mathrm{a}^3`$
$`2p\rho _\mathrm{a}\rho _\mathrm{i}^22h\rho _\mathrm{i}^2+2h\rho _\mathrm{e}^2=\mathrm{\hspace{0.33em}0}.`$
Using this equations together with the normalization allows for calculating the order parameter for zero-field ($`h=0`$) which yields the non-trivial ($`\rho _\mathrm{a}>0`$) solution

$$\rho _\mathrm{a}=\frac{83p^25p2\sqrt{2p(1p)(3p4)^2}}{9p^29p+8}.$$
(32)
This solution is valid below the mean field critical point $`p_\mathrm{c}=8/9`$ whereas the trivial solution $`\rho _\mathrm{a}=0`$ is valid for all $`p`$ but unstable above $`p_\mathrm{c}`$. Expanding (32) around the critical point leads to

$$\rho _\mathrm{a}=\frac{3}{8}\delta p+\mathrm{O}(\delta p^2)$$
(33)
with $`\delta p=(p_\mathrm{c}p)/p_\mathrm{c}`$. Thus the mean field exponent of the PCP is $`\beta =1`$ and the non-universal metric factor $`c_1^{\mathrm{PCP}}=3/8`$.

In order to obtain the order parameter in presence of an external field $`h`$ equations (5,5) are solved for $`\rho _\mathrm{i}`$ which yields

$`4h+4\rho _\mathrm{a}4h\rho _\mathrm{a}4p\rho _\mathrm{a}4\rho _\mathrm{a}^2`$ (34)
$`+\left\{12p^2\rho _\mathrm{a}^2+(2h+2\rho _\mathrm{a}2h\rho _\mathrm{a}2p\rho _\mathrm{a}2\rho _\mathrm{a}^22p\rho _\mathrm{a}^2)^2\right\}^{1/2}`$
$`\{(2h2\rho _\mathrm{a}+2p\rho _\mathrm{a}+2\rho _\mathrm{a}^22p\rho _\mathrm{a}^2)^2+\mathrm{}`$
$`\mathrm{}+4p\rho _\mathrm{a}(h+\rho _\mathrm{a}2h\rho _\mathrm{a}p\rho _\mathrm{a}2\rho _\mathrm{a}^2+h\rho _\mathrm{a}^2+2p\rho _\mathrm{a}^2+\rho _\mathrm{a}^3)\}^{1/2}`$ $`=`$ $`0.`$
To obtain the field dependence of the order parameter a series expansion around $`h=0`$ at $`p_\mathrm{c}=8/9`$ is performed which results in leading order

$$\rho _\mathrm{a}=\sqrt{\frac{3}{8}h},$$
(35)
i.e., the mean field values of the PCP are given by $`c_2^{\mathrm{PCP}}=\sqrt{3/8}`$ and $`\sigma =2`$.

Finally we derive the mean-field universal scaling function $`\stackrel{~}{R}`$ of the PCP. Therefore we write (34) as a function of the reduced control parameter $`\delta p`$ and perform the limits $`\rho _\mathrm{a}0`$, $`\delta \rho 0`$ and $`h0`$ with the constraint that $`\rho _\mathrm{a}/\sqrt{h}`$ and $`\rho _\mathrm{a}/\delta p`$ are finite. Thus we remains in leading order with

$$\frac{3}{8}\delta p\rho _\mathrm{a}\rho _\mathrm{a}^2+\frac{3}{8}h=\mathrm{\hspace{0.33em}0}.$$
(36)
Solving this equation yields

$$\rho _\mathrm{a}=\frac{1}{2}\frac{3}{8}\delta p+\sqrt{\frac{3}{8}h+\left(\frac{1}{2}\frac{3}{8}\delta p\right)^2}=\stackrel{~}{R}_{\mathrm{DP}}(c_1^{\mathrm{PCP}}\delta p,c_2^{\mathrm{PCP}}h).$$
(37)
Thus we have shown that both models are characterized by the same critical exponents and scaling function, i.e., similar to the $`1+1`$-dimensional case the mean-field solution of the PCP belongs to the mean-field universality class of directed percolation.

We would like to thank H. Hinrichsen for usefull discussions. This work was financially supported by the Minerva Foundation (Max Planck Gesellschaft).