text:
# Strong and electromagnetic interaction views on deuteron structure.

## I Introduction

The deuteron (bound state of proton and neutron) with spin $`𝒥`$ and parity $`P`$: $`𝒥^P=1^+`$, and with isospin $`I`$=0, has been investigated, theoretically and experimentally, since many decades. Few selected problems, which are very actual, will be discussed here according the following plan:

* The knowledge of the deuteron structure, itself, at the largest possible internal momenta (or at the shortest distances between the nucleons). The question which is addressed is the definition of the kinematical region where a description based on nucleons and mesons (basically impulse approximation (IA) with corrections due to meson exchange currents, isobars, relativistic effects..) has to be discarded in favor of a ’high energy view’ in terms of quarks and gluons. Here, in particular, pQCD gives predictions concerning the asymptotic behavior of the deuteron form factors and assumes hadron helicity conservation for the amplitudes.

The transition region between these two regimes should be the privileged domain of intermediate energy machines. High intensity is required by exclusive measurements and/or studies with secondary beams (polarization etc..).
* The deuteron as a probe to investigate the nucleon or the heavy nuclei structure.

+ it can be used to study the properties of the neutron, which is not available as a target;
    + it is an isoscalar probe, which can be very selective in exciting specific states in nucleons and nuclei.

In order to illustrate the different points listed above, we will review some of the last data obtained at Jefferson Laboratory and Saturne. We will show few results from elastic electron-deuteron scattering and proton-deuteron elastic and inelastic scattering. Particular attention will be devoted to polarization observables.

## II The deuteron electromagnetic form factors

### A Recent determination of the elastic deuteron electromagnetic form factors

The measurement of the differential cross section of elastic $`ed`$scattering, for a fixed value of $`Q^2`$, at different scattering angles, allows to determine the structure functions $`A(Q^2)`$ and $`B(Q^2)`$:

$$\frac{d\sigma }{d\mathrm{\Omega }}=\left(\frac{d\sigma }{d\mathrm{\Omega }}\right)_0𝒮,𝒮=A(Q^2)+B(Q^2)\mathrm{tan}^2(\theta _e/2)$$

with

$$\left(\frac{d\sigma }{d\mathrm{\Omega }}\right)_0=\frac{\alpha ^2cos^2(\theta _e/2)E^{}}{4E^3sin^4(\theta _e/2)},$$

where $`E`$ ($`E^{}`$) is the electron beam (the scattered electron) energy and $`\theta _e`$ the electron scattering angle in the Laboratory system. The structure functions $`A`$ and $`B`$ can be expressed in terms of the three form factors, $`G_c`$ (electric), $`G_m`$ (magnetic) and $`G_q`$ (quadrupole) as:

$$A(Q^2)=G_c^2(Q^2)+\frac{8}{9}\tau ^2G_Q^2(Q^2)+\frac{2}{3}\tau G_m^2(Q^2),B(Q^2)=\frac{4}{3}(1+\tau )\tau G_m^2(Q^2),\tau =\frac{Q^2}{4M^2},$$

where M is the deuteron mass. In case of unpolarized beam and target the outgoing deuteron is tensorially polarized and the components of the tensor polarization give useful combinations of form factors. In particular $`t_{20}`$ allows, together with $`A(Q^2)`$ and $`B(Q^2)`$, the determination of the three form factors:

$$t_{20}=\frac{1}{\sqrt{2}𝒮}\left[\frac{8}{3}\tau G_cG_q+\frac{8}{9}\tau ^2G_Q^2+\frac{1}{3}\tau \left(1+2(1+\tau )\mathrm{tan}^2(\theta _e/2)\right)\right]G_m^2,$$

At the Jefferson Laboratory (JLab), the elastic $`ed`$ cross section has been recently precisely measured up to large momentum transfer $`Q^26`$ (GeV/c)-2,  and $`t_{20}`$ has been measured for a momentum transfer up to $`Q^2=1.9`$ (GeV/c)-2 .

According to , the cross sections seems to scale as $`(Q^2)^{10}`$, as previously pointed out , and predicted by pQCD. However, from the $`t_{20}`$ data, it clearly appears that the pQCD limit is not yet reached, and that the data follow the trend suggested by IA. On the other hand, it is not possible, from these data, to constrain definitely different models or corrections. For a detailed comparison with theory, see, for example, .

The question is then, how to proceed further. On the experimental point of view, in the best presently achievable conditions of luminosity, about one month beam time was needed for the the polarization measurement and one week for the cross section at the largest values of $`Q^2`$. It seems hard to foresee more favorable experimental conditions, with the present technology related to electron machines. One possibility is to investigate the inelastic deuteron form factors through the reaction $`eded\pi ^0`$, where for the same momentum transfer for the electron, one acceeds to shortest distances in the deuteron, compared to elastic scattering .

In next section we will focus on a conceptual problem, which is the limit of the validity of the assumption of the one photon exchange mechanism in electron-hadron scattering at large momentum transfer. A related problem, that we will not discuss here, is that radiative corrections for polarization effects at large momentum transfer are not known.

### B Beyond the one-photon approximation

The formulas given above are valid if the momentum is transferred from the incident electron to the target by a virtual photon, with the underlying assumption that the possible two-photon contribution is small. The relative contribution of two photon exchange, from simple counting in $`\alpha `$, would be of the order of the fine structure constant, $`\alpha ={\displaystyle \frac{e^2}{4\pi }}{\displaystyle \frac{1}{137}}`$. However, more than 25 years ago it was observed that the relative role of two-photon exchange can increase significantly in the region of high momentum transfer . If the transferred momentum is equally shared between two virtual photons, due to the steep decrease of the deuteron form factors, the simple rule of $`\alpha `$-counting for the estimation of the relative role of two-photon contribution to the amplitude of elastic $`ed`$scattering does not hold anymore. This effect would manifest already at momentum transfer of the order of 1 GeV2, in particular in the region of diffractive minima.

In Ref.  the two-photon amplitude is purely imaginary, at least at small scattering angles, so it cannot interfere with the one-photon exchange amplitude in the differential cross section for unpolarized particles scattering. However, in this case, the polarization observables in elastic $`ed`$scattering have to be large, in particular the T-odd polarization observables. But the predicted increasing of the two-photon mechanism is so large that it may be observed even in the differential cross section of elastic $`ed`$scattering, at relatively large momentum transfer square, $`Q=812`$ fm-1. An evaluation of this contribution from the existing experimental data has been done in .

The crossing symmetry can provide a relation between the matrix elements $``$ of the the elastic $`e^{}+he^{}+h`$ scattering and the $`e^+e^{}`$-annihilation: $`e^++e^{}\overline{h}+h`$, in one-photon approximation.

$$\overline{|(eheh)|^2}=f(s,t)=\overline{|(e^+e^{}\overline{h}h)|^2}.$$
(1)
The line over $``$ denotes the sum over the polarizations of all particles (in initial and final states). The Mandelstam variable $`s`$ is the total energy square and $`t`$ is the momentum transfer square. They delimit different kinematical regions for the annihilation and the scattering channel.

The presence of a single virtual photon in the reaction $`e^++e^{}\gamma ^{}\overline{h}+h`$ constrains the total angular momentum $`𝒥`$ and the $`P`$-parity for the $`\overline{h}h`$system to take only one possible value, $`𝒥^P=1^{}`$, the quantum number of the photon. In the framework of the one-photon approximation, in the general case, $`\overline{|(e^+e^{}h\overline{h})|^2}`$ can be written (in CMS) as:

$$\overline{|(e^+e^{}h\overline{h})|^2}=a(t)+b(t)\mathrm{cos}^2\stackrel{~}{\theta },\text{and }\mathrm{cos}^2\stackrel{~}{\theta }=1+\frac{\mathrm{cot}^2\frac{\theta _e}{2}}{1+\tau },$$
(2)
where $`a(t)`$ and $`b(t)`$ are definite quadratic combinations of the electromagnetic form factors for the hadron $`h`$ and $`\stackrel{~}{\theta }`$ is the angle of the detected hadron.

In case of the presence of $`2\gamma `$ in the intermediate state, in the annihilation channel, any value of the total angular momentum and space parity is allowed, because the relative 3-momentum for the $`2\gamma `$-state is nonzero, contrary to the case of the one-photon mechanism. The $`\overline{h}h`$-system, produced through $`1\gamma `$\- and $`2\gamma `$-exchanges has different values of C-parity, because $`C(\gamma )=1`$ and $`C(2\gamma )=+1`$. Therefore the interference of one- and two-photon contribution must be an odd function of $`\mathrm{cos}\stackrel{~}{\theta }`$: $`\overline{e_1_2^{}}=\mathrm{cos}\stackrel{~}{\theta }(a_0+a_2\mathrm{cos}^2\stackrel{~}{\theta }+\mathrm{})`$.

The SLAC results  where obtained at a fixed electron scattering angle. The odd contribution can then be estimated from the cross sections at two angles, for the same $`Q`$, one value being given from a fit of the data of ref.  and the other by the recent data . We can then calculate separately the linear contribution or the cubic contribution in $`\mathrm{cos}\stackrel{~}{\theta }`$. It is expected that deviation from the linear $`\mathrm{cot}^2\frac{\theta _e}{2}`$ formula would not appear for $`Q`$ 5 fm-1. The resulting ratios C/A and D/A are reported in Fig. 1 and 2 as functions of $`Q`$, as open circles  and open squares . The data nicely agree: at large momentum transfer these ratios deviate from zero and show a dependence on the transferred momentum, which could result from two photon exchange. However the points had to be rescaled in order to have zero deviation at low $`Q`$ (corresponding solid symbols) due to systematic errors in the measurement of the cross section between the two sets of data  and .

While this can not be considered as a definite evidence for the presence of $`2\gamma `$-exchange in $`ed`$elastic scattering, it is the first attempt to obtain a quantitative upper limit of a possible $`2\gamma `$-contribution, using a parameterization of the $`1\gamma 2\gamma `$-interference and the existing experimental data.

The $`2\gamma `$-exchange in elastic hadron scattering can be experimentally searched in different ways: - through the comparison of the cross section for scattering of unpolarized electrons and positrons (by protons or deuterons) in the same kinematical conditions, - looking to the deviation from a straight line on the Rosenbluth plot or measuring specific properties of polarization phenomena: as nonzero $`Todd`$ polarization observables, and violation of definite relations between T-even polarization observables and the SF $`B(Q^2)`$. The measurement of cross section and polarization observables bring complementary and independent pieces of information, as they test the real and imaginary part of the $`2\gamma `$ contribution.

## III The deuteron structure from $`dp`$ break-up and backward elastic scattering

We showed that the measurement of cross section and analyzing powers in $`ed`$elastic scattering allows to determine all the three form factors (in one-photon approximation). These form factors, in IA, are integrals of the wave function over the radial coordinate. On the other hand, the reactions $`\stackrel{}{d}+p\stackrel{}{p}+d`$ (backward elastic scattering) and $`\stackrel{}{d}+pp+X`$ (deuteron break up) are directly related to the deuteron wave functions. In the IA, the tensor analyzing power $`T_{20}`$ and the polarization transfer $`\kappa _0`$ can be written as:

$$T_{20}=\frac{1}{\sqrt{2}}\frac{\sqrt{8}uww^2}{u^2+w^2},$$

$$\kappa _0=\frac{u^2w^212uw}{u^2+w^2},$$

from which one obtains a quadratic relation between $`T_{20}`$ and $`\kappa _0`$ :

$$\left(T_{20}+\frac{1}{2\sqrt{2}}\right)^2+\kappa _0^2=\frac{9}{8}$$

Measurements have been performed at Saturne and Dubna . The correlation between the two polarization observables is shown in Fig. 3. The full line is the IA prediction. The two sets of points correspond to backward elastic scattering (open circles) and deuteron break up (solid circles) and they show a very similar behavior. The deviation from IA, which gets larger at larger momenta has been interpreted in different models (see  and refs. herein). At this moment it seems fair to conclude that unambiguous signatures of quarks have not yet been found.

A promising way of looking to deuteron at very short distances seems to be the measurement of the tensor analyzing power of pions emitted at $`0^0`$, in the reaction $`\stackrel{}{d}+p\pi +X`$, in the cumulative region .

## IV  The deuteron as an isoscalar probe

### A The neutron electromagnetic form factors

Having high precision data on the differential cross section for $`ed`$ elastic scattering, and assuming a reliable model for their description, one can extract, in principle, the dependence of the electric neutron form factor $`G_{En}`$ on the momentum transfer $`Q^2`$. Such a procedure has been carried out in ref. , up to $`Q^2`$=0.7 (GeV/c)2. It can be extended at higher $`Q^2`$ using the elastic $`ed`$-scattering data mentioned above and recent data on the proton electric form factor . These data have been obtain by the recoil proton polarization measurement in $`\stackrel{}{e}+pe+\stackrel{}{p}`$, following an idea suggested more than 30 years ago  and extend up to $`Q^2`$=3.5 (GeV/c)2. The large sensitivity to the nucleon form factors of the models which describe the light nuclei structure, particularly the deuteron, was already carefully studied in , and it was pointed out that the disagreement between the relativistic impulse approximation and the data could be significantly reduced if $`G_{En}`$ were different from zero.

In the non relativistic IA, the deuteron form factors depend only on the deuteron wave function and on nucleon form factors:

$$G_c=G_{Es}C_E,G_q=G_{Es}C_Q,G_m=\frac{M_d}{M_p}\left(G_{Ms}C_S+\frac{1}{2}G_{Es}C_L\right),$$
(3)
where $`M_p`$ is the proton mass, $`G_{Es}`$=$`G_{Ep}`$+$`G_{En}`$ and $`G_{Ms}`$=$`G_{Mp}`$+$`G_{Mn}`$ are the charge and magnetic isoscalar nucleon form factors, respectively. The terms $`C_E`$, $`C_Q`$, $`C_S`$, and $`C_L`$ describe the deuteron structure and can be calculated from the deuteron $`S`$ and $`D`$ wave functions, $`u(r)`$ and $`w(r)`$  :

$$C_E=_0^{\mathrm{}}𝑑rj_0\left(\frac{Qr}{2}\right)\left[u^2\left(r\right)+w^2(r)\right],$$

$$C_Q=\frac{3}{\sqrt{2}\eta }_0^{\mathrm{}}𝑑rj_2\left(\frac{Qr}{2}\right)\left[u(r)\frac{w(r)}{\sqrt{8}}\right]w(r),$$

$$C_S=_0^{\mathrm{}}𝑑r\left[u^2(r)\frac{1}{2}w^2(r)\right]j_0\left(\frac{Qr}{2}\right)+\frac{1}{2}\left[\sqrt{2}u(r)w(r)+w^2(r)\right]j_2\left(\frac{Qr}{2}\right),$$
(4)
$$C_L=\frac{3}{2}_0^{\mathrm{}}𝑑rw^2(r)\left[j_0\left(\frac{Qr}{2}\right)+j_2\left(\frac{Qr}{2}\right)\right],$$

where $`j_0(x)`$ and $`j_2(x)`$ are the spherical Bessel functions. The normalization condition is $`_0^{\mathrm{}}𝑑r\left[u^2(r)+w^2(r)\right]=1.`$

With the help of expressions (3) and (4), the formula for $`A(Q)^2`$, can be inverted into a quadratic equation for $`G_{Es}`$. Then $`G_{Es}`$ is calculated using the experimental values for $`A(Q)^2`$, assuming, for the magnetic nucleon form factors $`G_{Mp}`$ and $`G_{Mn}`$ the usual dipole dependence, which is in agreement with the existing data at a 3% level, up to $`Q^2`$ 10 (GeV/c)2.

In Fig. 4 we illustrate the behavior of the different nucleon electric form factors: $`G_{Es}`$$`G_{Ep}`$ and $`G_{En}`$. The nucleon isoscalar electric form factor, derived from different sets of deuteron data, decreases when $`Q^2`$ increases. The solid line represents the Gari-Krümpelmann parametrization  for $`G_{Es}`$. The dipole behavior, which is generally assumed for the proton electric form factor is shown as a dotted line. The new $`G_{Ep}`$ data, which decrease faster than the dipole function, are also well reproduced by the Gari-Krümpelmann parametrization (thick dashed line).

The electric neutron form factor can be calculated from the isoscalar nucleon form factor, taking for $`G_{Ep}`$ a dipole behavior (solid stars) or a fit based on the new data (open stars). The last option leads to values for $`G_{En}`$ which are in very good agreement with the parametrization . These results shows that the neutron form factor is not going to vanish identically at large momentum transfer, but becomes more sizeable than predicted by other parametrizations, often used in the calculations  (thin dashed line). Starting from $`Q^22`$ (GeV/c)2 the form factor $`G_{En}`$ becomes even larger than $`G_{Ep}`$. Let us mention that a recent ’direct’ measurement  at $`Q^2=0.67`$ (GeV/c)2 finds $`G_{En}`$=$`0.052\pm 0.011\pm 0.005`$ in agreement with the present values.

Let us mention that the $`\gamma ^{}\pi ^\pm \rho ^{}`$-contribution, which is a good approximation for the isoscalar transition $`\gamma ^{}\pi ^+\pi ^{}\pi ^0`$ ($`\gamma ^{}`$ is a virtual photon), is typically considered as the main correction to IA, necessary, in particular, to improve the description of the SF $`A(Q^2)`$ . However the relative role of MEC is strongly model dependent  as the coupling constants for meson-NN-vertexes are not well known and arbitrary form factors are often added .

It should be pointed out that the $`\gamma ^{}\pi \rho `$ vertex is of magnetic nature and its contribution to $`A(Q^2)`$ has to be of the same order of magnitude as the relativistic corrections.

The forthcoming data about $`G_{En}`$, planned at JLab up to $`Q^2`$=2 (GeV/c)2,  will be crucial in this respect. The large sensitivity of the deuteron structure to the nucleon form factors shows the necessity to reconsider the role of meson exchange currents in the deuteron physics at large momentum transfer.

### B The nuclear structure

The measurement of polarization transfer for inelastic $`\stackrel{}{p}`$ and $`\stackrel{}{d}`$ scattering by nuclei allows the study of the nuclear response and to disentangle the spin and isospin components. In particular it has been shown that he spin-flip probability is a very good signature of the presence of $`\mathrm{\Delta }S=1`$ in the continuum, as well as for discrete states. A systematic work has been carried on at LAMPF and at SATURNE, on several nuclei, from $`{}_{}{}^{12}C`$ to $`{}_{}{}^{208}Pb`$. The nuclear response has been measured in different channels and compared to RPA calculations, in order to learn about correlations and collectivity in nuclei .

In order to go further, one should complete measurements at very small angles (with a zero degree facility) and with high statistics, and make a multipole analysis in the continuum.

### C The baryon resonances excitation in the inclusive $`\stackrel{}{d},p`$ scattering

An application of this method based on the measurement of the spin-flip probability, is the study of nucleon resonances, in the reaction $`\stackrel{}{d}+p\stackrel{}{d}+X`$. The selectivity of reactions such as $`p(d,d^{})X`$ or $`p(\alpha ,\alpha ^{})X`$  to the isoscalar part of the $`N^{}`$-electroexcitation makes these processes complementary to electron-nucleon inelastic scattering, for the study of the $`N^{}`$-structure. Particular attention has been devoted to the Roper resonance  In case of polarized deuteron beam, it has been shown  that the $`\omega `$-exchange model gives a natural and simple description of the polarization phenomena for $`\stackrel{}{d}+pd+X`$. The main ingredients of such model are the existing information about the deuteron electromagnetic form factors  and the ratio $`r`$ of the longitudinal and transversal isoscalar cross sections for the excitation of the $`N^{}`$-resonances .

The $`\omega `$-meson is preferred, among the isoscalar mesons as $`\sigma `$ or $`\eta `$, for several reasons. The $`\omega NN`$ coupling is large; the $`\omega `$-meson, being a spin 1 particle, can induce strong polarization effects and an energy-independent cross section. When it is considered as an $`isoscalarphoton`$, then the cross sections and the polarization observables can be calculated from the known electromagnetic properties of the deuteron and $`N^{}`$, through the vector dominance model.

The tensor analyzing power in $`d+pd+X`$, $`T_{20}`$, can be written in terms of the deuteron electromagnetic form factors as:

$$T_{20}=\sqrt{2}\frac{V_1^2+(2V_0V_2+V_2^2)r(t)}{4V_1^2+(3V_0^2+V_2^2+2V_0V_2)r(t)},$$
(5)
where $`V_0(t)`$, $`V_1(t)`$ and $`V_2(t)`$ are linear combinations of the standard electric, $`G_c`$, magnetic $`G_m`$ and quadrupole $`G_q`$ deuteron form factors. The ratio $`r`$ characterizes the relative role of longitudinal and transversal isoscalar excitations in the transition $`\omega +NX`$.

From Eq. (5) one can see that all information about the $`\omega NN^{}`$-vertex is contained in the function $`r`$ only. A zero value of $`r`$ results in a $`t`$ and $`w`$-independent value for $`T_{20}`$, namely $`T_{20}=1/\sqrt{8}`$, for any value of the deuteron electromagnetic form factors. The ratio $`r`$ is calculated using the collective string model in , assuming $`SU_{sf}(6)`$ symmetry, including the contributions of the following resonances: $`N_{11}(1440),S_{11}(1535),D_{13}(1520)\text{ and}S_{11}(1650)`$, which are overlapping in this energy region.

In Fig. 5 we report the theoretical predictions for $`T_{20}`$, in framework of $`\omega `$-model, together with the existing experimental data. In such approximation $`T_{20}`$ is a universal function of $`t`$ only, without any dependence on the initial deuteron momentum. The experimental values of $`T_{20}`$ for $`p(\stackrel{}{d},d)X`$ , for different momenta of the incident beam are shown as open symbols. These data show a scaling as a function of $`t`$, with a small dependence on the incident momentum, in the interval 3.7-9 GeV/c. On the same plot the data for the elastic scattering process $`e^{}+de^{}+d`$  are shown (filled stars).

These different data show a very similar behavior: negative values, with a minimum in the region $`|t|0.35GeV^2`$ and their value increase toward zero at larger $`|t|`$. The lines are the result of the $`\omega `$-exchange model for the $`d+pd+X`$ process: for $`r=0`$ (dashed-dotted line), the calculation based on  for the Roper excitation only is represented by the dotted line and for the excitation of all the four resonances by the the full line. The deuteron electromagnetic form factors have been taken from , a calculation based on relativistic impulse approximation, and they reproduce well the $`T_{20}`$data for $`ed`$ elastic scattering . When $`r0`$ or if the contribution of the deuteron magnetic form factor $`V_1(t)`$ is neglected, then $`T_{20}`$ does not depend on the ratio $`r`$, and coincides with $`t_{20}`$ for the elastic $`ed`$-scattering (with the same approximation).

From Fig. 5 it appears that the $`t`$behavior of $`T_{20}`$ is very sensitive to the value of $`r`$ especially at relatively small $`r`$, $`r0.5`$. The values of $`r`$, predicted by model , give a very good description of the data, when taking into account the contribution of all four resonances. These data, in any case, exclude a very small value of $`r,r0.1`$ as well as very large values of $`r`$. Only the Roper resonance has a nonzero isoscalar longitudinal form factor. Without excitation of the Roper resonance, $`r=0`$, (in the considered kinematical region) and the value for $`T_{20}`$ becomes $`t`$independent: $`T_{20}=1/2\sqrt{2}`$, in evident disagreement with existing data.

Due to their specific quark structure, the resonances lying in the concerned mass region, such as $`S_{11}(1535)`$, $`D_{13}(1520)`$ and $`S_{11}(1650)`$, are characterized by a pure isovector nature of longitudinal virtual photons absorbed by the nucleons. The isoscalar longitudinal amplitudes of $`S_{11}(1535)`$ and $`D_{13}(1520)`$ electroexcitation vanish due to a specific spin-flavor symmetry, while both isoscalar and isovector longitudinal couplings of $`S_{11}(1650),D_{15}(1675)`$ and $`D_{13}(1700)`$ vanish identically.

This behavior of the isoscalar form factors is essential for the correct description of the existing experimental data on the $`t`$dependence of $`T_{20}`$ for the process $`d+pd+X`$.

One could improve the model taking into account for example, other meson exchanges, or the effects of the strong interaction in initial and final states. However these corrections are strongly model- and parameter- dependent and the existing experimental data are not good enough to constrain the additional parameters which have to be added.

The successful description of the polarization observable $`T_{20}`$ can be considered as a strong indication that the $`\omega `$ exchange is the main mechanism for the considered process and that the Roper resonance is excited in this process.

## V Conclusions

Some of the recent results obtained in the study of the deuteron structure have been reviewed. Despite the large precision and the new region of internal momentum explored, questions arised many decades ago are still actual:

* The relative role of different possible components in the deuteron wave function: isobar configurations, six-quark components, etc..;
* The relativistic description of the deuteron structure: the number of independent components of the deuteron wave functions as well as the number and the nature of their arguments;
* The relative role of the different possible mechanisms in the simplest reactions with deuterons: for example the role of meson exchange current in the description of the electromagnetic form factors of the deuteron;
* The kinematical region corresponding to the transition regime in the deuteron structure: i.e. from the nucleon-meson picture of the deuteron to the quark-gluon description.

These problems are far from being solved. The experimental study of polarization phenomena remains the most fruitful way in this investigation. If, in the case of elastic $`ed`$-scattering we are near the limits of the experimental possibilities, in reactions such as the deuteron photodisintegration, $`\gamma +dn+p`$ or the coherent pion photoproduction on on the deuteron $`\gamma +dd+\pi ^0`$, a wide program of polarization experiments can be realized, due to the rich spin structure of the corresponding matrix elements. Moreover the experimental study of new reactions, such as the coherent neutral pion production on the deuteron, $`e+de+d+\pi ^0`$, in the near threshold region and in the $`\mathrm{\Delta }`$-region (for large momentum transfer) can bring useful information.

These conclusions apply to the hadronic sector, too. The complete experiment for the simplest case, the $`dp`$ backward elastic scattering, requires further steps which can be realized in principle in the nearest future. The study of the nucleon resonances through reactions induced by isoscalar probes as deuterons and $`\alpha `$ particles seems to be very promising, in particular the polarization observables.


format: text

text:
# A Simple Model for Plastic Dynamics of a Disordered Flux Line Lattice

## I Introduction

Collective transport in disordered media is a widespread and poorly understood phenomena. A great deal of experimental and theoretical effort in this area has been devoted to studying the nonlinear dynamics of the disordered flux line lattice (FLL) in type II superconductors. The FLL exhibits a threshold behavior due to the competition between pinning and flux line repulsion . In response to a force, such as that associated with a transport current, the three dimensional FLL can move smoothly via elastic deformations, maintaining its integrity and order. However, in another regime the FLL deforms plasticly. In that regime disorder becomes more important, and the moving FLL manifold tears as some flux lines (in two-dimensions, vortices) move while others do not. As a result, the flow pattern breaks up in a nonuniform way. It is generally believed that the microscopic structure of the FLL, or defects in it, is fundamental to transport behavior both in the plastic and elastic regimes.

Part of the attention to the dynamics of a moving FLL has been motivated by interests in possible, exotic phase transitions and glassy phases, melting and other complicated scenarios associated with structural order in the FLL. Most experiments, however, are essentially transport studies, and are, quoting Higgins and Bhattacharya, “notoriously ill-suited for the study of thermodynamic phase transitions. These experiments yield direct information only about the mobility of flux lines, i.e. on dynamics and pinning, which would then have to be connected, through highly model-dependent ways, to the structure of the state they pertain to.”

Here, we show that many of the empirical results found in transport studies of plastic flux flow in type II superconductors may be obtained with an extremely simple model. It has been recognized for many years that, in the presence of pinning, magnetic flux in type II superconductors forms a pile with an overall, global slope, akin to a sandpile. Penetration of magnetic flux into superconductors driven solely by the flux density gradient has been described using molecular dynamics (MD) simulations , and by the model used here . However, the flux gradient has not been taken into account in any previous numerical simulation studies of the current-voltage characteristic. One possible reason is that previous numerical studies of flux motion at the scale of the vortex cores have not been able to reach a sufficiently large system size. The Bassler-Paczuski (BP) model , on the other hand, is a coarse-grained model and describes the magnetic flux dynamics at the much larger scale of the London length, making the large system size limit much more accessible.

The numerical simulations presented here of the BP model  show that nonlinear behaviors, characteristic of experimental transport measurements of plastic flow in superconductors, arise as a result of vortex flow down a vortex density gradient. The effect of the transport current is modelled by a shift in boundary conditions, which leads to a generalized “tilt” of the vortex pile. Eventually the “tilt” is sufficiently great so that some vortices can flow down the pile in the steady state, leading to the onset of a finite voltage. The vortex flow forms a variety of river morphologies depending on the interaction strength between the vortices, compared to pinning, and the overall rate of flow. The flow patterns of magnetic flux are self-organized together with the magnetic flux profile, which is the substrate on which the flow takes place. Since the BP model does not contain any detailed information on the positions of the vortex cores at the micro-scale of the FLL, it cannot exhibit any structural ordering or disordering behavior.

Although the BP model can also be studied at finite temperatures, here we use the zero temperature limit where thermal fluctuations of the flux motion may be ignored. The zero temperature approximation seems reasonable to describe the plastic transport dynamics of the low temperature superconductors and, perhaps, some aspects of the high temperature superconductors as well. Also, we consider the limit where the depairing current density, $`j_0`$, is extremely large compared to the critical current density, $`j_c`$. This corresponds to the so called “weak pinning” regime.

We argue that many of the varieties of collective transport dynamics observed in superconductors may be generic to repulsive particle systems driven through a disordered media. These behaviors are directly related in our coarse-grained model to large scale morphologies of flow down a vortex density gradient and changes from filamentary strings, to a braided river, to uniform flow at high applied currents. Since the BP model arguably contains the essential physics of the disordered flux line system, at a coarse-grained level, and reproduces a wide variety of experimental results on transport properties, we propose that these self-organized, large scale flow morphologies are also governing the nonlinear, plastic dynamics in the actual physical system: flux lines driven through a superconductor with a disordered pinning landscape by an applied transport current.

### A Summary

In the next section, we focus attention on an extensive study of $`2HNbSe_2`$ as summarized by Higgins and Bhattacharya. The reader may already compare Fig. 1 and Fig. 2, which are schematic reproductions of experimental results, with Figs. 6 and 7, which are results from numerical simulations presented in this work. Section III presents a short review the standard theoretical picture, associated with the work of Shi and Berlinsky, on tearing of the FLL and the resultant nonlinear IV as found via molecular dynamics (MD) simulations. These MD simulations do not, however, take into account a vortex density gradient, and do not agree with some the varieties of plastic dynamics observed in experiments. Section IV describes the original coarse-grained model, first used to describe vortex avalanches and vortex rivers as seen in experiments , and explains the exact model used here to describe the transport experiments summarized in Section II. To eliminate all potentially spurious sources of noise, a completely deterministic variant of the original BP model is used. To describe the IV experiments, we apply a shift in boundary conditions or generalized “tilt” of the vortex pile to represent the effect of a transport current. The resultant vortex flow represents the measured voltage. Section V contains the main numerical results and comparisons with both experiments and MD simulations. The last section summarizes our conclusions.

## II Brief Summary of transport measurements

In a series of papers, Bhattacharya and Higgins describe experiments on the nonlinear transport properties of the FLL in the anisotropic superconductor $`2HNbSe_2`$. In this system, the pinning is extremely weak ($`j_c/j_010^6`$), the lattice is well-formed and a robust “peak effect” occurs slightly below $`H_{c2}`$. Since the London length, $`\lambda `$, is much less than the film thickness, the system operates in the three dimensional regime, where the flux line interactions decay exponentially for lengths larger than $`\lambda `$. They use the strong magnetic field dependence of the critical current (which they interpret in terms of a changing rigidity of the FLL) to explore the crossover between different type of dynamics including “elastic” flow, “plastic” flow, and “fluid” flow. In this regard, the material they study is an ideal experimental system, allowing the exploration of very different regimes in a well-controlled manner.

Near the upper critical magnetic field, they observe a pronounced peak in the pinning force. This is referred to as a “peak effect”. (It should be distinguished from the peak in the differential resistance to be described later.) Equivalently, the critical current, $`I_c`$, where some flux lines start to move, increases as the external magnetic field, $`H`$, increases. For $`I<I_c`$ the entire flux system is pinned, while for $`I>I_c`$, some magnetic flux lines flow across the sample leading to the onset of a finite voltage, $`V`$. We will focus mostly on results associated with the current-voltage (IV) relation as the external magnetic field, $`H`$, is varied.

As shown in Fig. 1b of Ref. , in the “peak regime” the IV characteristics of the superconductor vary enormously. For clarity, these experimental observations are illustrated schematically in Fig. 1. First, on increasing $`H`$, the critical current, $`I_c`$, increases. This is the “peak effect”. Second, below some threshold magnetic field, the IV curves always rise concave upward from $`I_c`$. This is the generic form of the IV for an FLL that is usually reported in the literature, and it is associated with an “elastic” regime. Third, when the external magnetic field enters the “peak regime”, the IV curves change drastically, starting as concave upwards but then bending over after a pronounced inflection point associated with a change of curvature. In the plastic regime, one obtains a characteristic S-shape IV curve. Close to onset, $`I_c`$, it is concave upward, but then bends over as $`I`$ increases further, saturating to a finite slope at large currents. Between the S-shaped and elastic IV curves may be a special curve which is always convex for $`I>I_c`$. This appears at approximately 5.8T in the experiments where the inflection point has moved to onset. Finally, above another magnetic field value, the inflection point is at currents larger than those used in the experiments, and there is no saturation in the slope of the IV curves that can be observed. The numerical simulations presented here reproduce the entire progression of these curves in the plastic regime and their changes as a parameter in the model, representing vortex interactions, is varied. It does not reproduce the behavior in the “elastic” regime, as explained later. The elastic regime is also not observed in two-dimensional MD simulations of driven vortices near the onset of flow.

Experimental measurements of the differential resistance, $`R=dV/dI`$, reveal a peak in $`R`$, corresponding to the inflection point for the S-shaped IV curves. The magnitude of this peak diminishes as the external magnetic field approaches $`H_{c2}`$, and the position of the peak shifts to higher currents. Except for the largest magnetic fields used, the differential resistance eventually saturates at constant values at large enough currents, indicating Ohmic or fluid like behavior for sufficiently high driving. Also, in the plastic regime corresponding to the S-shaped IV curves, the usual scaling ansatz associated with dynamic critical phenomena $`V(II_c)^\beta `$ does not appear to hold. These experimental observations are illustrated schematically in Fig. 2. The numerical simulations presented here reproduce this precise pattern of behavior for the differential resistance in the “peak” regime, including the changing position (to higher currents) and decreasing magnitude of the peak resistance as vortex interactions weaken, and the saturation to Ohmic or fluid-like behavior at high currents.

Over a narrow range of parameters in the plastic regime, Bhattacharya and Higgins observed jaggedness in the differential resistance. This corresponds to secondary peaks in addition to the main peak in the differential resistance. As the external magnetic field is varied, the peaks can be made to appear and disappear, but for a given value of magnetic field and for a given sample, the peaks are reproducible, and thus act as “fingerprints” of the underlying pinning disorder. Our model also exhibits such IV fingerprints.

## III Some Comments on Conventional Theory

The usual interpretation of the peak effect follows the original proposition of Pippard, that softer systems are pinned more strongly than more rigid ones. This picture has been put on a formal basis by Larkin and Ovchinnikov’s theory of collective pinning . It was found that critical current increases as the external magnetic field increases near $`H_{c2}`$ because the FLL elastic moduli soften. However, collective pinning theory doesn’t explain the shape of the IV curve outside of the elastic regime. Some advances have recently been made by Le Doussal and Giamarchi to probe the transition to plastic flow and other features .

MD simulations have indicated that, at least for the two dimensional case, the pinning forces and IV characteristics are determined by plastic deformations (tearing of the FLL) which fall outside the region of validity of collective pinning theory. It is believed that the two dimensional elastic system is always unstable at onset, when a finite voltage first appears, and flow must take place in the plastic regime. These simulations can also be interpreted as describing two dimensional cross-sections of the three dimensional FLL in the plastic regime. The scale of the MD simulations is that of the vortex cores, which is much smaller than the London length, $`\lambda `$.

Shi and Berlinsky presented a limited analytic treatment of the way in which “lattice defects alter the flow characteristics of the lattice under the influence of an external drive”. However both the density of lattice defects and the IV relation were determined by numerical MD simulations. These simulations showed that at high currents where the IV relation is linear, the defect density drops. It increases sharply at lower currents where the IV relation develops an S-shaped curve.

Further simulations have shown that the two dimensional plastic flow takes place in terms of rivers of moving vortices separated by islands where the vortices do not move. This channel flow behavior has been observed experimentally using Lorentz microscopy. Previous simulations of the coarse-grained BP model revealed results consistent with these and also showed that the channels can form a braided river which exhibit self-affine (multifractal) behavior similar to fluvial braided rivers (see for example Ref. ).

Nori and collaborators first studied, using MD simulations, flux driven into superconductors with random pinning, with the driving force solely due to the flux-density gradient. They elucidated many properties of the Bean state including the magnetic field profile, magnetization hysteresis loops, critical currents, vortex avalanches, and vortex rivers . None of these studies using an open system with an overall density gradient reported the IV characteristics, though.

Recently, Nori et al  have simulated the IV curve as the FLL softens by varying the vortex interaction parameter (see also Ref. ). As in all previous MD studies of IV behavior of the FLL , the vortices are contained in a periodic system, where they can neither enter nor leave. The initial condition is an ordered vortex lattice. Motion takes place via MD updates, with a uniform force applied to all vortices.

The most important difference, besides the scale of the model, with our description of the IV experiments, is that, in all of these MD studies of the current-voltage relation, no overall vortex density gradient can develop owing to the periodic boundary conditions. This makes the vortices travel perpetually around the system and forces the paths to circle, which is unphysical, and does not occur in the real system. The actual physical situation is an open system with flux pushed in and out, rather than a periodic one. The physical system adjusts its overall magnetic flux profile in response to applied forces. This is not possible in a periodic, closed system. Nevertheless, the basic result that the critical current increases as the vortex interactions weakens is obtained. However, all the IV curves measured via MD simulations of periodic systems, fall on top of each other, or overlap, at high currents, for different values of the vortex interaction strength. This is inconsistent with the experiments of Bhattacharya and Higgins (see Fig. 1), and reflects the fact that the artificially periodic system is not able to adjust its profile in response to applied forces. An even more significant, but related, difference is that the peak in the differential resistance grows monotonically and gets sharper as the resistance peak shifts to higher currents, the exact opposite of what happens in experiments.

In short, the systematically changing morphology of the IV curves and differential resistance for a soft FLL has not been theoretically described.

## IV The Model

### A Motivation

We take a different approach to modeling flux flow in superconductors. In order to describe collective transport in disordered media, it is reasonable to consider coarse-grained models where the details of the precise interactions between flux lines (or vortices) are lost but the general effects of repulsive interaction between granular or discrete objects, pinning, and over-damped motion leading to stick-slip dynamics (tearing) are preserved.

The BP model is an interacting sand-pile model of vortices in a type-II superconductor, where the “sand” grains, representing magnetic vortices, repel each other. It was originally motivated by the observation of (possibly) self-organized critical avalanches in field ramping experiments, where the distribution of flux packages falling into the interior coil of a hollow cylindrical superconductor were measured. Actually, the similarity between the Bean state (or vortex pile) and sandpiles was first pointed out by de Gennes. Later, Vinokur, Feigel’man, and Geshkenbein suggested that thermally induced flux creep would lead to a self-organized critical state in a type II superconductor, as did Tang. The key observation is that flux line flow always takes place on a flux pile which has an overall density gradient. This pile may be in a self-organized critical, or some other nontrivial state.

In addition to describing vortex avalanches in field ramping experiments, the BP vortex model has also been used to describe flux noise, vortex avalanches in the presence of a periodic, dense array of pinning centers, thermally activated flux creep, and magnetization loops.

#### 1 Coarse-graining

Consider a transverse two-dimensional slice of a superconducting slab at $`T=0`$. The BP model results from a coarse-grained description to the scale of the London length, $`\lambda `$, of the microscopic dynamics, and incorporates the features that are essential to produce the observed complex behavior: Repulsive interactions between vortices, variations in the pinning potential, and variations in the vortex density - all at the scale of $`\lambda `$. One can imagine imposing a grid of cells on the system; vortices in the model correspond to a vortex number in an extended region of the actual physical system. Pinning in the model corresponds to a number of point pins in an extended cell. Each lattice site in the model can hold many vortices, and can have a different, albeit quenched, pinning potential, due to the underlying randomness in the positions and strengths of the microscopic pinning centers. The model allows many vortices to interact with each other while maintain locality (at the scale of $`\lambda `$) in interactions. It is the only model of this sort that has been proposed to describe flux dynamics in superconductors. It enables numerically studies, using ordinary workstations, of the steady state and transient properties of systems larger than $`(500\lambda )^2`$ containing tens of millions of vortices.

### B Definition

The BP model is defined as follows (see Fig. 3). Consider a two-dimensional honeycomb lattice  where each cell, $`x`$, has three nearest neighbors, and is occupied by an integer number of vortices, $`m(x)`$. The total energy of the vortex system includes the repulsive pairwise interaction between the vortices and the attractive interaction of vortices with the pinning potential, $`\widehat{V}`$. For a given configuration of vortex number $`\{m(x)\}`$, the total energy of the system is:

$$H\left(\{m(x)\}\right)=\underset{i,j}{}J_{ij}m(i)m(j)\underset{i}{}\widehat{V}_{\text{pin}}(i)m(i).$$
(1)
Since the model describes a system coarse-grained to the scale of the London length, the repulsive interactions $`J_{ij}`$ are short-ranged. This usually includes an on-site interaction, and a weaker nearest-neighbor interaction.

As in the microscopic case described by MD simulations, the change in the total energy of the model when moving a unit vortex from one site to a nearest neighbor site is determined. This yields the force to move a unit vortex from $`x`$ to $`y`$, which is

$`F_{xy}`$ $`=`$ $`V_{\text{pin}}(y)V_{\text{pin}}(x)+[m(x)m(y)1]`$ (3)
$`+r[m(x1)+m(x2)m(y1)m(y2)].`$
As indicated in Fig. 3, the nearest neighbor cells of $`x`$ are $`y`$, $`x1`$, and $`x2`$, and the nearest neighbors cells of $`y`$ are $`x`$, $`y1`$, and $`y2`$, and $`0r<1`$. A slightly different implementation of the disorder is used than before. The normalized pinning potential $`V_{\text{pin}}(x)`$ is a random number taken from a uniform distribution in the interval between zero and $`V_{max}`$. In each time step, all cells of the lattice are updated in parallel. A single vortex moves from a cell to a neighboring cell if the force in that direction is positive, or equivalently if the total energy of the system is lowered. In Eq. 2, the units of force on a vortex have been normalized so that the on-site term is unity. Thus there are two dimensionless parameters remaining, $`V_{max}`$ and $`r`$.

Many alternatives exist to handle the situation when more than one unstable direction appears for a vortex to move. In the previous implementation of the model, one unstable direction was chosen at random. In order to simplify the model and eliminate all potentially spurious sources of noise, the most unstable direction that has the largest force is chosen and the vortex moves to that site. This represents an extremal process. In fact the entire model is now completely deterministic, corresponding to a $`T=0`$ limit of the dynamics.

### C The External Magnetic Field: Building a Vortex Pile

Flux lines enter the superconductor from the edges, pushed in by the external magnetic field. This is represented by putting all sites on the left edge of the model in contact with a reservoir of vortices at some potential, corresponding to the external magnetic field on the left side of the sample. The same is done for the right side of the sample. For simplicity, periodic boundary conditions are used for the top and bottom. If the two reservoirs are set equal, representing embedding the sample in an external magnetic field, vortices enter the system generating the classic V-shaped flux density curve as the external magnetic field is increased (see below).

#### 1 Details about Boundary Algorithm

In order to calculate the force on a vortex to move to or from a boundary site, a special algorithm must be used because one of the nearest neighbor sites of each boundary site is not on the lattice. The rule used here simply assumes that the “virtual” off lattice site neighboring each boundary site is occupied with an equal number of vortices as that boundary site. More specifically, at the beginning of each lattice update, all of the sites on the boundary are set to be occupied with the same number of vortices. The lattice update then proceeds, during which vortices can move off the boundary sites into the system, or from the system onto the boundary sites, thereby changing the number of vortices occupying a boundary site. However, at the beginning of the next lattice update all of the sites on the boundary are reset to their original value. Through this process, vortices can be removed or added to the system. In general, the left and right boundaries are held at different values. There is no pinning at the boundary sites.

#### 2 Details about Parallel Update

An artifact of parallel updating is the existence of local instabilities in which two, or more, vortices oscillate back and forth between neighboring sites. These local instabilities disappear if the model is coarse-grained, because then the neighboring sites are incorporated into a single one, and therefore are not important to the large scale behavior of the system. The instabilities can be eliminated by keeping track of the direction from which the last vortex moved onto each site and always forbidding a return movement backwards in that direction. Otherwise, backward jumps are rare and therefore disallowing them does not change the large scale behavior of the system. A similar rule applies to the boundary sites. The advantage of the parallel update is that it is numerically more efficient, than other update schemes and does not introduce any uncontrollable spurious effects.

### D The Transport Current: Shifting the Boundary Conditions

We consider an infinite slab of finite thickness in a parallel applied magnetic field, carrying a current perpendicular to the field. Depending on the direction of the applied current, the magnetic field on one side of the superconductor (e.g. the right hand side) will be decreased and the magnetic field on the other side will be increased (e.g. the left hand side). Assuming the applied magnetic field is sufficiently strong, this corresponds just to changing the heights of the magnetic flux pile on the left and right edge, or a general kind of “tilt” of the pile. See, for example, Ref. . Of course, the internal currents and forces inside the superconductor will readjust to accommodate this new boundary condition. The flux lines are considered to be perfectly stiff, and described by a two-dimensional slice of the three dimensional slab.

An applied magnetic field in the geometry described above will produce the well known “V profile” of the magnetic flux density discussed, for example, in Orlando and observed in many experiments, such as those of Behnia, et al. This is shown in Fig. 4. The actual magnetic profile depends on the history of the sample and how magnetic field has been applied in the past. Applying a finite current shifts the boundary conditions resulting in a hysteretic profile, also shown in Fig. 4.

#### 1 Onset

Eventually, as the applied current is increased further, a critical current is reached, where the shift of the boundary conditions is so large that steady vortex flow occurs down the gradient spannning the entire sample. A profile just above the critical current is also shown in Fig. 4. Note that this profile also retains some hysteretic properties. For example, it has a bump corresponding to where the increasing and decreasing portions of the magnetic field profile merged as the external current was increased above the threshold. This bump disappears if the boundary conditions on the sample are shifted further and then lowered back to the previous value.

### E Making IV Measurements

The IV characteristic is determined by the relation between applied current and the vortex flow, which induces a voltage. To our knowledge, this type of numerical measurement, made by shifting the boundary conditions on the Bean state, has not been investigated before. Here the IV characteristic is simply the relation between the magnitude of the shift (representing an applied current) and the average flow rate of vortices (representing a voltage) when the critical current (tilt) is exceeded. In order not to confuse the reader we use the term current to refer to the applied electrical transport current and the term flow to refer to the motion of magnetic vortices.

In general, the boundary sites in the model can be set to any real value, including non-integer values. This is because vortex number on the boundary sites describes the external magnetic field density at the boundary of the sample. However, only integer units of magnetic flux can enter the interior of the system, and thus only integer numbers of vortices occupy interior lattice sites. Obviously, the magnitude of the difference between the heights of the left and right boundaries, can also be set to non-integer values. However, if the boundary heights are shifted by less than whole integers, steps can appear in the IV data. These steps are caused by the fact that only whole numbers of vortices can occupy interior lattice sites. Since one vortex unit in the model represents many actual physical vortices, this is to some extent an unphysical artifact of the model, and the discreteness effect should be less apparent in experiments. To eliminate this effect, all of the IV data presented in this paper were calculated by shifting the boundary heights by only integer numbers of vortices.

In simulations of the model, there is almost no backward movement of vortices. Thus, the vortex flow can be determined by measuring the average number of vortices moving per lattice update (the average activity). Furthermore, the velocity of each moving vortex is one lattice site per update. Thus, the experimentally measured voltage, which is equal to the amount of moving flux times the velocity of that flux, is therefore proportional to the average vortex activity.

IV data from simulations of different size systems collapse nicely in a scaling plot, Fig. 5, if the voltage is measured as vortices moving per lattice update per lattice site (the average activity per site), and the current is measured as the average slope across the system (the magnitude of the height difference between the left and right boundaries divided by the length of the system). In the following sections, IV data is presented in these scaling units. In Fig. 5, the IV data was produced by repeatedly increasing the height of the left boundary by 1 vortex and lowering the height right boundary by 1 vortex. The vortex interaction strength was $`r=0.1`$.

The IV data presented in the following sections was calculated in a similar fashion. All of the data is for systems of size 200x400. Each IV data point was calculated by first shifting the left boundary height up 1 vortex and shifting the right boundary height down 1 vortex. As for Fig. 5, the lattice was updated 20000 times to eliminate transient behavior, and finally the lattice was updated another 20000 times during which the average vortex activity was measured.

## V Collective Transport Behavior of the Model

The IV relation was measured for different values of the parameters and for different system sizes. Our main result is shown in Fig. 6, where the parameter $`V_{max}=5`$ and the parameter $`r`$ is varied. The parameter $`r`$ represents the strength of repulsion between vortices at nearest neighbor cells (of size $`\lambda `$).

The first result is that as $`r`$ decreases, the critical current, which is the slope of the pile where vortices first start to flow, increases. Clearly, applying an increased $`r`$ to the pile in the steady state lowers the slope of the pile since formerly stable local slopes will now become unstable due to the increased repulsion between vortices at neighboring sites. Thus, fixing all other parameters, we can identify the parameter $`r`$ in the model as a way of controlling the critical current, or slope of the pile. In the superconductor, the critical current can be controlled by the applied magnetic field. In the peak effect regime, it turns out that increasing the applied magnetic field leads to a softer FLL and thus a higher critical current. Therefore, the regime where the critical current is an increasing function of the applied magnetic field is represented in our model by a critical current which is a decreasing function of the parameter $`r`$. This is made evident by comparing Figs. 1,2 and Figs. 6,7.

The second result is that, except for the largest $`r`$, all of the IV curves have a characteristic S shape. They start out at some $`I_c`$ increasing concave upward and then bend over saturating to a finite slope at large currents. The IV relations for a given realization of disorder, however, do not overlap at high currents, unlike the results obtained with previous MD simulations . Again, all of this agrees with the experimental results.

In Fig. 7, the differential resistance $`dV/dI`$ is calculated by taking numerical derivatives of the curves shown in Fig. 6. At the largest $`r`$ value, the resistance peak is very large. As the parameter $`r`$ decreases, the peak moves to higher currents and decreases in magnitude. In fact the same behavior is observed in experiments as the external magnetic field approaches $`H_{c2}`$. As mentioned before, the MD simulations give the opposite result of a peak that increases in magnitude as it shifts to higher currents.

In the actual experiments, at fields below the peak effect regime, the vortex flow is believed to be elastic and there is no observed peak in the differential resistance. As the magnetic field increases a peak in the differential resistance starts to develop, which then reaches a maximum, decreasing again for larger fields. Our simulations appear to describe the experiments once the peak in differential resistance has reached its maximum.

In order to obtain an elastic regime, we could consider a three dimensional coarse-grained model of repelling flux lines rather than point vortices representing a two-dimensional cross section of that system. Work is in progress along those lines.

The model discussed here does not describe the behavior of the superconductor in fields greater than that which gives the largest critical current, where the critical current decreases as the magnetic field increases. This may be due to the fact that we take the depairing current to be strictly infinite (i.e. $`j_c/j_0=0`$) and there is no transition to a non-superconducting state in the model presented here.

Note that the differential resistance curves for small $`r`$ contain secondary peaks, in addition to the main peak. This is similar to the jaggedness or “fingerprint” found in experiments. This jaggedness in our results occurs in the filamentary channel regime, discussed below, and is due to filaments opening and closing as the applied current increases.

### A Relation to River Morphology

So far we have only characterized the model using measurements analogous to those that experimentalists typically have available. However, numerical simulations can also provide a great deal of easily accessible information about the morphologies of flow patterns associated with different configurations. In fact, different morphologies of vortex flow patterns have been observed using MD simulations . To understand the nature of the differences in the shape of the IV curves, we have examined how the flow morphologies change as one increases the applied current for a system that has an S-shaped IV curve, and also for the singular IV curve that occurs at large $`r`$.

The paths that the vortices take as they cross the sample can be determined by measuring the average activity at each site. For example, Fig. 8 shows a series of grey scale images of the vortex flow patterns. These images represent “time-lapsed photographs” of the vortex activity. The different images in Fig. 8 show the vortex flow patterns for the case $`r=0.1`$ as a function of increasing external current, $`I`$. In these images, sites with no vortex activity are blank, sites with an activity greater than or equal to 0.5 (one vortex moves every other lattice update) are black, and sites with activity between 0 and 0.5 are indicated by grey dots with a darkness proportional to their activity. These patterns are fixed and do not change in time with fixed external driving conditions.

The nature of the flow can be quantified by the distribution of activity at the different lattice sites. Figure 9 shows histograms of the activity corresponding to the images of Fig. 8. These histograms are constructed with 1000 bins and normalized so that the area under the curve is equal to one. Note that there are peaks corresponding to zero activity not included in the figure. Those peaks at zero activity decrease in size as the current is increased.

As can be seen in Figs. 6 and 7, the IV curve for $`r=0.1`$ is S-shaped. Figure 8a shows the vortex flow pattern just above threshold for that case. The vortex flow takes place only on a single filamentary string with some small side branches. This behavior is also evident in the histogram Fig. 9a, which shows a small number of isolated peaks. As the tilt of the pile increases the number of filaments of the vortex flow increases, and they begin to merge. This process can be seen in Figs. 8b and 8c. The corresponding histograms of activity shown in Figs. 9b and 9c indicate an increasing number of peaks, and the development of a continuous distribution. During this merging process, while the vortex flow remains filamentary, the IV curve remains concave upward, and the differential resistance rises.

Eventually the filaments of vortex flow merge to form a braided river, as shown in Fig. 8d. This occurs around the peak in differential resistance, corresponding to a change in curvature of the IV relation. In this case, the corresponding histogram of activity (Fig. 9d) shows a continuous distribution of activity peaked at zero. Thus the peak in differential resistance signals a change in the underlying vortex flow morphology from filamentary strings to a braided river. This appears to be consistent with the observation of Kolton et al that at the peak in differential resistance “all the vortices are moving in a seemingly isotropic channel network with maximum interconnectivity” .

As the tilt of the pile is increased even further, the vortex flow becomes spatially more and more uniform, as the braided river floods. This uniform flow region corresponds to linear, or Ohmic behavior. These results can be seen in the river flow pictures in Figs. 8e and 8f. The corresponding histograms shows a peaked continuous distribution, which has separated from zero activity. As the tilt increases further the continuous distribution of activity obtains an increasing mean and narrowing width. The transition to Ohmic behavior comes about as the braided river floods at higher flow rates than can be supported by such a structure.

\[

\]

\[

\]

A similar progression of flow morphologies, from isolated filamentary channels, to a braided river, to a flooded river with uniform flow, occurs for other values of $`r`$ with S-shaped IV curves. However, a different scenerio occurs at large values of $`r`$ where the IV curve is not S-shaped, but is instead always concave downward from the onset of vortex motion. The IV curve for $`r=0.4`$ in Fig. 6 is an example. We will refer to this as the critical curve; its functional form is discussed in the next subsection. The differential resistance for the critical curve is discontinuous, as can be seen in Fig. 7.

For the critical curve, there is no filamentary region of vortex flow. Instead, right at the onset of vortex motion, the flow has the structure of a braided river. This can be seen in Fig. 10a, and in the corresponding activity histogram in Fig. 11a which shows a continuous distribution of activity. (Figures 10 and 11 were produced in the same manner as Figs. 8 and 9, respectively.) As the tilt of the pile is increased even further, the flow again increases and becomes more spatially uniform, similar to behavior of the S-shaped IV case after the peak in the differential resistance, where the braided river floods. This behavior can be seen in the river flow pictures in Figs. 10b, 10c and 10d, and the corresponding histograms in Figs. 11b, 11c and 11d. Since the critical curve has a diverging differential resistance at onset, this further supports our observation that the peak in the differential resistance is associated with a braided river structure.

\[

\]

\[

\]

### B Scaling of the Critical IV Curve

The critical IV curve for $`r=0.4`$ is well described by

$$V(II_c)^\beta ,$$

with $`\beta =0.6\pm 0.1`$, as shown in Fig. 12. For that measurement, the critical current, $`I_c,`$ was measured to be $`0.095`$, which was the largest current found with zero voltage. Since the current was sampled at values spaced by $`0.005`$, the uncertainty in $`I_c`$ is approximately $`\pm 0.005`$. Varying the value of $`I_c`$ over that range changes the value of $`\beta `$ that best fits the data near the onset of vortex motion, and allows an estimate of the error on $`\beta `$.

The exponent $`\beta `$ can be related via scaling arguments to the exponents characterizing the distribution of avalanches in the self-organized critical state. In the self-organized critical state, the average flow rate, $`V`$, of vortices is controlled rather than the overall slope of the system. The scaling argument is similar to that used in Corral and Paczuski  to describe the transition from avalanche to continuous flow in the one dimensional Oslo model . The excess slope above the critical slope at onset is $`\mathrm{\Delta }m=II_c=(\mathrm{\Delta }N)L^3`$, where $`\mathrm{\Delta }N`$ is the excess number of vortices in the system, and $`L`$ is the system size. If vortices are added very slowly then there will be distinct avalanches separated by intervals of no activity. In that regime, superposition applies and $`\mathrm{\Delta }NVL^z`$, where $`L^z`$ is the cutoff in the duration of the avalanches. In the rapidly driven regime, the avalanches overlap, and the excess slope becomes independent of system size, thus $`\mathrm{\Delta }mV^{1/\beta }`$. These two limits can be combined into a single cross-over scaling function:

$$\mathrm{\Delta }NL^af(VL^x).$$
(4)
The exponent $`x`$ in this expression measures the average duration of avalanches $`tL^x`$. Since avalanches in this case come about from adding an entire row of $`L`$ vortices to the system, the average size of avalanches, measured in terms of the number of topplings is $`sL^2`$, rather than $`sL`$, as in the Oslo model. Obviously, on average each site in the system topples once when a row of vortices is added in the self-organized critical state. Using this result together with conservation of probability gives a scaling relation between $`x`$ and $`z`$;

$$x=z+2D.$$

Combining all these results gives

$$\beta =\frac{x}{x+2z}=\frac{z+2D}{4D}.$$

Using the exponent values $`z=1.5`$ and $`D=2.7`$ obtained in Ref.  gives $`\beta =0.6`$, in good agreement with the numerical result presented in Fig. 12.

## VI Summary

We have studied the nonlinear current-voltage characteristics of flux flow in type II supeconductors with random pinning using a simple, cellular model. As in the physical system, vortices flow down a flux density gradient. Because the coarse-grained model does not include any degrees of freedom at scales much smaller than the London length, $`\lambda `$, it can not describe structural ordering or disordering of the flux line lattice. Despite this fact, simulations of the model reproduce many of the empirically observed transport characteristics of plastic flux flow in type II superconductors.

In particular, our results reproduce many of the features attributable to plastic flow of the FLL in the “peak regime”. By weakening the vortex interaction strength, we find an increase of the critical current and a falling of the magnitude of the peak in the differential resistance. The model also exhibits IV fingerprints, and crossover to Ohmic or linear behavior at high currents. Also, the IV curves for different vortex interaction strengths do not merge at large external currents. All of these features are completely consistent with experimental results. However, presumably due to the two-dimensionality of the model, the elastic behavior of the FLL is not reproduced.

The success of these efforts to describe the plastic transport behavior of magnetic flux in superconductors with a coarse-grained cellular model suggests a possibly generic explanation for and ubiquity of plastic flow phenomena observed in superconductors. It may not depend on the degree of disorder or defects in an underlying microscopic flux line lattice. Instead the behavior may be common to driven repulsive particle systems in a disordered media. The varieties of plastic flow behaviors in the model studied here result from the changing morphologies of the vortex flow pattern down a vortex density gradient, self-organizing into distinct large scale patterns. These include isolated filaments, which merge at higher flow rates to give a braided river, and lead to uniform flow, or Ohmic behavior, at the highest flow rates. The filamentary structure is associated with a concave IV characteristic; the braided river structure is associated with the peak in the differential resistance; the change to Ohmic behavior comes about as the braided river floods as it can not support a higher level of flow.

###### Acknowledgements.

This work was supported by the NSF, grant #DMR-0074613, and also by the Texas Center for Superconductivity. E. Altshuler thanks the World Laboratory for Pan-American Collaboration in Science and Technology for financial support. We thank Per Bak for comments on the manuscript, and C.S. Ting and G.F. Reiter for informative conversations on superconductivity.


text:
# From Superluminal Velocity To Time Machines?

Various experiments have shown superluminal group and signal velocities recently. Experiments were essentials carried out with microwave tunnelling , with frustrated total internal reflection , and with gain-assisted anomalous dispersion . According to text books a superluminal signal velocity violates Einstein causality implying that cause and effect can be changed and time machines known from science fiction could be constructed. This naive analysis, however, assumes a signal to be a point in the time dimension neglecting its finite duration. A signal is not presented by a point nor by its front, but by its total length. On the other hand a signal energy is finite thus its frequency band is limited, the latter is a fundamental physical property in consequence of field quantization with quantum $`h\nu `$. All superluminal experiments have been carried out with rather narrow frequency bands. The narrow band width is a condition sine qua non to avoid pulse reshaping of the signal due to the dispersion relation of the tunnelling barrier or of the excited gas, respectively . In consequence of the narrow frequency band width the time duration of the signal is long so that causality is preserved. However, superluminal signal velocity shortens the otherwise luminal time span between cause and effect.

Can a signal travel faster than light? If this happens, would it really violate the principle of causality stating that cause precedes effect ? The latter statement has been widely assumed as a matter of fact. It has been shown according to the theory of special relativity that a signal velocity faster than light allows to change the past. The line of arguments how to manipulate the past in this case is illustrated in Fig. 1 . There are two frames of reference displayed. In the first one at the time t = 0 lottery numbers are presented, whereas at t = -10 ps the counters were closed. Mary (A) sends the lottery numbers to her friend Susan (B) with a signal velocity twice the velocity of light. Susan moving in the second inertial system at a relative speed of 0.75c, sends the data back at an even faster speed of 4c, which arrives in the first system at t = -50 ps, thus in time to deliver the correct lottery numbers before the counters close at t = -10 ps.

The time shift of a point on the time coordinate into the past is given by the relation :

$$t_A=\frac{L}{c}\frac{(cc/Nc/N^{}+v_r/NN^{})}{(cv_r/N)},$$
(1)
where $`L`$ is the transmission length of the signal, $`v_r`$ is the relative velocity of the two inertial systems A and B, and Nc, N’c are the signal velocities in A and B, respectively. N and N’ are numbers assumed to be $`>`$ 1.

This is an example often encountered in the literature supposed to show that a superluminal signal velocity results in negative times and allows to manipulate the past. We show now that this simple model is not correct.

First we are going to recall the basic properties of a signal. Microwave pulses  and quite recently light pulses  of frequency $`\nu `$ and bandwidth $`\mathrm{\Delta }\nu `$ have been shown to travel at a velocity much faster than light. The pulses in the two experiments correspond to signals used nowadays in telephone as well as in inter-computer communication. Frequency band limitation of signals, the basis of the sampling theorem , is a backbone of digital communication technology discussed in detail in the literature (e.g. in Encyclopedia Britannica), but scarcely addressed to in textbooks of physics. The finite energy content of a signal actually implies the frequency band limitation . This fundamental physical property is in consequence of the energy of any frequency components of a signal to be $`nh\nu `$ where $`n`$ is a whole number, $`h`$ the Planck constant, and $`\nu `$ the frequency.

A pulse represents an amplitude modulated (AM) signal on a carrier frequency. The carrier frequency is in charge of the receivers address and the half-width of the pulse represents the number of digits, i.e. the information. In the case of modern fiber optics the relative frequency band width is $`10^3`$ approximately, in the superluminal microwave experiments the band width was $`10^1`$ and in the optical experiment mentioned above it was less than $`10^9`$. Due to the narrow frequency bands there was no significant pulse reshaping neither in the microwave tunnelling experiment nor in the gain-assisted light propagation experiment. The superluminal signals are shown together with the luminal reference signals in Fig. 2.

Thus in both experiments the signal travelled at a superluminal velocity, e.g. with 4.7c  or with –310c , respectively. Nevertheless, the principle of causality has not been violated in both experiments.

In the example with the lottery data the signal was assumed to be a point on the time coordinate. However, a signal has a finite duration as the pulse sketched along the time coordinate in Fig. 3.

(In the experiments in question 7.5$`\mu `$s and 5 ns, see Fig. 2.) Any information like a word has a finite extension on the time coordinate. In the two cited superluminal experiments the superluminal time shift compared with the pulse length is about 30$`\%`$ in the microwave experiment with the velocity 4.7c and about 1 $`\%`$ in the light experiment with the velocity –310c. Due to the signal’s finite duration of 200 ps the information is obtained only at positive times under the assumptions as illustrated in Fig. 3. The same holds a fortiori for the two discussed experiments. The finite duration of a signal is the reason that a superluminal velocity does not violate the principle of causality. On the other hand a shorter signal corresponds to a broader frequency band. In consequence of the dispersion relation of either a tunnelling barrier or of an excited atomic gas with an extremely narrow frequency regime of anomalous dispersion strong pulse reshaping would occur. Summing up, the principle of causality has not been violated by the experiments with superluminal signal velocities, but amazing the time span between cause and effect has been reduced compared with luminal propagation.

Acknowledgment
We gratefully acknowledge helpful discussions with P. Mittelstaedt.


text:
# Nuclear fragmentation: sampling the instabilities of binary systems

## Abstract

We derive stability conditions of Asymmetric Nuclear Matter ($`ANM`$) and discuss the relation to mechanical and chemical instabilities of general two-component systems. We show that the chemical instability may appear as an instability of the system against isoscalar-like rather than isovector-like fluctuations if the interaction between the two constituent species has an attractive character as in the case of $`ANM`$. This leads to a new kind of liquid-gas phase transition, of interest for fragmentation experiments with radioactive beams.

PACS numbers: 21.65.+f, 25.70.Pq, 71.10.Ay

A binary system manifests a richer thermodynamical behaviour as a consequence of a new parameter, the concentration, which is required for a complete description of its states. The phase transitions are more complex because they have to accommodate one more conservation law. A thermodynamical state can be unstable not only mechanically but also chemically.

Nuclear matter belongs to this class of systems. The process of multifragmentation following the collision of heavy nuclei at medium energies is expected to show features analogous to usual liquid-gas phase transitions of water . In particular, the spinodal decomposition appears to be an important mechanism leading to phase separation in symmetric nuclear matter ($`SNM`$) , as confirmed in recent experiments . The relevance of instabilities in $`ANM`$-fragmentation has been already discussed few years ago . It was shown that in these systems a kind of diffusive (or chemical) spinodal rather than the mechanical spinodal is significant. However the detailed nature of fluctuations responsible for such instability has not been clarified yet. This is a quite important aspect since, as we will see, it provides reliable information on the character of the interaction in the medium. We remark that the effect is leading to the ”isospin distillation” that can be experimentally observed, see the recent data .

In one-component systems the mechanical instability is related to instability against density fluctuations as the result of the strong attraction between constituents. In symmetric binary systems, like $`SNM`$, one encounters two kinds of density fluctuations: i) isoscalar, when the densities of the two components oscillate in phase with equal amplitude, ii) isovector when the two densities fluctuate still with equal amplitude but out of phase. Mechanical instability is associated with instability against isoscalar fluctuations leading to cluster formation while chemical instability is related to instability against isovector fluctuations, leading to species separation. We will show that in asymmetric binary systems, as in $`ANM`$, this direct correspondence between the nature of fluctuations and the thermodynamical instability is lost.

An appropriate framework for this investigation is provided by the Fermi liquid theory which has been already applied to the study of instabilities in symmetric binary systems as $`SNM`$ (the two components being protons and neutrons) and the liquid $`{}_{}{}^{3}He`$ (spin-up and spin-down components) .

We first investigate the thermodynamical stability of $`ANM`$ at $`T=0`$ extending to the asymmetric case the formalism introduced in . The distribution functions for protons and neutrons are:

$$f_q^{(0)}(ϵ_p^q)=\mathrm{\Theta }(\mu _qϵ_p^q),q=n,p$$
(1)
where $`\mu _q`$ are the corresponding chemical potentials. The nucleon interaction is characterized by the Landau parameters:

$$F^{q_1q_2}=N_{q_1}V^2\frac{\delta ^2}{\delta f_{q_1}\delta f_{q_2}},N_q(T)=\frac{2d𝐩}{(2\pi \mathrm{})^3}\frac{f_q(T)}{ϵ_p^q}$$
(2)
where $``$ is the energy density, V is the volume and $`N_q`$ is the single-particle level density at the Fermi energy. At $`T=0`$ this reduces to $`N_q(0)=mp_{F,q}/(\pi ^2\mathrm{}^3)`$, were $`p_{F,q}`$ is the Fermi momentum of the $`q`$-component. Thermodynamical stability for $`T=0`$ requires the ground state energy to be an absolute minimum for the undistorted distribution functions, such that the relation:

$$\delta \mu _p\delta \rho _p\mu _n\delta \rho _n>0$$
(3)
is satisfied when we deform proton and neutron Fermi seas. The distorted distribution functions can be written as $`f_q(𝐩)=\mathrm{\Theta }(ϵ_{F,q}(\theta )ϵ_p^q)`$, where $`ϵ_{F,q}(\theta )`$ is a direction dependent Fermi energy characterizing the distortion. We will follow the usual multipole expansion for the variation

$$[\delta ϵ_{F,q}(\theta )]=ϵ_{F,q}(\theta )\mu _q=\nu _q^lP_l(cos\theta )$$
(4)
and for the Landau parameters

$$F^{q_1q_2}(\theta _1\theta _2)=\underset{l}{}F_l^{q_1q_2}P_l(cos(\theta _1\theta _2)).$$
(5)
We limit ourselves to monopolar deformations considering here momentum independent interactions such that $`F_{l=0}^{q_1q_2}`$ are the only non-zero Landau parameters . Then, up to second order in the variations, Eq.(3) becomes

$$\delta \mu _p\delta \rho _p\mu _n\delta \rho _n=\frac{1}{2}(a\nu _{p}^{}{}_{}{}^{2}+b\nu _{n}^{}{}_{}{}^{2}+c\nu _p\nu _n)$$
(6)
where $`\nu _{n,p}\nu _{n,p}^0`$ and

$`a=N_p(0)(1+F_0^{pp});b=N_n(0)(1+F_0^{nn});`$ (7)
$`c=N_n(0)F_0^{pn}+N_p(0)F_0^{np}=2N_n(0)F_0^{pn}.`$ (8)
We diagonalize the r.h.s. of Eq.(6) by introducing the following transformation:

$`u`$ $`=`$ $`cos\beta \nu _p+sin\beta \nu _n,`$ (9)
$`v`$ $`=`$ $`sin\beta \nu _p+cos\beta \nu _n,`$ (10)
where the mixing angle $`0\beta \pi /2`$ is defined by

$$tg2\beta =\frac{c}{ab}=\frac{N_n(0)F_0^{pn}+N_p(0)F_0^{np}}{N_p(0)(1+F_0^{pp})N_n(0)(1+F_0^{nn})}.$$
(11)
Then Eq.(6) takes the form

$$\delta \mu _p\delta \rho _p\mu _n\delta \rho _n=Xu^2+Yv^2$$
(12)
where

$`X`$ $`=`$ $`{\displaystyle \frac{1}{2}}(a+b+sign(c)\sqrt{(ab)^2+c^2})`$ (13)
$``$ $`{\displaystyle \frac{N_p(0)+N_n(0)}{2}}(1+F_{0g}^s)`$ (14)
and

$`Y`$ $`=`$ $`{\displaystyle \frac{1}{2}}(a+bsign(c)\sqrt{(ab)^2+c^2})`$ (15)
$``$ $`{\displaystyle \frac{N_p(0)+N_n(0)}{2}}(1+F_{0g}^a).`$ (16)
Thanks to the rotation Eq.(10) we separate the total variation Eq.(3) into two independent contributions, the ”normal” modes, characterized by the ”mixing angle” $`\beta `$, which depends on the density of states and the details of the interaction.

In the symmetric case, $`N_p=N_nN`$, $`F_0^{nn}=F_0^{pp}`$ and $`F_0^{np}=F_0^{pn}`$, Eq.(6) reduces to

$`\delta \mu _p\delta \rho _p\mu _n\delta \rho _n`$ $`=`$ $`{\displaystyle \frac{N(0)}{2}}(1+F_0^s)(\nu _p+\nu _n)^2`$ (17)
$`+`$ $`{\displaystyle \frac{N(0)}{2}}(1+F_0^a)(\nu _p\nu _n)^2`$ (18)
where $`F_0^sF_0^{nn}+F_0^{np}`$ and $`F_0^aF_0^{nn}F_0^{np}`$ are symmetric and antisymmetric (or isoscalar and isovector) Landau parameters, and we recover the usual Pomeranchuk stability conditions for pure isoscalar/isovector fluctuations .

In the general case we interpret $`u`$\- and $`v`$-variations as new independent $`\mathrm{𝑖𝑠𝑜𝑠𝑐𝑎𝑙𝑎𝑟}`$-like and $`\mathrm{𝑖𝑠𝑜𝑣𝑒𝑐𝑡𝑜𝑟}`$-like fluctuations appropriate for asymmetric systems. The proton and neutron densities will fluctuate in phase for isoscalar-like variations and out of phase for isovector-like variations, see Eq.(10). Moreover, $`F_{0g}^s`$ and $`F_{0g}^a`$, defined by Eq. (14,16), can be considered as generalized symmetric and antisymmetric Landau parameters.

From Eq.(12) we see that thermodynamical stability requires $`X>0`$ and $`Y>0`$. Equivalently, the following conditions have to be fulfilled:

$$1+F_{0g}^s>0and1+F_{0g}^a>0,$$
(19)
They represent Pomeranchuk stability conditions extended to asymmetric binary systems.

The new stability conditions, Eq.(19), are equivalent to mechanical and chemical stability of a thermodynamical state, , i.e.

$$\left(\frac{P}{\rho }\right)_{T,y}>0and\left(\frac{\mu _p}{y}\right)_{T,P}>0$$
(20)
where $`P`$ is the pressure and $`y`$ the proton fraction, as it can be proved by observing that :

$`XY`$ $`=`$ $`N_p(0)N_n(0)[(1+F_0^{nn})(1+F_0^{pp})F_0^{np}F_0^{pn}]`$ (21)
(22)
$`=`$ $`{\displaystyle \frac{(N_p(0)N_n(0))^2}{(1y)\rho ^2}}\left({\displaystyle \frac{P}{\rho }}\right)_{T,y}\left({\displaystyle \frac{\mu _p}{y}}\right)_{T,P}`$ (23)
and:

$`\left({\displaystyle \frac{P}{\rho }}\right)_{T,y}={\displaystyle \frac{\rho y(1y)}{N_p(0)N_n(0)}}(ta+{\displaystyle \frac{1}{t}}b+c)`$ (24)
$`X(\sqrt{t}cos\beta +{\displaystyle \frac{1}{\sqrt{t}}}sin\beta )^2+Y(\sqrt{t}sin\beta {\displaystyle \frac{1}{\sqrt{t}}}cos\beta )^2`$ (25)
$`witht={\displaystyle \frac{y}{1y}}{\displaystyle \frac{N_n(0)}{N_p(0)}}.`$ (26)
From Eq.s(12,19), we are led to define as isoscalar instability the case when the state is unstable against isoscalar-like fluctuations, i.e. when $`1+F_{0g}^s<0`$ (or $`X<0`$). Analogously we deal with isovector instability when the system is unstable against isovector-like fluctuations i.e. when $`1+F_{0g}^a<0`$ (or $`Y<0`$).

One can easily see that in symmetric nuclear matter the isoscalar instability, ($`X<0,Y>0`$), appears as mechanical instability and the isovector instability, ($`X>0,Y<0`$), as chemical instability. Indeed now $`t=1`$, $`a=b`$, $`\beta =\pi /4`$ and so $`1+F_{0g}^s`$ (or $`X`$) and $`\left(\frac{P}{\rho }\right)_{T,y}`$ are proportional (see Eq.(26)). This is not true any longer for an asymmetric system.

To simplify the analysis let us assume for the moment that the quantities $`a`$ and $`b`$ remain positive in the density range we discuss here and so we can study the effect of the interaction between the two components, given by $`c`$, on the instabilities of the mixture. This is indeed the case of nuclear matter.

If $`c<0`$, i.e. for an attractive interaction between the two components, from Eq.(16) we see that the system is stable against isovector-like fluctuations but it becomes isoscalar unstable if $`c<2\sqrt{ab}`$ (see Eq.(14)). However thermodynamically this instability against isoscalar-like fluctuations will show up as a chemical instability if $`(tab/t)<c<2\sqrt{ab}`$ or as a mechanical instability if $`c<(tab/t)<2\sqrt{ab}`$ (see Eq.(26)). Therefore isoscalar-like instabilities are not necessarily equivalent to mechanical instabilities and may manifest instead as chemical instabilities.

If $`c>0`$, i.e. when the interaction between the components is repulsive, the thermodynamical state is always stable against isoscalar-like fluctuation, but can be isovector unstable if $`c>2\sqrt{ab}`$. Since the system is mechanically stable ($`a,b,c>0`$, see Eq. (26)), the isovector instability is now associated with chemical instability. Such situation will lead to a component separation of the liquid mixture.

Following this line a complete analysis of the instabilities of any binary system can be performed, in connection to signs, strengths and density dependence of the interactions.

We show now a quantitative calculation for asymmetric nuclear matter which illustrates the previous general discussion. We investigate the instabilities of $`ANM`$ characterized by a potential energy density of Skyrme-like type,

$`_{pot}(\rho _n,\rho _p)`$ $`=`$ $`{\displaystyle \frac{A}{2}}{\displaystyle \frac{(\rho _n+\rho _p)^2}{\rho _0}}+{\displaystyle \frac{B}{\alpha +2}}{\displaystyle \frac{(\rho _n+\rho _p)^{\alpha +2}}{\rho _0^{\alpha +1}}}`$ (27)
$`+`$ $`(C_1C_2({\displaystyle \frac{\rho }{\rho _0}})^\alpha ){\displaystyle \frac{(\rho _n\rho _p)^2}{\rho _0}}`$ (28)
where $`\rho _0=0.16\text{fm}^3`$ is the nuclear saturation density.

The values of the parameters $`A=356.8`$ MeV, $`B=303.9`$ MeV, $`\alpha =1/6`$, $`C_1=125`$ MeV, $`C_2=93.5`$ MeV are adjusted to reproduce the saturation properties of symmetric nuclear matter and the symmetry energy coefficient. We also extend the discussion to finite temperatures .

We first focus on the low density region, where phase transitions of liquid-gas type are expected to be seen in fragmentation events. Since $`a,b>0`$ and $`c<0`$, in the case of the liquid-gas phase transition in asymmetric nuclear matter we deal only with instability against isoscalar-like fluctuations as it happens for symmetric nuclear matter. However, at variance with the latter case, now the instability can manifest either as chemical or mechanical depending on the relative strength of the interactions in the system, as discussed previously.

In Figure 1 the circles represent the spinodal line corresponding to isoscalar-like instability, as defined above, for three values of the proton fraction. For $`y<0.5`$ under this border one encounters either chemical instability, in the region between the two lines, or mechanical instability, under the inner line (crosses). The latter is defined by the set of values ($`\rho ,T`$) for which $`\left(\frac{P}{\rho }\right)_{T,y}=0`$. We observe that the line defining chemical instability is more robust against the variation of the proton fraction in comparison to that defining mechanical instability: reducing the proton fraction it becomes energetically less and less favoured to break in clusters with the same initial asymmetry.

From Eq.s(23) and (26), we see that the quantity $`\left(\frac{\mu _p}{y}\right)_{T,P}`$ changes the sign passing through zero when we cross the spinodal line corresponding to chemical instability and passing through infinity when we cross the line associated with mechanical instability. However, from Eq.(23) (see also Fig. 3), we get that at the inner spinodal line, the ratio $`\left(\frac{\mu _p}{y}\right)_{T,P}/\left(\frac{\rho }{P}\right)_{T,y}XY`$ is a finite negative number. Therefore we conclude that the singularities of two quantities cancel by taking the ratio. Note also the smooth behaviour of the mixing angle $`\beta `$ through the mechanical instability line, Fig.2a.

Finally let us stress that, for $`ANM`$, even if we deal with instabilities against isoscalar-like fluctuations, we encounter a chemical effect, which gives raise to isospin distillation in equilibrium phase transitions . Indeed the variation of the asymmetry $`(I=12y)`$ is

$$\delta I=\nu _p[\frac{(1I_0)N_n}{(1+I_0)N_p}tg\beta 1]\nu _p[\chi 1]$$
(29)
where $`I_0`$ is the initial asymmetry. For $`y=0.5`$ we have $`\delta I=0`$ but for $`y<0.5`$, we find $`\chi <1`$ (see Fig.2b). Therefore $`\delta I<0`$ if $`\nu _p>0`$ and the opposite for $`\nu _p<0`$.

To be more general, we have extended our investigation also to the high density region. In Fig.3 we plot the density dependence of the generalized Landau parameters. We find that for the considered interaction the system exhibits another instability at high density, around $`1.5fm^3`$, where the quantity $`c`$ becomes positive and $`1+F_{0g}^a<0`$ (or $`Y<0`$).

Thermodynamically we see from Eq.(23) that this is again a chemical instability. However now it results from isovector-like fluctuations, in contrast to the low density instability. The reason is the change in the character of the interaction between the two components. Since the interaction becomes repulsive the nuclear phase can become unstable against proton-neutron separation. We also notice in Fig.3 that the generalized Landau parameters display a discontinuity where the quantity $`c`$ changes the sign. Other effective forces, with more repulsive symmetry terms, will not show this high density chemical instability , that actually could be of interest for other many body systems.

In conclusion, we have shown that in asymmetric binary systems the relevant instability regions are defined by the instabilities against isoscalar-like and isovector-like fluctuations. The kind of thermodynamical instability, chemical or mechanical, will depend on the relative strength between the various interactions acting in the system.

In particular, the liquid-gas phase transition in asymmetric nuclear matter results from instability against isoscalar-like modes rather than isovector-like due to the attractive character of the interaction between protons and neutrons. This is a qualitatively new effect which is leading to the observed ”isospin distillation” . More data are expected from the new Radioactive Beam facilities opening the possibility of direct studies of the charge dependent part of the nuclear interaction far from normal conditions.