text:
# Magneto-transport properties of dilute granular ferromagnets

## Abstract

We present magnetoresistance (MR) measurements performed on quench condensed granular Ni thin films which are on the verge of electric continuity. In these systems the electric conductivity is believed to be governed by the resistance between a very small number of grains. The films exihibit sharp resistance jumps as a function of magnetic field. We interpret these findings as being the result of magneto-mechanical distortions that occour in single grains which act as bottlenecks in the dilute percolation network. The observed features provide a unique measure of magnetostriction effects in nano-grain structures as well as being able to shed light on some of the properties of regular granular magnetic films.

Transport in granular metals, i.e. systems of metallic islands imbedded in an insulating matrix, has been an active area of research for a number of decades but many experimental findings are still far from being understood. Among these are the temperature dependence of the resistance ($`R\alpha \mathrm{exp}[\frac{T_0}{T}]^{0.5})`$ and the large noise associated with granular films. The case in which the metallic islands are ferromagnetic is of particular interest since these systems show magnetoresistance (MR) properties which are reminiscent of the giant magnetoresistance phenomena characteristic of magnetic multilayers.

A natural way to treat granular metals is by considering a Miller-Abrahams resistor network in a similar manner to that of disordered systems ambegaokar . Each pair of grains is represented by a resistor with conductance proportional to the tunneling probability between the grains adkins . In this representation the sample conductance is governed by ”critical resistances” which act as bottlenecks for the conductance and dominate the electric transport of the percolation network. Since the density of these critical resistors is determined by the distribution of resistances one can expect that as the distribution becomes very large the distance between critical resistors (i.e. the correlation length of the electric percolation network) grows until eventually it should approach the sample size. In this limit the sample is critical and the behavior of any bottleneck grain will have large influence on the total conductance. In this paper we describe conductance measurements performed on granular Ni films in which the transport is believed to be governed by a very small number of grains. The samples exhibit unique sharp resistance jumps as a function of magnetic field. We attribute these jumps to buckling type magnetostriction deformation of bottleneck grains leading to abrupt changes in the sample conductivity.

An established method to fabricate systems of granular metals in general strongin ; granular bob ; granular goldman ; granular rich  and granular ferromagnets in particular aviad1 ; aviad2 ; philmag  in a very controlled way is quench condensation. In this method thin films are grown by sequential evaporation on a cryogenically cold substrate under UHV conditions while monitoring the film thickness and resistance. If the samples are quench condensed on a non-passivated substrate such as SiO, they begin growing in a granular manner so that the film is constructed of separated islands. As more material is quench condensed, the average distance between the islands decreases and the resistance drops. In these samples there is a critical thickness $`d_C`$, below which no conductivity can be practically measured. Once the thickness, $`d`$, of the sample is larger than $`d_C,`$ the sheet resistance, R, drops exponentially with $`d`$ reflecting the tunneling nature of the transport in these systems.

The quench condensation method provides a very sensitive control on the sample growth process, allowing one to terminate the evaporation at any desired stage of the material deposition and ”freeze” the morphological configuration. In particular, it allows one to stop the film growth at a thickness for which a measurable conductivity first appears across the sample ($`dd_C`$). The current study was performed on quench condensed granular Ni films which were prepared under these conditions. Recently, we have studied a number of quench condensed granular materials in-which deposition was terminated at the verge of electric conductivity future . These show a number of features which differ from the usual 2D macroscopic systems and are characteristic of mesoscopic samples . Among these are the following future :

* Resistance switches and induced noise appear randomly (as a function of time) in many of these samples. There are also large sample to sample variations in the transport properties.
* Superconducting films show magnetoresistance oscillations similar to those which are observed in 1D superconducting granular wires herzog . The area associated with these oscillations is (700Å)2.
* As in macroscopic granular films, these samples exhibit an electric field dependence of the form I$`e^{\frac{F_0}{F}^{\frac{1}{3}}}`$for large electric field, F. However, the effective sample length extracted from the measured I-V curves is three orders of magnitude smaller than those of macroscopic systems.

These findings lead us to consider the possibility that despite the fact that these samples have macroscopic dimensions ($``$0.5 cm2), the percolation network is dilute enough so that the transport is governed by a very small area of the sample. This observation is supported by our measurements of the magneto-transport of granular ferromagnets described below. We name these samples ”dilute samples” as opposed to the regular 2D films.

The magnetoresistance of 2D quench condensed Ni has been extensively studied in the past aviad1 ; aviad2 ; philmag . A typical MR trace is shown in figure 1. The field in these experiments was applied perpendicular to the film. The MR curve shows a resistance peak centered at H=0 followed by negative MR up to fields of $`\pm 0.5T`$. This behavior is similar to that observed in insulating granular ferromagnets prepared by co-evaporation of a magnetic material and an insulation material gittelman ; gerber ; yang ; honda ; sankar  and is a result of spin dependent tunneling between grains which have randomly oriented magnetic moments abeles . Applying a magnetic field aligns these moments causing a resistance decrease. The magnitude of the MR in our films is always of the order of 1-2% at T=4K. A unique feature that appears in all of the high-resistance quench condensed films is the increase of the resistance at high fields (H$`>`$0.5T) rather than a saturation that is expected for large H where all the grain moments are assumed to be aligned. This finding has been a puzzling issue for several years. As we note below, the results presented in this current work may shed some light on this high-field positive MR trend.

The MR curves of the dilute samples show several unique features that are not seen in the regular granular systems. These are depicted in figures 2 and 3 where it is seen that at low fields the samples exhibit a negative MR similar to those of regular 2D films. However, at specific magnetic fields the resistance switches abruptly to a higher value. In some cases we observe a single resistance jump, such as that seen in the top frame of figure 2, after which the resistance saturates at a high resistance value. In other cases we observe two jumps (figures 2 bottom frame and figure 3). The first is a relatively small resistance increase followed by the continuation of the negative MR trend and the second is a large resistance jump leading to saturation at a high value. It is worth noting that although the qualitative behavior at low fields is similar to that of usual 2D granular systems, the magnitude of the negative MR in our samples can be much larger. While the MR reported in previous works on granular ferromagnets rarely acceded 2%, we commonly observe effects of the order of 10% and in some cases the MR is as large as 60% (see figure 3). Clearly, in these films there is a large sample-to-sample variation in the magnitude of the negative MR and of the sharp jumps as well as the in the magnetic field at which the jump occurs. This feature further reflects the mesoscopic nature of these geometries.

A number of properties of the MR curves should be noticed. First, although the traces are very reproducible, the jumps are hysteretic as seen in figure 2. Sweeping the field in opposite directions produces mirror images of MR. Second, raising the temperature results in two clear trends (figure 3): The magnitude of the jumps decreases with raising the temperature in a similar manner to that of the negative MR at low fields, and, at the same time the magnetic field of the jumps shifts to higher values. Finally, perhaps the most surprising finding is the fact that the saturated resistance value at high magnetic fields is larger than that of the value at H=0. At zero field the magnetic grain moments are believed to be totally random and the resistance is expected to be maximal. The high resistance value at high magnetic fields is therefore a very unexpected result.

In an attempt to understand the above behavior we observe that the electric properties of these samples may be dominated by a small number grains which are characterized by the largest inter-grain distance that participate in the conduction paths. We have considered the possibility that the jumps are a result of magnetic moment flip of such a bottleneck grain. Recently, Meilikhov grainflip  considered the behavior of small single-domain grains in an external magnetic field and showed that their magnetic moment may flip abruptly if the initial angle between the magnetic field and the easy axis of the grain is larger than $`\pi /2.`$If a bottleneck grain is subject to such magnetic moment flip the MR of the entire sample conductance would be subject to sharp features.

However, the above model seems to be inconsistence with a number of our observations. In the first place this scenario could naturally explain sharp resistance decreases. Our samples always exhibit sharp increases in the resistance. Furthermore, the resistance value after the jump is higher than that at H=0 where the grain moment orientations are assumed to be completely random. Secondly, a pre-requirement for a resistance jump is that the bottleneck grain moment orientation have a larger-than-$`\pi /2`$-angle with the magnetic field direction. Our grains are believed to be in the superparamagnetic state, and it is hard to see how this condition is fulfilled. But even if such a situation occurs, computer simulations on such a granular configuration show that if the grain moment begins with such an angle, the initial MR would be positive rather than negative. In addition to the above, according to the Meilikhov picture the field at which the jump occurs should decrease with increasing temperature since the thermal energy assists in overcoming the anisotropy energy. This is opposite to the trend observed in figure 3.

We have also considered the possibility that the resistance jumps are associated with switches between different current trajectories leading to sudden resistance changes. Such a process also doesn’t seem to explain our observations because of the reproducibility of the data. The jumps are always upwards, they are always of the order of magnitude of the negative MR itself and the curves trace each-other when sweeping the field many times. Random switches between different trajectories are expected to cause either increase or decrease of the resistance with the same probability, to be random in size and magnetic field and to be non-reproducible. Indeed, a number of our samples showed much larger resistance switches which were random both in field and in sign. An example is shown in figure 4. It is seen that after each switch the magnetoresistance curve follows a different trace having different MR magnitudes (ranging between 0.5-13 %). In particular, one of the traces shows a resistance jump similar to the effects described above. It seams reasonable to assume that this sample (unlike that of figures 1,2 and 3) is characterized by a number of bottleneck trajectories and the current switches randomly between these paths leading to the different MR fingerprints.

The above considerations lead us to suggest that the effects seen in our samples are due to magneto-mechanical deformation of the grains. Magnetostriction (MS) in bulk Ni is well known and widely studied but not a lot is known about such effects in nano-structures. Geometry changes in small grains are extremely difficult to measure and, in granular films which are achieved by co-evaporation, the grains are imbedded in an insulating matrix and any MS effects depend strongly on the elasticity of this insulator. The quench condensed systems provide us with a unique opportunity to study MS effects on small grains since they are free standing particles. We envision that as the grains start aligning parallel to each other, all perpendicular to the substrate, and dipole-dipole interactions cause magnetic repulsion between the grains. At some magnetic field, a buckling instability may occur and a small MS distortion of the grains will take place causing a conductance decrease. Since the conductance is achieved by tunneling, any small distortion (even by fractions of an Å) in the geometry of the grains that belong to a conduction bottleneck will have a large effect on the conductivity. We estimate the magnetostatic interactions between two neigbouring parallelly aligned grains to be approximately 70K (about an order of magnitude larger than the thermal energy of the system at the relevant measurements). Our experience shows that increasing the temperature of quench condensed grains above 40K causes annealing and irreversible morphology changes in the sample. It seems reasonable, therefore, that magnetic forces in our samples could also have a considerable effect on the morphology of the grains. According to this model, if an external field aligning the grains is applied perpendicular to the substrate, two opposing effect influence the resistance. One is the fact that the tunneling probability is higher for aligned grains. This acts to lower the resistance causing a negative MR background. The second is a magnetostatic repulsion causing the inter-grain distance to increase slightly, thus decreasing the tunneling probability at the buckling instability fields.

The model described above can shed light on some of the features observed in the MR of usual 2D quench condensed granular samples. As noted above, one of the puzzling features in these samples is the increase of the MR at high magnetic fields at which the resistance is expected to have reached saturation (see figure 1). This anomaly is observed mainly when the resistance of the granular films is relatively high. Such behavior can be explained by magnetic repulsion between grains causing minute geometrical changes which decrease the inter-grain tunneling probabilities causing resistance increases. In the dilute systems, this manifests itself as sharp resistance jumps while in regular samples one observes a continuous mild increase in resistance with increasing field. Such an effect can be particularly expected in quench condensed samples in which the grains are free standing on a substrate separated by vacuum. In granular samples, which are achieved by co-evaporation of metal and insulating materials, geometrical modifications are less probable.

In summary, we have fabricated and measured dilute granular Ni films in which the transport is believed to be governed by a very small number of grains. These samples exhibit sharp resistance increases at discrete magnetic fields. In order to explain these findings we proposed a qualitative model which involves buckling of bottleneck grains due to magnetostatic repulsion. The unique geometry of our samples enables one to detect the properties of a single grain even if the sample has macroscopic dimensions. Our granular samples also provide a unique opportunity to study magneto-mechanical effects on nanometer sized magnetic structures. The results obtained on these extreme samples shed some light on the behavior of the ferromagnetic granular samples in general and they may also be relevant for the development of new types of magneto-machanically based single-grain devices.

We gratefully acknowledge illuminating discussions with R.C. Dynes, Z. Ovadyahu and M. Pollak. This research was supported by the Binational USA-Israel fund grant number and by the Israeli academy of science.


format: text

text:
# Incoherent 𝜌⁰ electroproduction off nuclei

## I Introduction

High energy meson electroproduction off complex nuclei offers a promising tool to study the physics of hadron formation. The nuclear environment makes it possible to investigate the properties of hadrons immediately after their creation. In addition one can vary the energy and virtuality of the exchanged photon to examine the phenomenon of color transparency .

In the present paper we study exclusive incoherent $`\rho ^0`$ electroproduction off 14N and 84Kr in the energy regime of the HERMES experiment . In this kinematic regime the coherence length of the photon’s $`\rho ^0`$ component is large and quantum mechanical interference effects play an important role. In  we presented a method to account for this so-called shadowing effect in photoproduction within a semi-classical BUU transport model which we now apply for electron induced reactions. We distinguish between the ’initial state interactions’ of the photon that give rise to shadowing and the final state interactions (FSI) of the produced particles. The coupled channel treatment of the FSI within the transport model allows for a broader spectrum of FSI than the usually used Glauber theory . In addition we account for Fermi motion and Pauli blocking during the reaction and take the finite lifetime of unstable reaction products into account.

The time that the products of an elementary photon-nucleon interaction need to evolve to physical hadrons is called formation time. Any change in the reaction probability of a hadron $`h`$ during its formation time affects the so-called nuclear transparency

$$T_A=\frac{\sigma _{\gamma ^{}AhX}}{A\sigma _{\gamma ^{}NhX}}$$
(1)
if the formation length, i.e. the length that $`h`$ travels during its formation time, exceeds the average internucleon distance ($`1.7`$ fm) inside the nucleus.

In the case of diffractive vector meson electroproduction ($`\gamma ^{}NVN`$) the virtual photon initially produces a colorless $`q\overline{q}`$-pair whose size is expected to decrease with increasing photon virtuality $`Q^2`$. As long as the $`q\overline{q}`$-pair is very small, i.e. in the early stage of its evolution into the physical vector meson, it mainly reacts via its color dipole moment. This leads to a cross section that is quadratic in the $`q\overline{q}`$-pair’s size. At large enough energies the $`q\overline{q}`$-pair is frozen in this small sized configuration over a distance that, because of time dilatation, can exceed the diameter of a nucleus. This effect leads to a large nuclear transparency for vector meson production at large photon energy $`\nu `$ and virtuality $`Q^2`$ (color transparency). For a detailed review of color transparency see e.g.  and references therein.

It is not clear whether the energy and the virtuality of the photon in the HERMES experiment is large enough to see an onset of color transparency. In the present paper we do not find evidence for a finite formation time of the produced vector meson. The 14N data for the transparency ratio as function of the coherence length is compatible with the assumption that the diffractively produced $`\rho ^0`$ starts interacting with a hadronic cross section right after its production. Kopeliovich et al.  point out that this might be an accidental consequence of the specific correlation between the $`Q^2`$ of the photon and the coherence length of the $`\rho ^0`$ in the HERMES data. In their work a light-cone QCD formalism was used to incorporate formation and coherence length effects in coherent and incoherent vector meson electroproduction. In addition the effect of gluon shadowing was studied, but found to be negligible in the case of incoherent $`\rho `$ production. In Ref.  the decrease of the transparency ratio due to the finiteness of the $`\rho ^0`$ lifetime was discussed as well. This effect is automatically included within our transport model. Our coupled channel calculation shows in addition that the effects of side-feeding in the FSI are unimportant because of the kinematic cuts of the HERMES experiment. Therefore Glauber theory seems to be appropriate for these reactions. However, we find that for heavy nuclei the kinematic cuts lead to an additional reduction of the transparency ratio due to elastic FSI which scatter the produced particle out of the acceptance window. One must therefore be careful which cross section one uses to describe the FSI in Glauber theory.

Our paper is structured in the following way: In Section II A we show how we describe the electron-nucleon interaction and how we account for shadowing within our model. In the limit of purely absorptive FSI our result for the incoherent vector meson production cross section turns out to be equivalent to the known Glauber result of Huefner et al. . However, in contrast to the result of Ref.  our expression allows for an intuitive physical interpretation. The transport model itself is sketched in Section II B. The results of our calculations for the exclusive incoherent $`\rho ^0`$ production cross section within the coupled channel treatment are presented in Section III in comparison with experimental data from the HERMES experiment. We close with a short summary and outlook in Section IV.

## II Model

### A Shadowing

In the one photon exchange picture the interaction of an electron and a nucleon (or nucleus) can be reduced to the reaction of a virtual photon on the hadronic target. The exchange of more than one photon is suppressed by further factors of the fine structure constant $`\alpha _{em}`$. We follow the method of Friberg and Sjostrand  and use the event generator PYTHIA v6.2  to describe the interaction of the (virtual) photon and a nucleon. The basic idea is that in a photon hadron collision the photon not necessarily interacts as a point particle (direct interaction) but might fluctuate into a vector meson $`V=\rho ^0,\omega ,\varphi ,J/\mathrm{\Psi }`$ (vector meson dominance) or perturbatively branch into a $`q\overline{q}`$ pair before the interaction (generalized vector meson dominance, GVMD). Within the PYTHIA model it turns out that the latter is very unlikely in the kinematic regime of the HERMES experiment (photon energy $`\nu 1020`$ GeV, $`Q^20.55`$ GeV2) as can be seen from Fig. 1 where we show the contribution of the different photon components to the total $`\gamma ^{}N`$ cross section at an invariant mass $`W=`$5 GeV. In most cases the photon will therefore scatter deep inelastically from a parton of the target nucleon or fluctuate into a vector meson before it reaches the nucleon. In the latter case the vector meson might either scatter diffractively from the nucleon or a hard scattering between the constituents of the vector meson and the nucleon might take place. The hard scattering leads to the excitation of one or two hadronic strings which finally fragment into hadrons. We assume that the hadrons that emerge from the edges of these strings, i.e. those that contain quarks not originating from the string fragmentation, can interact with a hadronic cross section right after the photon nucleon reaction. Consequently, for the diffractively produced vector meson we do not use a formation time since none of its constituents arise from a string fragmentation.

If the struck nucleon is embedded in a nucleus one has to account for its Fermi motion and binding energy as well as Pauli blocking of final state nucleons. In addition one has to be aware that the nuclear environment influences the VMD part of the photon nucleon interaction, since the vector meson components get modified on their way through the nuclear medium to the interaction point.

We express the physical photon state $`|\gamma `$ in terms of vector meson states $`|V`$ and a state $`|\gamma _0`$ which consists of the point-like photon and the GVMD part. This makes sense if one assumes that the GVMD part does not get shadowed or as in our case is unimportant for kinematic reasons:

$$|\gamma =\left(1\underset{V=\rho ,\omega ,\varphi ,J/\mathrm{\Psi }}{}\frac{e^2}{2g_V^2}F_V^2\right)|\gamma _0+\underset{V=\rho ,\omega ,\varphi ,J/\mathrm{\Psi }}{}\frac{e}{g_V}F_V|V.$$
(2)
The formfactor

$$F_V^2=\left(\frac{W^2}{Q^2+W^2}\right)^3\left[1+0.5\frac{4m_V^2Q^2}{(m_V^2+Q^2)^2}\right]\left(\frac{m_V^2}{m_V^2+Q^2}\right)^2$$
(3)
is taken from  and also accounts for contributions of longitudinal photons. In Eq. (3) $`\nu `$ denotes the energy of the photon, $`Q^2`$ its virtuality and $`m_V`$ the mass of the vector meson $`V`$. $`W`$ is the invariant mass of the photon nucleon system.

For a photon with momentum $`k\stackrel{}{e}_z`$ the wave function of the component $`V`$ can be written as

$$\psi _V^{(0)}(\stackrel{}{b},z)=\chi (\stackrel{}{b})e^{ikz}.$$
(4)
On the way through the nucleus the vector meson components get modified. Consider an ensemble of $`A`$ nucleons at positions $`\{\stackrel{}{r}_i\}=\{(\stackrel{}{s}_i,z_i)\}`$ which are labeled in such a way that $`z_1<z_2<\mathrm{}<z_A`$. According to Glauber theory  the wave function $`\psi _V`$ behind the first nucleon looks like

$`\psi _V^{(1)}(\stackrel{}{b},z)`$ $`=`$ $`\psi _V^{(0)}(\stackrel{}{b},z)\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_1)\psi _V^{(0)}(\stackrel{}{b},z_1)e^{ik_V(zz_1)}`$ (5)
$`=`$ $`\chi (\stackrel{}{b})\left(e^{ikz}\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_1)e^{iq_Vz_1}e^{ik_Vz}\right)`$ (6)
if one neglects off-diagonal scattering from one vector meson component into another (diagonal approximation). The phase factor in the second term arises from putting the vector meson component on its mass shell, i.e. $`k_V=\sqrt{\nu ^2m_V^2}`$ and $`q_V=kk_V`$ is the corresponding momentum transfer. The profile function $`\mathrm{\Gamma }_V`$ is related to the elastic vector meson nucleon scattering amplitude $`f_V`$ in the following way:

$$\mathrm{\Gamma }_V(\stackrel{}{b})=\frac{1}{2\pi i\stackrel{}{k}_V}e^{i\stackrel{}{k}_t\stackrel{}{b}}f_V(\stackrel{}{k}_t)d^2k_t.$$
(7)
After passing the second nucleon the wave function takes on the form

$$\psi _V^{(2)}(\stackrel{}{b},z)=\chi (\stackrel{}{b})(e^{ikz}\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_2)e^{iq_Vz_2}e^{ik_Vz}\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_1)(1\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_2)e^{iq_Vz_1}e^{ik_Vz}).$$
(8)
At position $`z`$ the whole ensemble of $`A`$ nucleons has led to the modification

$`\psi _V^{(A)}(\stackrel{}{b},z)`$ $`=`$ $`\chi (\stackrel{}{b})\left(e^{ikz}{\displaystyle \underset{j=1}{\overset{A}{}}}\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_j)e^{iq_Vz_j}{\displaystyle \underset{kj}{\overset{A}{}}}\left[1\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_k)\mathrm{\Theta }(z_kz_j)\right]\mathrm{\Theta }(zz_k)e^{ik_Vz}\right)`$ (9)
$`=`$ $`\left\{1{\displaystyle \underset{j=1}{\overset{A}{}}}\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_j)e^{iq_V(z_iz)}{\displaystyle \underset{kj}{\overset{A}{}}}\left[1\mathrm{\Gamma }_V(\stackrel{}{b}\stackrel{}{s}_k)\mathrm{\Theta }(z_kz_j)\right]\mathrm{\Theta }(zz_k)\right\}\psi _V^{(0)}(\stackrel{}{b},z)`$ (10)
$`=`$ $`\left\{1\mathrm{\Gamma }_V^{(A)}(\stackrel{}{b},z;\{\stackrel{}{r}_i\})\right\}\psi _V^{(0)}(\stackrel{}{b},z).`$ (11)
Hence, the physical photon state changes on its way through the nucleus to position $`\stackrel{}{r}=(\stackrel{}{b},z)`$ according to

$$|\gamma (\stackrel{}{r})=\left(1\underset{V=\rho ,\omega ,\varphi ,J/\mathrm{\Psi }}{}\frac{e^2}{2g_V^2}F_V^2\right)|\gamma _0+\underset{V=\rho ,\omega ,\varphi ,J/\mathrm{\Psi }}{}\frac{e}{g_V}F_V\left(1\overline{\mathrm{\Gamma }_V^{(A)}}(\stackrel{}{r})\right)|V$$
(12)
where the nuclear profile function from (11) must be averaged over the positions of the nucleons in the nucleus

$`\overline{\mathrm{\Gamma }_V^{(A)}}(\stackrel{}{b},z)`$ $`=`$ $`0|\mathrm{\Gamma }_V^{(A)}(\stackrel{}{b},z;\{\stackrel{}{r}_i\})|0`$ (13)
$`=`$ $`{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}𝑑z_in(\stackrel{}{b},z_i){\displaystyle \frac{\sigma _{VN}}{2}}(1i\alpha _V)e^{iq_V(z_iz)}\mathrm{exp}\left[{\displaystyle \frac{1}{2}}\sigma _{VN}(1i\alpha _V){\displaystyle \underset{z_i}{\overset{z}{}}}𝑑z_kn(\stackrel{}{b},z_k)\right].`$ (14)
To get to the second line of Eq. (14) we have followed the way of , i.e. making use of the independent particle model and the large $`A`$ limit and used the optical theorem to express the elastic forward scattering amplitude in terms of the total $`VN`$ cross section $`\sigma _V`$:

$$f_V(\stackrel{}{0})=\frac{ik_V}{4\pi }\sigma _V(1i\alpha _V).$$
(15)
The nucleon number density is denoted by $`n(\stackrel{}{r})`$ and normalized to $`A`$.

Neglecting any influence of the FSI, the reaction amplitude for the process $`\gamma Nf`$ on a nucleon at position $`\stackrel{}{r}`$ inside a nucleus changes compared to the vacuum due to these ’initial state interactions’ of the photon:

$$f|\widehat{T}|\gamma f|\widehat{T}|\gamma (\stackrel{}{r}).$$
(16)
By comparing Eq. (2) and (12) one sees that to account for shadowing in incoherent photoproduction one only has to multiply the amplitude of each vector meson states $`|V`$ by the corresponding factor $`\left(1\overline{\mathrm{\Gamma }_V^{(A)}}(\stackrel{}{r})\right)`$. This is equivalent to the method we used in .

We now compare our method of describing shadowing of incoherent electroproduction with conventional Glauber results. We therefore calculate the incoherent vector meson electroproduction cross section off a nucleus and assume purely absorptive FSI for the moment. After averaging over all positions $`\stackrel{}{r}=(\stackrel{}{b},z)`$ we get:

$`\sigma _{\gamma AVA^{}}`$ $`=`$ $`\sigma _{\gamma NVN}{\displaystyle d^2b\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑zn(\stackrel{}{b},z)}`$ (19)
$`\times \left|1{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}𝑑z_in(\stackrel{}{b},z_i){\displaystyle \frac{\sigma _V}{2}}(1i\alpha _V)e^{iq_V(z_iz)}\mathrm{exp}\left[{\displaystyle \frac{1}{2}}\sigma _V(1i\alpha _V){\displaystyle \underset{z_i}{\overset{z}{}}}𝑑z_kn(\stackrel{}{b},z_k)\right]\right|^2`$
$`\times \mathrm{exp}\left[\sigma _V^{inel}{\displaystyle \underset{z}{\overset{\mathrm{}}{}}}𝑑z^{}n(\stackrel{}{b},z^{})\right].`$
In (19) we again neglect off-diagonal scattering, i.e. $`V`$ production is triggered by the $`V`$ component of the photon only. The last factor in (19) includes the inelastic $`VN`$ cross section $`\sigma _V^{inel}`$ and accounts for the FSI.

The production process described by Eq.(19) can be understood as follows: The incoherent vector meson production takes place on a nucleon at position $`(\stackrel{}{b},z)`$. The production is triggered either directly by the photon or via an intermediate vector meson which was earlier produced on a nucleon at position $`(\stackrel{}{b},z_i)`$ without excitation of the nucleus. The interactions of this vector meson on its way from $`z_i`$ to position $`z`$ are of optical potential type and leave the nucleus in its ground state. This propagation is described by the second exponential in (19). The interference between the direct and the indirect process leads to shadowing. The nucleus gets excited in the incoherent reaction at position $`z`$. The possibility that the vector meson is lost on its way from position $`z`$ out of the nucleus is taken care of by the last exponential.

Expression (19) formally differs only slightly from the approximate result for the incoherent vector meson production cross section which is given without derivation in Ref. . The latter can be obtained by making the following replacements in (19):

$$e^{iq_V(z_iz)}e^{iq_Vz_i},\sigma _V^{inel}\sigma _V.$$
(20)
The last replacement simply means that one neglects the possibility of elastic $`VN`$ scattering whereas the occurrence of the different phase factor is unclear. The Authors of Ref.  have given a seemingly different result for incoherent $`V`$ electroproduction:

$`\sigma _{\gamma AVA^{}}`$ $`=`$ $`\sigma _{\gamma NVN}{\displaystyle }d^2b\{{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}dzn(\stackrel{}{b},z)e^{\sigma _V^{inel}T_z(\stackrel{}{b})}`$ (24)
$`+{\displaystyle \frac{1}{2}}{\displaystyle \frac{\sigma _V}{\sigma _V^{el}}}(\sigma _V^{inel}\sigma _V^{el}){\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}𝑑z_1n(\stackrel{}{b},z_1){\displaystyle \underset{z_1}{\overset{+\mathrm{}}{}}}𝑑z_2n(\stackrel{}{b},z_2)`$
$`\times \mathrm{cos}\left[q_V(z_1z_2)\right]e^{\frac{1}{2}(\sigma _V^{inel}\sigma _V^{el})T_{z_2}(\stackrel{}{b})\frac{1}{2}\sigma _VT_{z_1}(\stackrel{}{b})}`$
$`{\displaystyle \frac{1}{4}}{\displaystyle \frac{(\sigma _V)^2}{\sigma _V^{el}}}\left|{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}dzn(\stackrel{}{b},z)e^{iq_Vz}e^{\frac{1}{2}\sigma _VT_z(\stackrel{}{b})}|^2\right\}`$
with $`T_z(\stackrel{}{b})=_z^{\mathrm{}}n(\stackrel{}{b},z^{})𝑑z^{}`$ and $`T(\stackrel{}{b})=T_{z=\mathrm{}}(\stackrel{}{b})`$. Expression (24) can be interpreted as follows: The last term represents the coherent $`V`$ photoproduction cross section which is subtracted from the inclusive $`V`$ photoproduction cross section to yield the incoherent part. In Appendix A we show that the physically more transparent expression (19) is mathematically identical to (24).

### B Transport model

Up to now we have ignored the effects of Fermi motion, binding energies and Pauli blocking of final state nucleons in the process described by (16). In addition we neglected the finite life time of the $`\rho ^0`$ and took only absorptive FSI into account. These shortcomings are avoided when using the model of Ref. . The propagation of the final state $`|f`$ through the nucleus is treated within a semi-classical transport model based on the Boltzmann-Uehling-Uhlenbeck (BUU) equation. The BUU equation describes the time evolution of the phase space density $`f_i(\stackrel{}{r},\stackrel{}{p},t)`$ of particles of type $`i`$ that can interact via binary reactions. Besides the nucleons these particles involve baryonic resonances and mesons ($`\pi `$, $`\eta `$, $`\rho `$, $`K`$, …) that are produced either in the primary reaction or during the FSI. For a particle species $`i`$ the BUU equation can be written as:

$$\left(\frac{}{t}+\frac{H}{\stackrel{}{r}}\frac{}{\stackrel{}{r}}\frac{H}{\stackrel{}{r}}\frac{}{\stackrel{}{p}}\right)f_i(\stackrel{}{r},\stackrel{}{p},t)=I_{coll}[f_1,\mathrm{}f_i,\mathrm{},f_M].$$
(25)
For baryons the Hamilton function $`H`$ includes a mean field potential which depends on the particle position and momentum. The collision integral on the right hand side accounts for the creation and annihilation of particles of type $`i`$ in a collision as well as elastic scattering from one position in phase space into another. For fermions Pauli blocking is taken into account in $`I_{coll}`$ via blocking factors. The BUU equations of each particle species $`i`$ are coupled via the mean field and the collision integral. The resulting system of coupled differential-integral equations is solved via a test particle ansatz for the phase space density. For details of the transport model see Ref. .

The classes of FSI that are included in the transport model goes far beyond what can be achieved within Glauber theory. As a result the finally observed $`\rho ^0`$ does not need to be created in the primary reaction but might be produced during the FSI via side feeding. It is therefore clear that a purely absorptive treatment of the FSI as in Glauber theory can only be used if one is sure that one has eliminated the possibility of side feeding by applying enough constraints on the observable (see  for details).

We also stress that within our model instable particles might decay during their propagation through the nucleus. In case of the $`\rho ^0`$ this means that both pions have to escape the nucleus without further rescattering to make the identification of the $`\rho ^0`$ still possible.

## III Results

Before we turn to nuclear targets we first verify that the input of our model is reasonable. We therefore look at exclusive $`\rho ^0`$ production off hydrogen. We use the same kinematical cuts as in the HERMES experiment , i.e. the final state has to consist of two oppositely charged pions with invariant mass between $`0.6`$ GeV and $`1`$ GeV. The four momentum transfer $`|tt_{max}|`$ between the virtual photon and the $`\pi ^+\pi ^{}`$-pair has to be smaller than 0.4 GeV2 and we apply the exclusivity measure

$$\mathrm{\Delta }E=\frac{p_Y^2m_N^2}{2m_N}<0.4\text{ GeV},$$
(26)
where $`m_N`$ denotes the nucleon mass and

$$p_Y=p_N+p_\gamma p_\rho $$
(27)
the 4-momentum of the undetected final state. In Eq. (27) $`p_\gamma `$ and $`p_\rho `$ denote the 4-momenta of the incoming photon and the detected $`\pi ^+\pi ^{}`$ pair and $`p_N`$ is the 4-momentum of the struck nucleon which, for the calculation of $`p_Y`$, is assumed to be at rest. In Fig. 2 we compare our calculation of the exclusive $`\rho ^0`$ production cross section off hydrogen with experimental data . In the whole $`Q^2`$ region covered by the HERMES experiment we find very good agreement for a broad range of the invariant mass $`W`$ of the photon nucleon system. Also the slope of the differential production cross section is reproduced very well by the PYTHIA model as can be seen from the solid line in Fig. 3 where we show our calculation of $`\frac{d\sigma }{dt}(\gamma ^{}p\rho ^0p)`$ together with the HERMES data .

For our calculations on exclusive $`\rho ^0`$ electroproduction off nuclei we again use the kinematic cuts of the HERMES collaboration . This means that we restrict our exclusivity measure to the region

$$2\text{ GeV}<\mathrm{\Delta }E<0.6\text{ GeV},$$
(28)
and introduce a lower boundary for the four-momentum transfer $`|tt_{max}|>`$0.09 GeV2 as imposed by the HERMES collaboration to get rid of coherently produced $`\rho ^0`$. From the dashed line in Fig. 3 one sees that the differential $`\rho ^0`$ electroproduction cross section off 14N is again in excellent agreement with the HERMES data . Throughout our calculations the effect of the nucleon potential turns out to be negligible. This is reasonable since the involved energies are much larger than the typical binding energies which are in the order of a few MeV. A combined effect of Fermi motion and Pauli blocking on the incoherent differential production cross section is visible at $`|tt_{max}|<0.1`$ GeV2. At $`|tt_{max}|0.05`$ GeV2 the differential cross section is about 25% smaller than the calculation without Fermi motion and Pauli blocking (dotted curve). The reason for this reduction is that in the case that the bound nucleon moves towards the incoming photon the outgoing nucleon might be Pauli blocked for small momentum transfers. Note that this Pauli blocking just means that the whole nucleus absorbs the transferred momentum, i.e. this event contributes to the coherent production cross section which we do not consider here.

In Fig. 4 we show the transparency ratio $`T_A`$ for exclusive $`\rho ^0`$ production as a function of the coherence length $`l_\rho =q_\rho ^1`$. The solid line is the result that one gets if one uses Eq. (19) and accounts for two-body correlations by making the substitution

$$n(\stackrel{}{b},z_i)n(\stackrel{}{b},z_i)(1j_0(q_c|z_iz|))$$
(29)
in (19) with $`q_c=0.78`$ GeV. This is necessary to avoid unphysical contributions from processes where $`z_iz`$ which would contribute for small values of the coherence length $`l_\rho `$. In Ref.  we showed that this Bessel function parameterization yields a good description of shadowing in photoabsorption. In Eq. (19) we use for the total $`\rho ^0N`$ cross section $`\sigma _{\rho ^0}=25`$ mb and for the elastic part $`\sigma _{el}=3`$ mb. These two values correspond to the $`\rho N`$ cross sections used within the transport model for the involved $`\rho ^0`$ momenta.

The result of the transport model is represented by the open squares. For each HERMES data point  we have made a separate calculation with the corresponding $`\nu `$ and $`Q^2`$. In the case of 14N the Glauber and the transport calculation are in perfect agreement with each other and the experimental data. This demonstrates that, as we have discussed in Ref. , Glauber theory can be used for the FSI if the right kinematic constraints are applied. For comparison we also show the result of Huefner et al.  (dashed curve) as well as the result that one gets when using the approximate expression by Yennie  including two-body correlations (dotted curve). The somewhat larger transparency ratio of Huefner et al. arises from the different density distribution that the authors of Ref.  use which reduces the effect of FSI. The difference of the Yennie result (dotted curve) and our calculation arises mainly from the different cross section in the FSI (second replacement in (20)). The different phase factor leads to the change in the transparency ratio at small coherence lengths, since $`q_V=l_V^1`$.

After applying all of the above cuts, nearly all of the detected $`\rho ^0`$ stem from diffractive $`\rho ^0`$ production for which the formation time is zero. The 14N data seems to support the assumption that the time needed to put the preformed $`\rho ^0`$ fluctuation on its mass shell and let the wave function evolve to that of a physical $`\rho ^0`$ is small for the considered values of $`Q^2`$. Furthermore, the photon energy is too low to yield a large enough $`\gamma `$ factor to make the formation length exceed the internucleon distance and make color transparency visible. This conclusion is at variance with that reached in Ref. .

We now turn to 84Kr where we expect a stronger effect of the FSI. Unfortunately there is yet no data available to compare with. As can be seen from Fig. 4 the transport calculation for 84Kr gives a slightly smaller transparency ratio than the Glauber calculation, especially at low values of the coherence length, i.e. small momenta of the produced $`\rho ^0`$. There are two reasons for this: About 10% of the difference arises from the fact that within the transport model the $`\rho ^0`$ is allowed to decay into two pions. The probability that at least one of the pions interacts on its way out of the nucleus is about twice as large as that of the $`\rho ^0`$. The other reason is that in the Glauber calculation (19) only the inelastic part of the $`\rho ^0N`$ cross section enters whereas the transport calculation contains the elastic part as well. Thus all elastic scattering events out of the experimentally imposed $`t`$-window are neglected in the Glauber description. It is because of this $`t`$-window that also elastic $`\rho ^0N`$ scattering reduces the transport transparency ratio shown in Fig. 4. Both effects are more enhanced at lower energies and become negligible for the much smaller 14N nucleus.

## IV Summary and Outlook

We have developed a method to account for coherence length effects in incoherent electroproduction off nuclei which allows us to distinguish between the initial state interactions of the photon (shadowing) and the FSI of the reaction products. We have shown that our result is equivalent to the exact result of Glauber theory if one treats the FSI as a purely absorptive one. We have then performed a coupled channel treatment of the FSI within a semi-classical transport model and calculated the transparency ratio for exclusive incoherent $`\rho ^0`$ photoproduction off 14N and 84Kr. The result for 14N is in agreement with experimental data and with the Glauber prediction. The latter shows that in the case of 14N Glauber theory is applicable after the kinematic cuts of the HERMES experiment are applied. Since we do not use a formation time for diffractively produced vector mesons we deduce that one cannot see an onset of color transparency in the nitrogen data. For the 84Kr target no experimental data is available to compare with. However, we find deviations from the simple Glauber model because of the finite life time of the $`\rho ^0`$ and elastic scattering out of the kinematically allowed $`|t|`$-region. These effects should be taken into account when evaluating the 84Kr data in search of color transparency. As discussed by Kopeliovich et al.  one might see an onset of color transparency when investigating the transparency ratio as a function of $`Q^2`$ for fixed coherence length . Including a formation time for exclusive $`V`$ production is not straightforward within our model since a change of the $`VN`$ cross section during the vector meson formation will influence both the initial as well as the final state interactions. This is planned for future work, when hopefully also the 84Kr data has become available.

## Acknowledgments

The authors want to thank A. Borissov and W. Cassing for useful discussions. This work was supported by DFG.

## A

In the following we show the equality of Eqs. (19) and (24). Since the real part of the $`VN`$ scattering amplitude was neglected in the derivation of (24) we also set $`\alpha _V=0`$ in Equation (19) and end up with:

$`\sigma _{\gamma AVA^{}}`$ $`=`$ $`\sigma _{\gamma NVN}{\displaystyle d^2b\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑zn(\stackrel{}{b},z)e^{\sigma _V^{inel}_z^{\mathrm{}}𝑑z^{}n(\stackrel{}{b},z^{})}}`$ (31)
$`\times \left|1{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}𝑑z_in(\stackrel{}{b},z_i){\displaystyle \frac{\sigma _V}{2}}e^{iq_V(z_iz)}\mathrm{exp}\left[{\displaystyle \frac{1}{2}}\sigma _V{\displaystyle \underset{z_i}{\overset{z}{}}}𝑑z_kn(\stackrel{}{b},z_k)\right]\right|^2`$
$`=`$ $`\sigma _{\gamma NVN}{\displaystyle d^2b\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑zn(\stackrel{}{b},z)e^{\sigma _V^{inel}_z^{\mathrm{}}𝑑z^{}n(\stackrel{}{b},z^{})}}`$ (34)
$`\times \{1\sigma _V{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}dz_in(\stackrel{}{b},z_i)\mathrm{cos}\left[q_V(z_iz)\right]\mathrm{exp}[{\displaystyle \frac{1}{2}}\sigma _V{\displaystyle \underset{z_i}{\overset{z}{}}}dz_kn(\stackrel{}{b},z_k)]`$
$`+{\displaystyle \frac{\sigma _V^2}{4}}\left|{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}dz_in(\stackrel{}{b},z_i)e^{iq_Vz_i}\mathrm{exp}[{\displaystyle \frac{1}{2}}\sigma _V{\displaystyle \underset{z_i}{\overset{z}{}}}dz_kn(\stackrel{}{b},z_k)]|^2\right\}`$
$`=`$ $`\sigma _{\gamma NVN}{\displaystyle }d^2b\{{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}dzn(\stackrel{}{b},z)e^{\sigma _V^{inel}T_z(\stackrel{}{b})}`$ (39)
$`\sigma _V{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_1n(\stackrel{}{b},z_1){\displaystyle \underset{z_1}{\overset{\mathrm{}}{}}}𝑑z_2n(\stackrel{}{b},z_2)\mathrm{cos}\left[q_V(z_1z_2)\right]`$
$`\times \mathrm{exp}\left[\sigma _V^{inel}{\displaystyle \underset{z_2}{\overset{\mathrm{}}{}}}𝑑z^{}n(\stackrel{}{b},z^{}){\displaystyle \frac{1}{2}}\sigma _V{\displaystyle \underset{z_1}{\overset{z_2}{}}}𝑑z^{}n(\stackrel{}{b},z^{})\right]`$
$`+{\displaystyle \frac{\sigma _V^2}{4}}{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑zn(\stackrel{}{b},z)`$
$`\times |{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}dz_in(\stackrel{}{b},z_i)e^{iq_Vz_i}\mathrm{exp}[{\displaystyle \frac{1}{2}}\sigma _V{\displaystyle \underset{z_i}{\overset{z}{}}}dz^{}n(\stackrel{}{b},z^{}){\displaystyle \frac{1}{2}}\sigma _V^{inel}{\displaystyle \underset{z}{\overset{\mathrm{}}{}}}dz^{}n(\stackrel{}{b},z^{})]|^2\}.`$
In the second step we have renamed the integration variables and rewritten the integral limits. The first term of (31) already equals that the first term in (24). The exponent in the second term of (31) yields

$$\sigma _V^{inel}T_{z_2}(\stackrel{}{b})\frac{1}{2}\sigma _V\left(T_{z_1}(\stackrel{}{b})T_{z_2}(\stackrel{}{b})\right)=\frac{1}{2}\left(\sigma _V^{inel}\sigma _V^{el}\right)T_{z_2}(\stackrel{}{b})\frac{1}{2}\sigma _VT_{z_1}(\stackrel{}{b})$$
(40)
and the exponent of the third term in (31) is

$$\frac{1}{2}\left(\sigma _VT_{z_i}(\stackrel{}{b})T_z(\stackrel{}{b})\right)\frac{1}{2}\sigma _V^{inel}T_z(\stackrel{}{b})=\frac{1}{2}\sigma _VT_{z_i}(\stackrel{}{b})+\frac{1}{2}\sigma _V^{el}T_z(\stackrel{}{b}).$$
(41)
We now further manipulate the last term of (31):

$`{\displaystyle \frac{\sigma _V^2}{4}}`$ $`{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑zn(\stackrel{}{b},z)e^{\sigma _V^{el}T_z(\stackrel{}{b})}\left|{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}𝑑z_in(\stackrel{}{b},z_i)e^{iq_Vz_i}e^{\frac{1}{2}\sigma _VT_{z_i}(\stackrel{}{b})}\right|^2`$ (43)
$`={\displaystyle \frac{\sigma _V^2}{4}}{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}𝑑z_1{\displaystyle \underset{\mathrm{}}{\overset{z}{}}}𝑑z_2n(\stackrel{}{b},z)n(\stackrel{}{b},z_1)n(\stackrel{}{b},z_2)e^{\sigma _V^{el}T_z(\stackrel{}{b})}\mathrm{cos}\left[q_V(z_1z_2)\right]e^{\frac{1}{2}\sigma _V(T_{z_1}(\stackrel{}{b})+T_{z_2}(\stackrel{}{b}))}.`$
Rewriting the integral in the following form

$$\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑z\underset{\mathrm{}}{\overset{z}{}}𝑑z_1\underset{\mathrm{}}{\overset{z}{}}𝑑z_2=\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑z_2\underset{\mathrm{}}{\overset{z_2}{}}𝑑z_1\underset{z_2}{\overset{\mathrm{}}{}}𝑑z+\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑z_1\underset{\mathrm{}}{\overset{z_1}{}}𝑑z_2\underset{z_1}{\overset{\mathrm{}}{}}𝑑z$$
(44)
and using the symmetry of the integrand with respect to the variables $`z_1`$ and $`z_2`$ we can rewrite (43) as

$`{\displaystyle \frac{\sigma _V^2}{2}}`$ $`{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_2n(\stackrel{}{b},z_2){\displaystyle \underset{\mathrm{}}{\overset{z_2}{}}}𝑑z_1n(\stackrel{}{b},z_1)\mathrm{cos}\left[q_V(z_1z_2)\right]e^{\frac{1}{2}\sigma _V(T_{z_1}(\stackrel{}{b})+T_{z_2}(\stackrel{}{b}))}{\displaystyle \underset{z_2}{\overset{\mathrm{}}{}}}𝑑zn(\stackrel{}{b},z)e^{\sigma _V^{el}T_z(\stackrel{}{b})}`$ (48)
$`={\displaystyle \frac{1}{2}}{\displaystyle \frac{\sigma _V^2}{\sigma _V^{el}}}{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_2n(\stackrel{}{b},z_1){\displaystyle \underset{\mathrm{}}{\overset{z_2}{}}}𝑑z_1n(\stackrel{}{b},z_2)\mathrm{cos}\left[q_V(z_1z_2)\right]e^{\frac{1}{2}\sigma _V(T_{z_1}(\stackrel{}{b})+T_{z_2}(\stackrel{}{b}))}\left(1e^{\sigma _V^{el}T_{z_2}}\right)`$
$`={\displaystyle \frac{1}{2}}{\displaystyle \frac{\sigma _V^2}{\sigma _V^{el}}}{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_1n(\stackrel{}{b},z_1){\displaystyle \underset{z_1}{\overset{\mathrm{}}{}}}𝑑z_2n(\stackrel{}{b},z_2)\mathrm{cos}\left[q_V(z_1z_2)\right]e^{\frac{1}{2}\sigma _V(T_{z_1}(\stackrel{}{b})+T_{z_2}(\stackrel{}{b}))}`$
$`+{\displaystyle \frac{1}{2}}{\displaystyle \frac{\sigma _V^2}{\sigma _V^{el}}}{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_1n(\stackrel{}{b},z_1){\displaystyle \underset{z_1}{\overset{\mathrm{}}{}}}𝑑z_2n(\stackrel{}{b},z_2)\mathrm{cos}\left[q_V(z_1z_2)\right]e^{\frac{1}{2}(\sigma _V^{inel}\sigma _V^{el})T_{z_2}(\stackrel{}{b})\frac{1}{2}\sigma _VT_{z_1}(\stackrel{}{b}))}.`$
In the first step we have performed the integral over $`z`$. In the last step we have again rewritten the integral over $`z_1`$ and $`z_2`$ to be able to combine the second term of (48) with the second term of (31) to get the second term in (24). Now we only have to show the equality of the last term in (24) and the first one in (48). Therefore we rewrite

$`{\displaystyle \frac{1}{4}}`$ $`{\displaystyle \frac{\sigma _V^2}{\sigma _V^{el}}}\left|{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑zn(\stackrel{}{b},z)e^{iq_Vz}e^{\frac{1}{2}\sigma _VT_z(\stackrel{}{b})}\right|^2`$ (50)
$`={\displaystyle \frac{1}{4}}{\displaystyle \frac{\sigma _V^2}{\sigma _V^{el}}}{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_1{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_2n(\stackrel{}{b},z_1)n(\stackrel{}{b},z_2)\mathrm{cos}\left[q_V(z_1z_2)\right]e^{\frac{1}{2}\sigma _V(T_{z_1}(\stackrel{}{b})+T_{z_2}(\stackrel{}{b}))}`$
and reformulate the integral as

$`{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_1{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_2`$ $`=`$ $`{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_1{\displaystyle \underset{z_1}{\overset{\mathrm{}}{}}}𝑑z_2+{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_2{\displaystyle \underset{z_2}{\overset{\mathrm{}}{}}}𝑑z_1`$ (51)
$`=`$ $`2{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}𝑑z_1{\displaystyle \underset{z_1}{\overset{\mathrm{}}{}}}𝑑z_2`$ (52)
where the second equality again follows from the symmetry of the integrand with respect to the exchange of $`z_1`$ and $`z_2`$. The two expressions (19) and (24) are therefore equivalent.


text:
# Exact two-qubit universal quantum circuit

## Preliminary

We first briefly introduce some basic facts about Cartan decomposition of $`\text{SU}(4)`$ and local equivalence of two-qubit gates. Any two-qubit unitary operation $`U\text{SU}(4)`$ can be decomposed as Khaneja et al. (2001); Kraus and Cirac (2001); Zhang et al. (2003)

$$U=k_1e^{c_1\frac{i}{2}\sigma _x^1\sigma _x^2}e^{c_2\frac{i}{2}\sigma _y^1\sigma _y^2}e^{c_3\frac{i}{2}\sigma _z^1\sigma _z^2}k_2,$$
(1)
where $`\sigma _\alpha ^1\sigma _\alpha ^2=\sigma _\alpha \sigma _\alpha `$, $`\sigma _\alpha `$ are the Pauli matrices, and $`k_1`$, $`k_2\text{SU}(2)\text{SU}(2)`$ are local gates. The tetrahedral representation of nonlocal gates in Zhang et al. (2003) defines a unique set of coefficients $`c_j`$ satisfying:

$$\pi c_2c_1c_2c_30.$$
(2)
Also from Zhang et al. (2003), we know that local gates $`U\text{SU}(2)\text{SU}(2)`$ correspond to the case when $`c_1=c_2=c_3=0`$, or $`c_1=\pi `$, $`c_2=c_3=0`$; the SWAP gate to $`c_1=c_2=c_3=\frac{\pi }{2}`$; and the CNOT gate to $`c_1=\frac{\pi }{2}`$, $`c_2=c_3=0`$.

Two unitary transformations $`U`$, $`U_1\text{SU}(4)`$ are called *locally equivalent* if they differ only by local operations: $`U=k_1U_1k_2`$, where $`k_1`$, $`k_2\text{SU}(2)\text{SU}(2)`$. It was shown in Zhang et al. (2003) that any nonlocal two-qubit operation that is not locally equivalent to the SWAP gate is entangling. The SWAP gate and its local equivalence class are thus the only nonlocal two-qubit operations that transform unentangled states to unentangled states, i.e., that do not introduce any entanglement.

## Universal quantum circuit

We now present an analytic way to implement any arbitrary two-qubit gate $`U\text{SU}(4)`$ by constructing a closed form solution for a universal quantum circuit that is composed of a small number of repetitions of a given entangling operation $`U_g`$ together with local gates. Local gates are assumed to be implementable at ease, as is the case in many of the current proposed physical implementations of quantum computation Nielsen and Chuang (2000).

An arbitrary two-qubit operation $`U\text{SU}(4)`$ can be written as in Eq. (1). Letting $`k_x=e^{\frac{\pi }{4}i\sigma _y}e^{\frac{\pi }{4}i\sigma _y}`$ and $`k_y=e^{\frac{\pi }{4}i\sigma _x}e^{\frac{\pi }{4}i\sigma _x}`$, we have

$`\sigma _x^1\sigma _x^2=k_x^{}\sigma _z^1\sigma _z^2k_x,\sigma _y^1\sigma _y^2=k_y^{}\sigma _z^1\sigma _z^2k_y.`$ (3)
Substituting Eq. (3) into Eq. (1), we find that an arbitrary two-qubit gate $`U\text{SU}(4)`$ can be written as

$`U=(k_1k_x^{})e^{c_1\frac{i}{2}\sigma _z^1\sigma _z^2}(k_xk_y^{})e^{c_2\frac{i}{2}\sigma _z^1\sigma _z^2}(k_y)e^{c_3\frac{i}{2}\sigma _z^1\sigma _z^2}(k_2),`$ (4)
where $`k_1`$, $`k_2`$, $`k_x`$ and $`k_y`$ are all local gates. Since we have all the local gates at our full disposal, it is evident that we only need to implement the nonlocal block $`e^{c_j\frac{i}{2}\sigma _z^1\sigma _z^2}`$ from the given entangling gate $`U_g`$ together with local gates, for general values of $`c_j`$ between 0 and $`\pi `$. We have found the following analytic construction of this general block:

Apply $`U_g`$ at most twice to build a gate $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ with $`\gamma (0,\frac{\pi }{2}]`$ (Proposition 1). The value of $`\gamma `$ obtained depends on the starting $`U_g`$.

Apply $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ $`n`$ times, until $`n\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$.

Apply $`e^{n\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ with $`n\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$ twice, to simulate the nonlocal block $`e^{c_j\frac{i}{2}\sigma _z^1\sigma _z^2}`$ (Proposition 2).

Build the quantum circuit according to Eq. (4).

We now describe this construction in more detail.

###### Proposition 1

Any arbitrary given entangling gate $`U_g`$ can simulate a gate $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$, where $`\gamma (0,\frac{\pi }{2}]`$, by a quantum circuit that applies $`U_g`$ at most twice.

###### Proof.

From Sec. a, an arbitrary entangling gate $`U_g`$ can be uniquely written as $`U_g=k_lAk_r`$, where

$`A=e^{\gamma _1\frac{i}{2}\sigma _x^1\sigma _x^2}e^{\gamma _2\frac{i}{2}\sigma _y^1\sigma _y^2}e^{\gamma _3\frac{i}{2}\sigma _z^1\sigma _z^2},`$
$`k_l`$, $`k_r\text{SU}(2)\text{SU}(2)`$, and $`\pi \gamma _2\gamma _1\gamma _2\gamma _30`$. Since $`U_g`$ is an entangling gate, i.e., it is neither local nor locally equivalent to the SWAP gate, we do not need to take the following cases into account and may ignore them: (1) $`\gamma _1=\gamma _2=\gamma _3=0`$; (2) $`\gamma _1=\pi `$, $`\gamma _2=\gamma _3=0`$; and (3) $`\gamma _1=\gamma _2=\gamma _3=\frac{\pi }{2}`$. The first two cases correspond to local gates, and the third one to the SWAP gate. For all the remaining possibilities of $`\gamma _j`$, we distinguish the following four cases:

*Case 1:* If $`\gamma _3=\gamma _2=0`$ and $`\gamma _1(0,\pi )`$, then $`k_xAk_x^{}=e^{\gamma _1\frac{i}{2}\sigma _z^1\sigma _z^2}`$ is in the required form.

*Case 2:* If $`\gamma _3=0`$ and $`\gamma _1=\gamma _2=\frac{\pi }{2}`$, we construct the quantum circuit

$`e^{\frac{\pi }{4}i\sigma _y^1}Ae^{\frac{\pi }{4}i\sigma _y^2}e^{\frac{\pi }{4}i\sigma _z^1\frac{\pi }{4}i\sigma _z^2}Ae^{\frac{\pi }{4}i\sigma _z^1+\frac{\pi }{4}i\sigma _z^2}e^{\frac{\pi }{4}i\sigma _y^1}=e^{\frac{\pi }{4}i\sigma _z^1\sigma _z^2}.`$
*Case 3:* If $`\gamma _3=0`$, $`\gamma _1(0,\pi )`$, and $`\gamma _2(0,\frac{\pi }{2})`$, we have the quantum circuit as $`Ae^{\frac{i\pi }{2}\sigma _x^1}Ae^{\frac{i\pi }{2}\sigma _x^1}=e^{2\gamma _1\frac{i}{2}\sigma _z^1\sigma _z^2}`$.

*Case 4:* If $`\gamma _3(0,\frac{\pi }{2})`$, we use the quantum circuit $`Ae^{\frac{i\pi }{2}\sigma _z^1}Ae^{\frac{i\pi }{2}\sigma _z^1}=e^{2\gamma _3\frac{i}{2}\sigma _z^1\sigma _z^2}`$.

These four cases exhaust all the possibilities of $`\gamma _j`$. We have thus obtained a closed form solution to construction of the gate $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ from $`U_g`$ with local gates, where $`\gamma (0,2\pi )`$ and $`\gamma \pi `$. We can always choose $`\gamma (0,\pi )`$ since

$`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}=ie^{\frac{\pi i}{2}\sigma _z^1}e^{(\pi +\gamma )\frac{i}{2}\sigma _z^1\sigma _z^2}e^{\frac{\pi i}{2}\sigma _z^2}.`$ (5)
Furthermore, if $`\gamma [\frac{\pi }{2},\pi )`$, we can construct the following quantum circuit to bring $`\gamma `$ into the interval $`(0,\frac{\pi }{2}]`$:

$`ie^{\frac{\pi i}{2}\sigma _z^1}e^{\frac{\pi i}{2}\sigma _y^1}e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}e^{\frac{\pi i}{2}\sigma _y^1}e^{\frac{\pi i}{2}\sigma _z^2}=e^{(\pi \gamma )\frac{i}{2}\sigma _z^1\sigma _z^2}.`$ (6)
The proof of the Proposition is thereby complete. ∎

We have thus shown that given an entangling gate $`U_g`$ together with local gates, we can construct a gate $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ with $`\gamma (0,\frac{\pi }{2}]`$ in an analytic form (Step 1). Now it is evident that if the constructed gate has $`\gamma (0,\frac{\pi }{4}]`$, then it can be applied for $`n`$ times until $`n\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$ (Step 2). In the next Proposition, we will use the resulting gate $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ with $`\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$, as a basic building block to simulate any generic nonlocal block $`e^{c\frac{i}{2}\sigma _z^1\sigma _z^2}`$ (Step 3). From Eqs. (5) and (6), we only need to consider the case when $`c(0,\frac{\pi }{2}]`$.

###### Proposition 2

Given a gate $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$, where $`\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$, together with local gates, the following quantum circuit can simulate the gate $`e^{c\frac{i}{2}\sigma _z^1\sigma _z^2}`$ for any $`c(0,\frac{\pi }{2}]`$:

In the above quantum circuit, we have

$`U_1=\left(\begin{array}{cc}ip& iq\\ q& p\end{array}\right),U_2=\left(\begin{array}{cc}ip& q\\ iq& p\end{array}\right),`$ (7)
$`b=\mathrm{cos}^1\left((\mathrm{cos}c\mathrm{cos}^2\gamma )/\mathrm{sin}^2\gamma \right),`$ (8)
where

$`p=\sqrt{{\displaystyle \frac{1}{2}}\left(1+{\displaystyle \frac{\mathrm{tan}\frac{c}{2}}{\mathrm{tan}\gamma }}\right)},q=\sqrt{{\displaystyle \frac{1}{2}}\left(1{\displaystyle \frac{\mathrm{tan}\frac{c}{2}}{\mathrm{tan}\gamma }}\right)}.`$ (9)
###### Proof.

We first justify the condition $`\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$. From Eq. (8), we have $`\mathrm{cos}c=\mathrm{sin}^2\gamma \mathrm{cos}b+\mathrm{cos}^2\gamma `$. Therefore, $`\mathrm{cos}2\gamma \mathrm{cos}c1`$, which yields that $`0c2\gamma `$. To cover the full range $`(0,\frac{\pi }{2}]`$ of $`c`$, we therefore require that $`\gamma \frac{\pi }{4}`$.

We now derive a few formulas required for the proof below. It is straightforward to show that $`p^2+q^2=1`$. The identity $`\mathrm{sin}\frac{c}{2}=\mathrm{sin}\gamma \mathrm{sin}\frac{b}{2}`$ follows from Eq. (8) by direct derivations. This yields

$$\mathrm{tan}^2\gamma =\frac{\mathrm{sin}^2\frac{c}{2}}{\mathrm{sin}^2\frac{b}{2}\mathrm{sin}^2\frac{c}{2}},$$

whence

$$pq=\frac{1}{2}\sqrt{1\frac{\mathrm{tan}^2\frac{c}{2}}{\mathrm{tan}^2\gamma }}=\frac{\mathrm{cos}\frac{b}{2}}{2\mathrm{cos}\frac{c}{2}}.$$
(10)
The Proposition can now be proved. Since $`p^2+q^2=1`$, it is easy to see that $`U_1U_1^{}=U_2U_2^{}=I`$. Hence, $`U_1`$ and $`U_2`$ are indeed single-qubit gates. The quantum circuit can be rewritten as:

$`(IU_1)e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}(Ie^{(b+\pi )\frac{i}{2}\sigma _y})e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}(IU_2)`$ (11)
$`=`$ $`\left(\begin{array}{cc}W& \\ & V\end{array}\right),`$
where

$`W`$ $`=`$ $`U_1e^{i\frac{\gamma }{2}\sigma _z}e^{(b+\pi )\frac{i}{2}\sigma _y}e^{i\frac{\gamma }{2}\sigma _z}U_2,`$ (12)
$`V`$ $`=`$ $`U_1e^{i\frac{\gamma }{2}\sigma _z}e^{(b+\pi )\frac{i}{2}\sigma _y}e^{i\frac{\gamma }{2}\sigma _z}U_2.`$ (13)
After substituting Eq. (7) into Eq. (12) and applying the identities $`\mathrm{sin}\frac{c}{2}=\mathrm{sin}\gamma \mathrm{sin}\frac{b}{2}`$ and Eq. (10), we obtain

$`W_{11}`$ $`=`$ $`p^2\mathrm{sin}{\displaystyle \frac{b}{2}}e^{i\gamma }+2pq\mathrm{cos}{\displaystyle \frac{b}{2}}q^2\mathrm{sin}{\displaystyle \frac{b}{2}}e^{i\gamma }=e^{c\frac{i}{2}},`$
$`W_{22}`$ $`=`$ $`p^2\mathrm{sin}{\displaystyle \frac{b}{2}}e^{i\gamma }+2pq\mathrm{cos}{\displaystyle \frac{b}{2}}q^2\mathrm{sin}{\displaystyle \frac{b}{2}}e^{i\gamma }=e^{c\frac{i}{2}},`$
$`W_{12}`$ $`=`$ $`W_{21}=2ipq\mathrm{sin}{\displaystyle \frac{b}{2}}\mathrm{cos}\gamma i(p^2q^2)\mathrm{cos}{\displaystyle \frac{b}{2}}=0.`$
Hence $`W=e^{c\frac{i}{2}\sigma _z}`$. Similarly, we find $`V=e^{c\frac{i}{2}\sigma _z}`$. Eq. (11) now becomes

$`(IU_1)e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}(Ie^{(b+\pi )\frac{i}{2}\sigma _y})e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}(IU_2)`$
$`=`$ $`\left(\begin{array}{cc}e^{c\frac{i}{2}\sigma _z}& \\ & e^{c\frac{i}{2}\sigma _z}\end{array}\right)=e^{c\frac{i}{2}\sigma _z^1\sigma _z^2},`$
which completes the proof. ∎

Note that in the above Proposition, for the extreme case when $`\gamma =\frac{\pi }{2}`$, corresponding to starting from a CNOT gate, we have $`b=c`$, and

$`U_1={\displaystyle \frac{1}{\sqrt{2}}}\left(\begin{array}{cc}i& i\\ 1& 1\end{array}\right),U_2={\displaystyle \frac{1}{\sqrt{2}}}\left(\begin{array}{cc}i& 1\\ i& 1\end{array}\right).`$
As a physical example, let us consider neutral atoms in an optical lattice as a simulator for a solid state many-body spin system. The simulation objective may for instance be implementation of $`\sqrt{\text{SWAP}}`$. While this is readily generated in spin systems, from the isotropic exchange Hamiltonian Burkard et al. (1999), it is not directly accessible for neutral atoms in optical lattices. A convenient experimentally accessible nonlocal transformation in this setting is the Controlled-PHASE gate $`C_\varphi `$, where the PHASE gate is $`\left(\begin{array}{cc}1& \\ & e^{i\varphi }\end{array}\right)`$. From the Cartan decomposition, we have

$`C_\varphi =e^{i\frac{\varphi }{4}}e^{i\frac{\varphi }{4}\sigma _z}e^{i\frac{\varphi }{4}\sigma _z}e^{\frac{\varphi }{2}\frac{i}{2}\sigma _z^1\sigma _z^2}.`$ (14)
On the other hand, the $`\sqrt{\text{SWAP}}`$ gate can be written as:

$`\sqrt{\text{SWAP}}=e^{i\frac{\pi }{8}}e^{\frac{\pi }{4}\frac{i}{2}\sigma _x^1\sigma _x^2}e^{\frac{\pi }{4}\frac{i}{2}\sigma _y^1\sigma _y^2}e^{\frac{\pi }{4}\frac{i}{2}\sigma _z^1\sigma _z^2}.`$
From the quantum circuit in Eq. (4), since $`c_1=c_2=c_3=\frac{\pi }{4}`$, we need only to implement the nonlocal gate $`e^{\frac{\pi }{4}\frac{i}{2}\sigma _z^1\sigma _z^2}`$. For $`\varphi [\frac{\pi }{2},\pi ]`$, from Proposition 2 and Eq. (14), we get

$`e^{\frac{\pi }{4}\frac{i}{2}\sigma _z^1\sigma _z^2}`$ $`=`$ $`e^{i\frac{\varphi }{2}}(e^{i\frac{\varphi }{4}\sigma _z}U_1e^{i\frac{\varphi }{4}\sigma _z})C_\varphi (e^{i\frac{\varphi }{4}\sigma _z}`$
$`e^{\frac{i}{2}(b+\pi )\sigma _y}e^{i\frac{\varphi }{4}\sigma _z})C_\varphi (IU_2),`$
where $`b=\mathrm{cos}^1\left((\frac{1}{\sqrt{2}}\mathrm{cos}^2\frac{\varphi }{2})/\mathrm{sin}^2\frac{\varphi }{2}\right)`$, $`U_1`$ and $`U_2`$ are given as in Eq. (7), and

$`p=\sqrt{{\displaystyle \frac{1}{2}}\left(1+{\displaystyle \frac{\sqrt{2}1}{\mathrm{tan}\frac{\varphi }{2}}}\right)},q=\sqrt{{\displaystyle \frac{1}{2}}\left(1{\displaystyle \frac{\sqrt{2}1}{\mathrm{tan}\frac{\varphi }{2}}}\right)}.`$
Thus a spin-spin interaction can be simulated in an optical lattice with only 2 repetitions of a Controlled-PHASE gate $`C_\varphi `$ having $`\varphi [\frac{\pi }{2},\pi ]`$.

## Uniform upper bound

One often desires to simulate arbitrary two-qubit operation by applying the given entangling two-qubit operation as infrequently as possible. From the construction procedure described above, we first use the given entangling gate $`U_g`$ to implement a gate $`U_f=e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ with $`\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$ (Proposition 1), and then apply $`U_f`$ twice, to implement a generic nonlocal gate $`e^{c\frac{i}{2}\sigma _z^1\sigma _z^2}`$ (Proposition 2). From the decomposition of $`\text{SU}(4)`$ in Eq. (4), any arbitrary two-qubit unitary operation contains at most three such nonlocal blocks, resulting in the quantum circuit

where each nonlocal block $`e^{c_j\frac{i}{2}\sigma _z^1\sigma _z^2}`$ is simulated as shown in the circuit of Proposition 2. It is clear that overall we only need to apply the gate $`U_f`$ at most $`6`$ times, in order to simulate an arbitrary two-qubit operation. We thereby obtain an upper bound for the applications of the given entangling gate $`U_g`$ to construct an exact universal quantum circuit. The value of this upper bound depends only on the nonlocal part of the given gate. For example, for a Controlled-PHASE gate $`C_\varphi `$ with parameter $`\varphi [\frac{\pi }{2},\pi ]`$, it takes 6 applications of $`C_\varphi `$ and 7 local gates to simulate any arbitrary two-qubit operation. However, when $`\varphi (0,\frac{\pi }{2})`$, we first need to apply $`C_\varphi `$ $`n`$ times until $`n\varphi \frac{\pi }{2}`$. Consequently, in this case it takes $`6n`$ applications of $`C_\varphi `$ and $`6n+1`$ local gates to implement any arbitrary two-qubit operation. Furthermore, this upper bound is uniform in the sense that no matter which two-qubit unitary operation is to be implemented, we can always construct a quantum circuit to simulate this operation with applications of the given entangling gate $`U_g`$ not exceeding the upper bound.

It is instructive to compare these results with the numerical solution for construction of CNOT obtained in Bremner et al. (2002). For the gate $`U_g=e^{\frac{\pi }{3}\frac{i}{2}\sigma _z^1\sigma _z^2}`$, both procedures need only two applications to obtain the CNOT gate. When $`U_g=e^{\frac{\pi }{5}\frac{i}{2}\sigma _z^1\sigma _z^2}`$, our uniform construction requires four applications, whereas the procedure of Bremner et al. (2002) only needs three applications to get CNOT. This difference derives from the fact that our procedure provides a uniform solution and is not optimized for any specific gate, whereas the procedure of Bremner et al. (2002) is near optimal for CNOT. This comparison reveals that the uniform property and optimality cannot necessarily be satisfied simultaneously.

Our final analysis concerns the efficiency of these analytic circuits. We first show that our basic building blocks of the quantum circuit, namely the gates $`e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$, are locally equivalent to the Controlled-Unitary (Controlled-$`U`$) gates.

###### Proposition 3

Consider an arbitrary single-qubit gate $`U=\mathrm{exp}\{i\gamma \widehat{n}\stackrel{}{\sigma }\}`$, where $`\gamma ^+`$, $`\widehat{n}=(n_x,n_y,n_z)`$ is a unit vector in $`^3`$, and $`\stackrel{}{\sigma }`$ denotes the vector $`(\sigma _x,\sigma _y,\sigma _z)`$ of Pauli matrices. The corresponding Controlled-$`U`$ gate can be simulated by the following quantum circuit:

where

$$U_1=\{\begin{array}{cc}\left(\begin{array}{cc}i\sqrt{\frac{1n_z}{2}}& \sqrt{\frac{1+n_z}{2}}\\ \frac{n_yn_xi}{\sqrt{2\left(1n_z\right)}}& \frac{n_x+n_yi}{\sqrt{2\left(1+n_z\right)}}\end{array}\right),& \text{for }n_z\pm 1;\hfill \\ \sigma _x,& \text{for }n_z=1;\hfill \\ I,& \text{for }n_z=1.\hfill \end{array}$$

This proposition can be proved by substitution of $`U_1`$, followed by direct algebraic computation 111See EPAPS Document No. E-PRLXXX-00-000000 for proof. A direct link to this document may be found in the online article’s HTML reference section. This document may also be reached via the EPAPS homepage (http://www.aip.org/pubservs/epaps.html) or from ftp.aip.org in the directory /epaps/. See the EPAPS homepage for more information..

We saw above that all the gates $`U_f=e^{\gamma \frac{i}{2}\sigma _z^1\sigma _z^2}`$ with $`\gamma [\frac{\pi }{4},\frac{\pi }{2}]`$ have the same upper bounds. When $`\gamma =\frac{\pi }{2}`$, $`U_f`$ is locally equivalent to the CNOT gate. From Zhang et al. (2003), there exists a one-to-one map from the local equivalence classes of Controlled-$`U`$ gates to the points in the interval $`[0,\frac{\pi }{2}]`$. Furthermore, those gates having the same upper bounds as the CNOT gate constitute half of this interval, namely $`[\frac{\pi }{4},\frac{\pi }{2}]`$. Therefore, using the length of the interval as a measure, we conclude that exactly half of the Controlled-$`U`$ gates can be used to construct universal quantum circuits that satisfy the same upper bound, i.e., they can be used just as efficiently as the CNOT gate. (Note that this is true despite the fact that CNOT is the only Controlled-$`U`$ gate providing perfect entanglement Zhang et al. (2003).) Consequently, there is no need to restrict practical studies of physical implementation of quantum circuits for universal computation or for quantum simulations to the standard model of CNOT with local gates.

In summary, we have provided an analytic approach to construct a universal quantum circuit that can simulate any arbitrary two-qubit operation given any entangling gate $`U_g`$ supplemented with local gates. Closed form solutions have been derived for each step in this explicit construction procedure. The procedure was illustrated on a physical example of simulation of a solid state spin system with neutral atoms in an optical lattice. Our approach provides a uniform upper bound for the applications of the given entangling gate $`U_g`$. It was found that precisely half of all the Controlled-$`U`$ gates have the same uniform upper bounds as the CNOT gate, i.e., they are equally efficient as the CNOT gate for providing realizable implementation of arbitrary two-qubit operation. This offers new options for realization of interactions in simulation of one quantum many-body system by another, as well as for efficient implementation of quantum computation.

We thank the NSF for financial support under ITR Grant No. EIA-0205641 (SS and KBW). KBW thanks the Miller Institute for Basic Research for a Miller Research Professorship 2002-2003.


text:
# Thermodynamics of ionization and dissociation in hydrogen plasmas including fluctuations and magnetic fields

## 1 Introduction

In plasmas and in particular in hydrogen systems, ionization and dissociation processes play an important role. We find the equilibrium values of the degrees of ionization and dissociation by minimization of the free energy at fixed total mass ebkrkr76 ; ebetal91 ; beuleal99 . Further the ionization and dissociation rates are derived from kinetic equations ebetal91 ; klim . Here we will start with a study of the geometry of the free energy landscape ebfree90 ; ebhi00 . The idea to apply geometrical methods in thermodynamical research goes back to pioneering papers of Gibbs gibbs1892 . We apply Gibb’s method here to plasmas. We will study different approximations to the thermodynamic functions and the influence of several effects, our approach is limited to the region $`T>10000K`$. In a preceding work we have developed a similar approach for ionization processes described by the Saha equation ebhi00 . We include here also effects of magnetic fields which are of increasing interest not only for astrophysical applications but also for plasmas confined in magnetic traps. A preliminary work on magnetic field effects was published earlier ebstor00 . The magnetic field introduces an anisotropy into the system. This changes the ideal gas contribution, the bound state energies and the contributions due to interactions of the free particles. Here we will take into account the ideal gas corrections in an exact way, and the bound state and interaction corrections in approximations which include the quadratic orders in the field $`O(B^2)`$. Practically this restricts our approach to fields $`B<10^5T`$.

We start with exploring the geometry of the free energy landscape as a function of the degrees of ionization and dissociation. On this basis we are able to describe not only the chemical equilibria but also fluctuations around the equilibrium values ebfree90 .

## 2 The free energy density of hydrogen

### 2.1 Defining the basic variables of the system

We consider in the following hydrogen plasmas at fixed temperature $`T`$, and proton density $`n`$. Our aim is to study ionization processes

$$Hp+e^{}$$
(1)
and dissociation processes

$$H_2H+H.$$
(2)
For simplicity, the formation of $`H_2^+`$ and $`H^{}`$ species will be neglected.

In order to define the thermodynamics of the system we introduce the free energy

$$F(T,V,\left\{N\right\}),$$
(3)
the density of the free energy

$$f(T,n)=\frac{F}{V}=f(T,n_e,n_i,n_a,n_m),$$
(4)
the degrees of ionization and dissociation beuleal99 . We will use here the following definitions:

$$\begin{array}{cccc}\hfill \text{}\alpha =& \frac{n_i}{n_i+n_a+2n_m}\hfill & \hfill \beta =& \frac{n_a}{n_a+2n_m}\hfill \\ \hfill \text{}\beta _a=& \frac{n_a}{n_i+n_a+2n_m}\hfill & \hfill \beta _m=& \frac{2n_m}{n_i+n_a+2n_m}.\hfill \end{array}$$
(5)
We note that $`\beta _a`$ is the relative amount of protons bound in atoms and $`\beta _m`$ the relative amount bound in molecules. Due to the balance relation for the total proton density

$$n=n_i+n_a+2n_m.$$
(6)
we find the relations

$$\beta _a=(1\alpha )\beta ,\beta _m=(1\alpha )(1\beta ).$$
(7)
The condition of neutrality requires that electron and ion densities are always equal:

$$n_e=n_i.$$
(8)
Therefore the free energy density depends only on 4 independent parameters

$$ff(T,n,\alpha ,\beta ).$$
(9)
If we include magnetic field effects we have 5 independent parameters

$$ff(T,n,B,\alpha ,\beta ).$$
(10)
The equilibrium composition is found by the minimization procedure

$$\frac{\delta f}{\delta \alpha }=0,\frac{\delta f}{\delta \beta }=0.$$
(11)
Instead of $`\beta `$ we may use $`\beta _a`$ or $`\beta _m`$ as variational parameters. The relation

$$\alpha +\beta _a+\beta _m=1$$
(12)
gives the possibility for visualization of the free energy surface on a simplex.

### 2.2 The different contributions to the free energy

We assume the following structure ebetal91 ; beuleal99 :

$$f=f_e^{id}+f_i^{id}+f_a^{id}+f_m^{id}+f_{ee}+f_{ii}+f_{ie}+f_{hs}.$$
(13)
The first 3 terms describe the ideal contributions of the electrons, ions, atoms and molecules. The next 3 terms represent the Coulomb interactions electron-electron, ion-ion, and ion-electron. For these terms different approximations exists as e.g. the quantum Debye-Hückel-approximation (QDHA) and Padé approximation ebkrkr76 . Basically we will use here the QDHA:

$$f_{ee}+f_{ii}+f_{ie}f_{DH}.$$
(14)
In the QDHA screening and quantum effects are both taken into account in a first approximation. The inclusion of further effects as Hartree-Fock etc. by means of Padé approximations makes no principal difficulties but increases the numerical efforts. We will show that in the region of partial ionization and dissociation the QDHA yields a reasonable description of ionization/dissociation phenomena. The interactions of the neutral components are taken into account by hard core approximation including effects of reduced volume ebetal91 ; beuleal99 .

We will discuss now these contributions in more detail:

* the ideal free energy of electrons $`f_e^{id}`$. We use Fermi-Dirac statistics including excluded volume corrections representing the free energy in the following way ebetal91 :

$$f_e=n_ek_BTz\left(\frac{1}{2}\frac{n_e\mathrm{\Lambda }_e^3}{1\eta _0}\right).$$
(15)
The function $`z(x)`$ which is due to Zimmermann ebetal91  interpolates between the low density limit: $`z(x)\mathrm{ln}x1`$ and the high density limit: $`z\left(x\right)x^{\frac{2}{3}}`$. Further $`\mathrm{\Lambda }_e=h/\sqrt{2\pi m_ek_BT}`$ is the thermal wave length of the electrons, correspondingly we will use next the thermal wave lengths of ions, atoms and molecules ($`\mathrm{\Lambda }_i`$, $`\mathrm{\Lambda }_a`$, $`\mathrm{\Lambda }_m`$). The excluded volume factor $`\eta _0<1`$ expresses that some part of the total volume $`V`$ which is occupied by atoms and molecules is not accessible to the electrons. The accessible volume $`(1\eta _0)V`$ is smaller than the total volume and the effective densities are higher. A strict theory for $`\eta _0`$ is not yet available. We will use here

$$\eta _0=\frac{4}{3}\pi R_0^3(n_a+2n_m)$$
(16)
as the definition for the packing parameter with an estimate for the atomic radius $`R_0=0.78\text{Å}`$. Effects of a magnetic field are included by the replacement ebstor00

$$n_e\mathrm{\Lambda }_e^3n_e\mathrm{\Lambda }_e^3\frac{\mathrm{tanh}(x_e)}{x_e},$$
(17)
where $`x_e=\mathrm{}\omega _c^e/2k_BT`$ and $`\omega _c^e=eB/m_e`$.
* the ideal free energy of the ions $`f_i^{id}`$ which is given by the classical term

$$f_i^{id}=n_ik_BT\left[\mathrm{ln}\left(n_i\mathrm{\Lambda }_i^3\right)1\right]n_ik_BT\mathrm{ln}(1\eta _0).$$
(18)
The second term describes the excluded volume effect for the classical bare protons ebetal91 . The parameter $`x_i=\mathrm{}\omega _c^i/2k_BT`$ ($`\omega _c^i=eB/m_i`$) is very small $`x_i1`$ for small $`B`$. Therefore the magnetic field effects can be neglected for protons.
* the ideal atom contribution: $`f_a^{id}`$. We assume Boltzmann statistics i.e.

$$f_a^{id}=n_ak_BT\left[\mathrm{ln}\left(\frac{n_a\mathrm{\Lambda }_a^3}{\sigma _a}\right)1\right]$$
(19)
including the internal states by a Brillouin-Planck-Larkin partition function ebetal91

$$\sigma _a(\xi )=\underset{n=1}{\overset{\mathrm{}}{}}n^2\left\{e^{\left(\frac{\xi }{2}\right)^2\frac{1}{n^2}}1\left(\frac{\xi }{2}\right)^2\frac{1}{n^2}\right\},$$
(20)
where $`\xi =2\sqrt{\mathrm{Ry}/k_BT}`$ with the ionization energy of hydrogen $`\mathrm{Ry}13.598\mathrm{eV}`$. Including the magnetic field we have to replace the atomic partition function $`\sigma _a`$ by an effective partition function,

$$\sigma _a(T)\sigma _{\mathrm{eff}}(T,B),$$
(21)
which reads ebstor00

$$\sigma _{\mathrm{eff}}(T,B)=\sigma _a(\xi _e)+\frac{x_e^2}{24}\left[\sigma _B(\xi _e)+\frac{\xi _e^4}{192}\left(1+\frac{\pi ^2}{3}\right)\right]$$
(22)
with the interaction parameter $`\xi _e=e^2/4\pi ϵ_0k_BT\lambda _e`$ ($`\lambda _e=\mathrm{}/\sqrt{2m_ek_BT}`$ is the thermal wave length for the relative motion of electrons and protons (with an infinit proton mass)) and ebstor00

$`\sigma _B(\xi )`$ $`=`$ $`{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}2n^2(1+n^2)[e^{\left(\frac{\xi }{2}\right)^2\frac{1}{n^2}}1\left({\displaystyle \frac{\xi }{2}}\right)^2{\displaystyle \frac{1}{n^2}}`$ (23)
$`{\displaystyle \frac{1}{2!}}\left({\displaystyle \frac{\xi }{2}}\right)^4{\displaystyle \frac{1}{n^4}}]{\displaystyle }_{n=1}^{\mathrm{}}n^4(5+7n^2)({\displaystyle \frac{2}{\xi }})^2`$
$`\times [e^{\left(\frac{\xi }{2}\right)^2\frac{1}{n^2}}1\left({\displaystyle \frac{\xi }{2}}\right)^2{\displaystyle \frac{1}{n^2}}{\displaystyle \frac{1}{2!}}\left({\displaystyle \frac{\xi }{2}}\right)^4{\displaystyle \frac{1}{n^4}}`$
$`+{\displaystyle \frac{1}{3!}}\left({\displaystyle \frac{\xi }{2}}\right)^6{\displaystyle \frac{1}{n^6}}].`$
* the molecular contribution: $`f_m^{id}`$. We will use the representation beuleal99

$$f_m^{id}=n_mk_BT\left[\mathrm{ln}\left(\frac{n_m\mathrm{\Lambda }_m^3}{\sigma _m}\right)1\right]+f_m^{vib}+f_m^{rot}$$
(24)
including the internal molecular and atomic states,

$$\sigma _m=\sigma _a^2\mathrm{exp}\left(\frac{4.746\mathrm{eV}}{k_BT}\right),$$
(25)
and the contributions of vibrational and rotational states fowler

$`f_m^{vib}`$ $`=`$ $`n_mk_BT\left[\mathrm{ln}\left(1e^{\frac{T_v}{T}}\right)+{\displaystyle \frac{T_\nu }{2T}}\right],`$
$`T_v=6210\mathrm{K}`$
$`f_m^{rot}`$ $`=`$ $`n_mk_BT\left[\mathrm{ln}{\displaystyle \frac{T_r}{T}}{\displaystyle \frac{T_r}{3T}}{\displaystyle \frac{1}{90}}\left({\displaystyle \frac{T_r}{T}}\right)^2\right],`$
$`T_r=85\mathrm{K}.`$
* the Debye-Hückel contribution $`f_{DH}`$. Using here the so-called $`\frac{\mathrm{\Lambda }}{8}`$-approximation ebetal91  we arrive at the formula

$$f_{DH}=k_BT\frac{\kappa ^3}{12\pi }\tau \left(\kappa a(T)\right).$$
(28)
Here the inverse Debye-radius is defined by

$$\kappa =\sqrt{\frac{2n_ie^2}{\epsilon _0k_BT}}$$
(29)
and the effective thermal electron radius by

$$a(T)=\frac{\mathrm{\Lambda }_e}{8}$$
(30)
where $`\mathrm{\Lambda }_e`$ is the electron thermal wave length defined above. Further the $`\tau `$-function has the standard Debye-Hückel form ebkrkr76 :

$$\tau (x)=\frac{3}{x^3}\left[\mathrm{ln}\left(1+x\right)x+\frac{x^2}{2}\right].$$
(31)
The advantage of the QDHA is that magnetic field effects may easily taken into account by a factor ebstor00

$$a(T)a(T,x_e)=a(T)\left[1\frac{1}{48}x_e^2\right].$$
(32)
The magnetic field reduces the effective diameter. This reflects the fact that the magnetic field localizes the particles perpendicular to the field.
* the hard-core contribution $`f_{hs}`$. For this term we use the Carnahan-Starling approximation ebetal91 :

$$f_{hs}=k_BT\frac{\left(n_a+2n_m\right)\left(4\eta 3\eta ^2\right)}{\left(1\eta \right)^2}.$$
(33)
The volume fraction of the neutrals is here defined by the effective packing parameter

$$\eta =\frac{4\pi }{3}R^3\left(n_a+2n_m\right).$$
(34)
Following the work of Juranek and Redmer jure00  we fix the effective packing radius of atoms at a mean value $`R=0.37\text{Å}`$. This is of course a rough approximation which might be refined by taking into account different radii for the atoms and for the molecules and by introducing temperature-dependent radii based on fluid variational theory jure00 .

## 3 Discussion of the approximations and influence of parameters

For this section we neglect the magnetic field, i.e. we assume $`B=0`$. We are interested in partial ionization and dissociation, i.e. in the region of the density-temperature plane where electrons, ions, atoms, and molecules exist. In this region the electrons may show quantum effects but the ions, atoms, and molecules are still classical.

The excluded volume effect, i.e. the effect that the electrons and ions cannot penetrate into the neutral particles, are taken into account in the ideal contributions. The excluded volume factor $`\eta _0`$ is defined above with a fixed radius $`R_0`$. If this effect is not considered we find a lower ionization rate $`\alpha `$ at higher densities, shown in figure 1.

A small difference of the radius $`R_0`$ has an influence on the degree of ionization and dissociation, see figure 1. For larger radii, $`\alpha `$ begins to increase at smaller densities since the effective volume $`V_{\mathrm{eff}}=(1\eta _0)V`$, which is available for the electrons, decreases with increasing $`\eta _0`$ and thus increases the effective electron density $`n_e^{\mathrm{eff}}=n_e/V_{\mathrm{eff}}`$. Therefore the system tends to higher degrees of ionization. The transition between low and complete ionization takes place in a very narrow density range. For high densities the excluded volume effect has to be considered and has still to be refined maybe by taking into account temperature dependencies of the radius jure00  and by introducing occupation probabilities sach91 ; po96 ; chpo98 .

The nonideal effects of the free charged particles can be described more correctly by means of Padé approximation (see e. g. ebri85 ; steb98 ; stbl00 ). Figure 2 shows the ionization and dissociation rate as a function of the proton density calculated by using the QDHA and by the Padé formulae given in steb98 . For the density range of $`10^{22}10^{25}\mathrm{cm}^3`$ and the temperature range $`110^4510^4\mathrm{K}`$ we find only small deviations of the QDHA from the Padé approach. For higher temperatures above $`10^5\mathrm{K}`$ the deviations are much larger and the QDHA is less accurate. In the following we will restrict our study to temperatures much below $`10^5\mathrm{K}`$, therefore it is justified to use further the QDHA to simplify the formulae and to reduce calculation efforts.

## 4 The geometry of the free energy landscape

On the basis of the formulae given in section 2 we represented the geometry with respect to $`f(\alpha ,\beta )`$ at fixed $`T`$, $`n`$ and starting with $`B=0`$ (figures 3, 4). We used the variables $`\alpha `$, $`\beta _a`$, and $`\beta _m`$ which allow a representation on a simplex according to eq. (12). The minimum of the free energy corresponds to the stationary degrees of ionization and dissociation. The form of the landscape around the minimum determines the thermal fluctuations.

Depending on the values of temperature and proton density we find different locations of the minimum. At $`T=300\mathrm{K}`$, $`n=210^{19}\mathrm{cm}^3`$ we observe a minimum in the corner $`\alpha =0`$, $`\beta _a=0`$, $`\beta _m=1`$, meaning that only molecules exist (figure 3). At different temperature and proton density, e.g. $`T=410^4\mathrm{K}`$, $`n=10^{22}\mathrm{cm}^3`$ we observe a minimum in the center of the simplex. This corresponds to a situation where free electrons and ions and bound states (atoms and molecules) exist under the same conditions (figure 4).

With the procedure described above we calculated the stationary degrees of ionization and dissociation for a lattice of points in the temperature-density-plane. In this way we obtained the degree of ionization as a function of density and temperature (figure 5).

We observe full ionization ($`\alpha =1`$ on a plateau in the upper right corner and at the left edge of the diagram) at very high densities and in the region of lower densities but high temperatures. At sufficiently low density we find only bound states when the temperature is low ($`\alpha =0`$ in a region beginning at the lower edge of the diagram). In the intermediate region of proton densities we find a valley of bound states stretching over the whole temperature range with small ionization taking place only at very high temperatures ebkrkr76 .

Now we will look at the formation of molecules (figure 6). We observe similar features as above in the corresponding parameter regions. At very high densities (left edge) no molecules exist as well as in the region of lower densities and high temperatures (upper right corner). When the density is sufficiently low and the temperature is low as well we find that the bound states observed in figure 5 are practically only made up of molecules ($`\beta _m1`$ at the lower edge). In the intermediate density region we find a mountain of molecules ($`\beta _m=1`$) and a decline of the degree of dissociation with rising temperature.

## 5 Fluctuations according to Boltzmann-Gibbs

Let us first consider the general schema of the theory of isothermal fluctuations at given temperature $`T`$ ll ; klim95 . We follow Klimontovich who developed the general fluctuation theory of Boltzmann and Gibbs klim95 . In the general case we study a set of intrinsic parameters $`A(X,a)=(A_1(X,a),\mathrm{},A_n(X,a))`$ depending on dynamical variables $`X`$ and external parameters $`a`$ and define a conditional free energy $`F(a,T|A)`$ which is the free energy at fixed values $`A=A(X,a)`$. The intrinsic parameters $`A(x,a)`$ do not need to have any thermodynamic functions which correspond to them. The equilibrium value of the free energy is $`F(a,T)`$. Then according to the Boltzmann-Gibbs principle the probabilities of fluctuations of $`A`$ are given by

$$p(A|a,T)=\mathrm{exp}\left[\frac{F(a,T)F(a,T|A)}{k_BT}\right].$$
(35)
In our case the dynamical variables are $`A=\{\alpha ,\beta \}`$, for $`F`$ we choose the free energy per proton $`\varphi =f/n`$ with the equilibrium value $`\varphi _0=f_0/n`$ at fixed temperature $`T`$.

Then according to the Boltzmann-Gibbs principle the probability distribution of isothermal fluctuations reads

$$p(\alpha ,\beta |n,T)=\frac{e^{\frac{\varphi _0\varphi (\alpha ,\beta )}{k_BT}}}{_{\alpha ,\beta }e^{\frac{\varphi _0\varphi (\alpha ,\beta )}{k_BT}}}.$$
(36)
An illustration of the probability distribution

calculated for a lattice of 50x50 points on the $`\alpha ,\beta `$plane is shown in figure 7.

Further we show in figure 5 the degree of ionization together with the dispersion of the degree of ionization coded in gray levels. The diagram was obtained with the following procedure: For a given temperature and proton density the stationary degrees of ionization and dissociation ($`\overline{\alpha },\overline{\beta }`$) were calculated. Then we additionally fixed the degree of dissociation at its stationary value $`\beta =\overline{\beta }`$ and calculated the probability $`p(\alpha )`$ to find arbitrary values for the degree of ionization. The dispersion of the degree of ionization $`\overline{(\delta \alpha )^2}=\overline{(\alpha \overline{\alpha })^2}`$ was then calculated as the root of the mean quadratic deviation from the stationary value $`\overline{\alpha }`$.

We observe low dispersion (light) at very high densities due to the hard-sphere contribution. High dispersion (dark) is observed in two transition regions, one when ionization processes start taking place at intermediate densities and temperatures (just above $`\alpha =0`$), the other at the transition to complete ionization (just below $`\alpha =1`$) where the remaining bound states break up.

Inside the simplex we assume that the probability distribution function $`p(\alpha ,\beta =\overline{\beta })`$ is an one dimensional gaussian distribution for values $`\alpha `$ inside the borders and for small fluctuations $`\delta \alpha =\alpha \overline{\alpha }`$

$$p(\alpha )=C\mathrm{exp}\left[\frac{1}{2k_BT}\left(\frac{^2\varphi }{\alpha ^2}\right)_{\alpha =\overline{\alpha }}\left(\alpha \overline{\alpha }\right)^2\right]$$
(37)
where $`\overline{\alpha }`$ is the equilibrium value and therefore the most probable value. The constant $`C`$ is determined from the normalization. The dispersion in the fluctuations is given by

$$\overline{(\delta \alpha )^2}=k_BT\left[\left(\frac{^2\varphi }{\alpha ^2}\right)_{\alpha =\overline{\alpha }}\right]^1.$$
(38)
We see that the second deviation of $`\varphi (\alpha ,\beta =\overline{\beta })`$, which is the curvation of the free energy landscape determines the dispersion $`\overline{(\delta \alpha )^2}`$. Near to the corner of the simplex the distribution is highly non-gaussian and highly unsymmetrical. In this case the border of frequent fluctuations (the range of the dispersion) may be estimated by finding the roots of $`p(\alpha )p(\alpha =\overline{\alpha })/\sqrt{e}`$. For gaussian distributions this criterion agrees with the one given above.

In density-temperature regions where the minimum of the free energy is developed very well e.g. if the minimum is located in a corner of the simplex (see e.g. figure 3) we have small dispersions. Is the free energy landscape flat e.g. if the minimum is located inside the simplex (see e.g. figure 4) we have large dispersions.

To further illustrate this, we show the dispersion of the degree of ionization at fixed temperature, corresponding to a cut through the diagram in figure 5 parallel to the $`n`$-axis. The dispersion is indicated by error bars, which are cut at the physically impossible values $`\alpha >1`$ and $`\alpha <0`$ (figure 8). Again the dispersion is small at high and low densities and is high near the transition $`\alpha =1\alpha <1`$ where bound states emerge.

## 6 Magnetic field effects

We want to discuss now the influence of a constant uniform magnetic field on the ionization and dissociation equilibrium. We have to carry out the replacements (17), (21), and (32). These formulae equations are valid for weak magnetic fields, i.e.

$$\frac{x_e^2}{48}1.$$
(39)
In figure 9 the degree of ionization and dissociation for various magnetic field strengths are plotted over the temperature. We find at a fixed density that the degree of ionization of a plasma in a magnetic field is higher compared to the field free case and increases with the field strength. For temperatures higher than $`10^5\mathrm{K}`$ and weak magnetic fields there are no differences in the ionization degree. Further we see in figure 10 that the curves deviate in the mid-density range where atoms and molecules exist up to the density where the transition to complete ionization takes place.

## 7 Conclusions

Based on a detailed study of the free energy landscape we developed a new description of ionization and dissociation phenomena. The relatively simple approximation QDHA used for the interaction term in the free energy allows the study of the dispersion of the degrees of ionization and dissociation and of magnetic field effects. The QDHA may be replaced by more refined approximations as Padé approximation. We have shown that the study of free energy landscapes in combination with the fluctuation theory yields relevant information on dispersion of the degrees of ionization and dissociation. The natural dispersion is usually in the range of $`510\%`$ except on the borders. This has consequences for many experimental quantities as e.g. conductivity and optical properties. Since the conductivity is in first approximation proportional to $`\alpha `$, i.e. $`\sigma \alpha \sigma _0`$, we may conclude that the thermal dispersion of the plasma conductivity $`\sigma `$ has nearly the same shape as the dispersion of $`\alpha `$. We note that the thermal fluctuations of $`\alpha `$, $`\beta `$, $`\sigma `$ etc. are in general not of Gaussian character what is due to the existence of borders as e.g. $`0\alpha 1`$. Studying the influence of magnetic fields we have shown that a remarkable influence on the degree of ionization starts at $`B10^4\mathrm{T}`$, the effects increases with the magnetic field strength. At fields $`B510^4\mathrm{T}`$ the degree of ionization and also the conductivity and other related quantities are substantially higher than for the field free case, the present approach however is no more valid at these extremely high fields.
Recent experimental work on dense hydrogen concentrates mainly on Hugoniots knetal01 . The approach presented here is restricted to temperatures $`T>10000K`$, therefore it cannot be used for the calculation of the low-temperature parts of Hugoniots. An extension taking into account recent results on the low temperature EOS of dense hydrogen ebetal91  is in progress.

###### Acknowledgements.

The authors thank Stefan Hilbert, Hauke Juranek, Ronald Redmer, Gerd Röpke, Manfred Schlanges and Werner Stolzmann for many helpful discussions. The authors W. Ebeling and H. Hache profited very much from a stay at the Rostock University.