text:
# Lévy Statistics for Random Single–Molecule Line Shapes in a Glass

## Abstract

We demonstrate that the statistical behavior of random line shapes of single tetra-tert-butylterrylene chromophores embedded in an amorphous polyisobutylene matrix at $`T=2`$ K, is described by Lévy statistics as predicted theoretically by Barkai, Silbey, and Zumofen \[Phys. Rev. Lett. 84 5339 (2000)\]. This behavior is a manifestation of the long range interaction between two level systems in the glass and the single molecule. A universal amplitude ratio is investigated which shows that the standard tunneling model assumptions are compatible with the experimental data.

Experimental advances MO  have made it possible to measure the fluorescence of a single molecule (SM) embedded in a glass. Because each individual molecule is in a unique static and dynamic environment, the fluorescence of chemically identical SMs varies from molecule to molecule. In this way the molecules serve as local reporters on the dynamics and statics of the host glass. SM experiments have been performed both in low-temperature glasses Zumb1 ; Zumb2 ; Orrit ; Elisa ; Donley2  and recently close to the glass transition temperature Desch1 . For low-temperature glasses the fundamental question is: Is the standard tunneling model valid for glasses? Related questions are how to analyze the complex line shape behaviors of SMs in glasses, and what do their random behaviors teach us on the SM-glass system.

The standard tunneling model AHV  was suggested in the early seventies to explain universal features of glasses; for example many glasses show a heat capacity which is nearly linear in temperature. At the center of this phenomenological model is the concept of the two-level system (TLS). It is assumed that at cryogenic temperatures excitations in glassy materials are two-level tunneling systems whose energies and tunneling matrix elements are randomly distributed. More recently, Geva and Skinner Geva  modeled behaviors of SMs in glasses based on the standard tunneling model. Orrit and co-workers Zumb1 ; Zumb2  used the fluorescence of single terrylene molecules in the polymer polyethylene to obtain the first direct experimental proof that two level systems actually exist in an amorphous material. More recent experiments Orrit  revealed behaviors not compatible with the standard model for $`21`$ out of $`70`$ molecules; e.g., a molecule coupled to a three level system REM . We note that a fundamental first-principle understanding of TLSs is still missing Wol , although numerical simulations HSilbey  give some evidence on the microscopic nature of a few of these entities.

There are many open questions concerning the standard model. For example: Are the TLSs uniformly distributed in space, or do they preferentially appear at boundaries of clusters of atoms/molecules, as suggested in Philips . Another open question is the nature of the interaction between the TLSs and the SM e.g., is it dipolar as suggested in Geva ? A method to obtain this important information using SM spectroscopy was suggested in Barkai  (see details below).

A theoretical investigation of the distribution of random line shapes of SMs in glass was carried out in Barkai ; BZS  based on the standard model approach Geva . Interestingly, the theoretical results obtained in Barkai  showed that the problem of random line shapes of SMs in glasses is related to Lévy statistics, thus the generalized Central Limit Theorem Feller  applies to this problem. This connection to Lévy statistics is a manifestation of long–range interactions between the TLSs in the glass and the SM (see details below). We note that Lévy statistics is known to describe several other long-range interaction models in diverse fields such as turbulence Min , and random magnetic systems Maass . Stoneham’s theory Stoneham  of inhomogeneous line broadening in crystals with defects is based on long-range forces, and it can be interpreted in terms of Lévy statistics BZS .

In this paper we analyze the statistical properties of random line shapes in a glass and compare it to the theoretical predictions in Barkai . In Fig. 1 we show eight lines of single tetra-tert-butylterrylene chromophores embedded in an amorphous polyisobutylene matrix at $`T=2`$ K. For experimental details see Naumov . The lines are typically multi-peaked similar to the numerical predictions in Geva . Each line is different from any other line, since each individual SM is in a unique environment. The multi-peaked behavior of the line shapes in Fig. 1 can be qualitatively explained using standard model arguments. If a SM is coupled to a single slow flipping TLS embedded in its vicinity, we expect that when the TLS flips from its up state to its down state or vice versa, the SM absorption frequency will shift. In this case the line shape of the SM is a doublet. It follows that the frequency of a SM coupled to $`N`$ such independent two-level systems will jump between $`2^N`$ states. Hence as shown in Fig. 1, the single molecule’s line shape will be composed of $`2^N`$ peaks. The width of these peaks depends on the dynamics of TLSs situated far from the molecule, fast dephasing processes, and the life time of the electronic transition Geva ; Barkai ; ADV ; Brown .

To obtain the lines in Fig. 1 we used the spectral trail technique introduced by Moerner and co-workers MT , following the spectral activity of the molecule during a scan time which in our experiments was fixed to be $`120`$ Sec. Following the jump history of the SM enables us to identify the peaks in Fig. 1 as originating from a SM. Without using this method it is practically impossible to say if the lines in Fig. 1 are due to contribution from several SMs or originate from a SM.

Random multi-peaked lines as shown in Fig. 1 are a novel feature of SM spectroscopy, not observed in ensemble averaging techniques. To obtain information on the SM-glass system we investigate the distribution of the line shapes of SMs in a glass. Mathematically we are dealing with the question of the distribution function of a function \[i.e., the line $`I(\omega _L)`$ is a function of laser frequency\]. Recently, Refs. Barkai ; BZS , suggested to characterize the line shape of each molecule by its cumulants $`\kappa _1,\kappa _2,\kappa _3,\mathrm{}`$. Thus each line is characterized by an infinite set of cumulants, which are random variables that vary from one molecule to the other, thus reflecting the disordered nature of the glass. The distribution functions $`P(\kappa _1),P(\kappa _2),P(\kappa _3),\mathrm{}`$ completely characterizes the statistical properties of the line shapes of SMs in a glass. The cumulants are obtained from the moments of the line shape Naumov2 , $`m_n=\omega _L^nI\left(\omega _L\right)d\omega _L`$, according to the well-known relations $`\kappa _1=m_1`$, $`\kappa _2=m_2(m_1)^2`$ etc. From $`244`$ SM spectra, we obtained the histograms of the first two cumulants, $`P(\kappa _1),P(\kappa _2),`$ which yield important information on the glass–SM system. We consider here the distribution of the cumulants and not the distribution of the moments, since the cumulants were predicted to be described by Lévy statistics Barkai .

To better understand the meaning of our data analysis and the relationship to Lévy statistics, it is useful to recall five main assumptions of the model used in Geva ; Barkai .
(i) The absorption frequency of the SM follows a stochastic trajectory described by $`\omega (t)=\omega _0+_{n=1}^N\xi _n(t)\nu _n(r)`$, where $`\omega _0`$ is the bare frequency of the molecule. $`\xi _n(t)`$ are random functions of time which follow a two-state process, $`\xi _n(t)=1`$ when the $`n`$-th TLS is in state up or $`\xi _n(t)=0`$ when it is in state down. Thus the flips of the TLSs induce spectral diffusion. The flipping rate between the up and down states is determined by a rate $`R`$ which varies from TLS to TLS \[the distribution of jumping rates of TLSs in a glass spans many orders of magnitude from nano seconds to (at least) days\]. (ii) The SM frequency shifts are: $`\nu _n=2\pi \alpha \mathrm{\Psi }\left(\mathrm{\Omega }\right)\frac{A}{E}\frac{1}{(r_n)^3}`$. The most important ingredient of the theory is the long-range interaction between the TLS and SM, $`\nu 1/r^3`$, reflecting the assumption of dipolar interaction. This long-range type of interaction is the first important ingredient in the relationship between SM spectroscopy in glasses and Lévy statistics. Other parameters controlling the frequency shifts are: $`\alpha `$ the SM-TLS coupling constant, $`E=\sqrt{A^2+J^2}`$ the energy splitting of the TLS, $`\mathrm{\Psi }\left(\mathrm{\Omega }\right)`$ a dimensionless function of order unity describing the orientation of the TLS and SM, and finally the random parameters of the TLS: $`A`$ (asymmetry parameter) and $`J`$ (tunneling matrix element). (iii) The TLSs are uniformly distributed is space and are non-interacting. This assumption is the second important condition for Lévy statistics to be valid. (iv) The standard tunneling model is valid; this model determines the distribution of parameters describing $`\nu _n`$ and $`R_n`$ as well as the density of TLSs. Note, however, that the Lévy statistics and results in Barkai  are not limited to this model. (v) The stochastic Kubo–Anderson theory of line shapes is applicable, implying weak laser fields. Under these conditions the following two limiting behaviors were found Barkai .

The first case corresponds to the fast–modulation limit $`\nu _nR_n`$ for all TLSs in the vicinity of the molecule. In this case, also called motional narrowing limit, the lines of individual molecules are Lorentzian in shape. Then the lines are characterized by two parameters only, e.g., the width at half maximum and the center location. The distributions of these two parameters are Lévy stable laws Barkai . From Fig. 1 it is clear that the fast modulation limit does not describe our experimental results.

The second case corresponds to the slow-modulation limit. If $`R_n\nu _n`$ for all the TLSs in the vicinity of the molecule, the shape of the line is random and typically multi–peaked Barkai . In this slow–modulation limit the distributions of line shape cumulants, $`P(\kappa _1),P(\kappa _2),P(\kappa _3),`$ etc, are all Lévy stable. Specifically, the probability density function of the first cumulant $`\kappa _1`$ is given by the symmetrical Lévy density, $`P(\kappa _1)=l_{1,0}(\kappa _1)`$, namely the Lorentzian

$$P(\kappa _1)=\frac{1}{\pi }\frac{z_1}{\kappa _1^2+z_1^2},$$
(1)
where $`z_1`$ is a scaling parameter which can be calculated from the theory in Barkai . In Fig. 2 we show that our experimental results are compatible with the theoretical prediction. We also fitted our results to a Gaussian probability density (not shown) and found that the distribution of the first cumulant is definitely not Gaussian.

We note that the reference frequency, determining the laser detuning in Fig. 1, was choosen on the maximum intensity of the spectrum of the SM $`\omega _{max}`$. While the theoretical reference frequency in Barkai  was the bare frequency of the molecule $`\omega _0`$. Using numerical simulation based on the approach in Geva , we observed that also when the reference frequency is choosen as $`\omega _{max}`$ distribution of first cumulant is well fitted by a Lorentzian, in agreement with the experimental results in Fig. 2. The value of $`z_1`$ in the two approaches slightly differs, as we discuss below. We now discuss the distribution of the second cumulant (i.e. the variance) which is not sensitive to definition of the reference frequency.

The distribution of the second cumulant $`\kappa _2`$ is given by the one-sided Lévy stable law, $`P(\kappa _2)=l_{1/2,1}(\kappa _2)`$, namely Smirnov’s one-sided probability density

$$P(\kappa _2)=\frac{1}{\left(z_{1/2}\right)^2}\frac{2}{\sqrt{\pi }}\left(\frac{2\kappa _2}{z_{1/2}^2}\right)^{3/2}\mathrm{exp}\left(\frac{z_{1/2}^2}{2\kappa _2}\right),$$
(2)
where the scaling parameter $`z_{1/2}`$ was derived in Barkai . As shown in Fig. 3, the experimental data for the distribution of $`\kappa _2`$ are compatible with the theoretical prediction; the long tail of the Lévy stable law is visible. Yet, data of a larger number of molecules is needed to improve the statistical fluctuations. The analysis of a larger number of molecules will also enable us to compare theory and experiments for the higher-order cumulants $`\kappa _3`$, $`\kappa _4`$ etc.

As mentioned, the Lévy behavior is due to the long-range dipolar interaction of a SM with many TLSs, hence the Gauss or Lévy Central Limit Theorem arguments are expected to hold. The Lévy central limit theorem applies since averaged frequency shifts diverge; $`\nu _0^{\mathrm{}}r^{d1}/r^3dr=\mathrm{}`$ where $`d=3`$ is the dimensionality of the problem (for mathematical details see Barkai ). Hence, the Lévy behavior obtained for the distributions $`P(\kappa _1)`$ and $`P(\kappa _2)`$ is used to test the assumptions of dipolar interactions and uniform-distribution of the TLSs in space. The information about the random distribution of the parameters of the glass, i.e. $`A`$ and $`J`$, are contained in the values of $`z_{1/2}`$ and $`z_1`$ which we now discuss.

There are two distinct types of parameters describing the SM-glass system: Those describing the bath of TLSs and the coupling constant $`\alpha `$ which depends on the properties of the SM probe. An important unsolved problem is how to extract information from SM experiments which is sensitive only to the statistical properties of the TLSs and is not affected by $`\alpha `$ Elisa . The scaling parameters $`z_1`$ and $`z_{1/2}`$ depend on $`\alpha `$, the magnitude of these two parameters depends also on the precise modeling of the orientation function $`\mathrm{\Psi }(\mathrm{\Omega })`$ Barkai . Thus $`z_1`$ and $`z_{1/2}`$ are not universal functions, in the sense that they depend on properties of the SM under investigation and not on the properties of the glassy state (which are supposed to be universal, according to standard tunneling model). In fact, Donley et al Elisa  suggested that the coupling constant $`\alpha `$ itself should be a random variable. This may seem to limit our ability to investigate low-temperature glasses with SM spectroscopy. If the statistical analysis of line shapes depends on an additional unknown distribution function of the coupling constant (besides the standard distributions of the glass parameters), fitting of data to the theory becomes rather arbitrary.

However, based on Eq. 9 in Barkai  one can show that the ratio $`z_{1/2}/z_1`$ depends only on the statistical properties of the glass and not on the distribution of the coupling constant $`\alpha `$. More precisely,

$$\frac{z_{1/2}}{z_1}=\frac{1}{\sqrt{2\pi }}\frac{\frac{A}{E}\text{Sech}\left(\frac{E}{2k_bT}\right)_{AJ}}{\frac{A}{E}\frac{1}{1+\mathrm{exp}\left(E/k_bT\right)}_{AJ}},$$
(3)
where the averaging is performed over the TLS parameters $`A`$ and $`J`$ Barkai . Since Eq. (3) is independent of the exact distribution of $`\alpha `$, it is a useful tool for describing the behavior of glasses. To derive Eq. (3) we assumed that the random variable $`\alpha `$ is independent of the glass parameters $`A`$ and $`J`$.

From our fits we find $`z_{1/2}=0.175`$ GHz and $`z_1=0.0485`$ GHz, which yields $`z_{1/2}/z_1=3.6`$. The theoretical prediction, based on Eq. (3) yields $`z_{1/2}/z_1=2.4`$. The deviation between theory and experiment is now explained. As mentioned, the theoretical prediction is based on the assumption that the bare frequency of the molecule $`\omega _0`$ is the reference frequency for measurement. Using numerical simulations Geva  with the parameter set relevant for our experiment (details will be published) we find that the ratio of $`z_1`$ when $`\omega _0`$ is the reference frequency and $`z_1`$ when $`\omega _{\mathrm{max}}`$ is the reference frequency, is $`z_1^{\omega _0}/z_1^{\omega _{max}}=1.6`$. Varying the value of $`\alpha `$ in our simulations, in the interval $`10\text{GHz}\text{nm}^3<\alpha <40\text{GHz}\text{nm}^3`$ we observed that the ratio $`z_1^{\omega _0}/z_1^{\omega _{max}}`$ does not depend on $`\alpha `$. As expected $`z_1^{\omega _0}>z_1^{\omega _{max}}`$ since the value of $`\kappa _1`$ becomes smaller (in statistical sense) if we assign the origin to the maximum of the spectrum. Using the correction factor $`z_1^{\omega _0}/z_1^{\omega _{max}}=1.6`$ we find that the theory yields $`z_{1/2}/z_1=3.8`$. Taking into account that the standard tunneling model does not address the chemical composition of the disordered system or the chemical and geometric details of the SM under investigation, we believe that the theoretical result is in surprisingly reasonable agreement with experiment. Measurements of the ratio $`z_{1/2}/z_1`$ for other types of SMs and glasses and for wider range of control parameters (i.e., temperature and scan time) will show whether SM data are compatible with the universal predictions of the standard model.

Our approach based on the analysis of the distributions of line shape cumulants is different from a second approach used in other studies Geva ; Zumb2 ; Elisa . Previously, the line was characterized by its width at half maximum and the distribution of the random line shapes was described by distribution of line widths. Clearly such a method does not capture the multi–peaked behavior of the SM lines shown in Fig. 1. An approach based on the distribution of line widths is useful in the fast-modulation limit, which as mentioned, (i) obeys Lévy statistics but (ii) is not relevant to the experiment under consideration.

To summarize SM spectroscopy is an excellent method to investigate disordered systems by removing the ensemble averaging found in other techniques. We have shown that our experimental results are compatible with Lévy statistics and with standard tunneling model predictions. In particular the following two assumptions are reasonable: (i) The two-level systems are uniformly distributed in space, (ii) The frequency shifts are caused by dipolar interactions $`\nu \frac{1}{r^3}`$. We introduced the universal ratio, $`z_{1/2}/z_1`$ which is sensitive to details of the standard model, but not to the coupling of the SM to the TLSs in the glass (i.e., not to $`\alpha `$). The comparison between the theoretical and the experimental value of this ratio can be used to test the validity of the standard model predictions.


format: text

text:
# Controlled creation of maximally entangled states using a SQUID ring coupled to an electromagnetic field

## Abstract

We present solutions of the time dependent Schrödinger equation for a SQUID ring coupled to an electromagnetic field, both treated quantum mechanically. We that show the SQUID ring can be used to create a maximally entangled state with the em field that without dissipation remains constant in time. Using methods familiar in quantum optics, we extend the model to include the effects of coupling this system to a dissipative environment. With this model we show that although such an environment makes a noticeable difference to the time evolution of the system, it need not destroy the the entanglement of this coupled system over time scales required for quantum technologies.

entanglement SQUID rings quantum information

Experiments to probe the quantum properties of superconducting Josephson weak link (SQUID) circuits Friedman et al. (2000); van der Wal et al. (2000); Rouse et al. (1995); Silvestrini et al. (2000); Nakamura et al. (1999) (in this work a superconducting ring enclosing a single weak link) have highlighted the possible use of such devices in emerging quantum technologies such as quantum computation, communication and information processing Orlando et al. (1999); Makhlin et al. (1999); Averin et al. (1990). In an earlier paper, motivated in part by experiment we solved the time dependent Schrödinger equation (TDSE) for an isolated SQUID ring modulated by a classical, monochromatic, microwave field Clark et al. (1998). Recent studies Everitt et al. (2001); Al-Saidi and Stroud (2002) of a fully quantum mechanical model of electromagnetic (em) fields coupled to SQUID rings have, from a theoretical standpoint, furthered our understanding of the manner in which such circuit may be utilised.

An essential feature revealed by these studies is that a quantum mechanical SQUID ring is highly non-perturbative in nature, a property that can be used to affect very strongly the non-linear dynamical behaviour of SQUID ring based quantum circuit systems. Furthermore, these effects can be controlled through the external (bias) magnetic flux $`\mathrm{\Phi }_x`$ applied to the ring. Even so, although several interesting new results came from the investigation of coupled ring-em field systems, in each case we chose, for ease of computation, a SQUID ring Hamiltonian description which was time independent. However, there are many circumstances where time dependence needs to be introduced into this description. In the work described here this arises because we wish to entangle a SQUID ring with an em field mode. As we will show, this can be achieved by applying a time dependent magnetic flux to the SQUID ring which changes slowly enough to maintain the ring-field system in its current state but fast enough to ensure strong entanglement. In this paper we solve the TDSE or, more generally, the master equation where dissipation involved, for tensor products of SQUID ring-em field systems with a particular time dependent potential. We show that we can use the flux dependent, non-perturbative nature of the SQUID ring to great utility in controlling its level of entanglement with an em field. As we shall also demonstrate, the level of this entanglement can decay with time or remain constant depending on whether or not dissipation (decoherence) is present.

In dealing with a SQUID ring coupled to an em field as quantum mechanical objects we represent the field as an equivalent circuit, a parallel, inductor-capacitor (LC), oscillator Everitt et al. (2001). We also assume that the coupling between the SQUID ring and the em field is inductive Likharev (1986). The schematic for this coupled circuit system is shown inset in 1.

With reference to the inset in figure 1, we can write down a Hamiltonian for the coupled two mode (SQUID ring + em field oscillator mode) system of the form

$$H=H_e+H_sH_{es}$$
(1)
where $`H_e`$ and $`H_s`$ are, respectively, the Hamiltonians for the em field and the SQUID ring and the interaction energy between them (for our purposes purely inductive) is $`H_{es}`$. For these we adopt the convention that operators with the subscripts $`s`$ or $`e`$ represent, respectively, those quantities associated with the SQUID ring or the em field. We also note that we use the first four energy eigenstates of each Hamiltonian to represent the individual circuit component of the system. This number was found to be sufficient to compute accurately the results presented in this paper.

The Hamiltonians for these quantum circuits can be expressed in terms of the conjugate variables of magnetic flux$`\mathrm{\Phi }`$ (coordinate) and electric displacement flux$`Q\left(i\mathrm{}/\mathrm{\Phi }\right)`$and the appropriate circuit parameters. Thus, for a SQUID ring (inductance $`\mathrm{\Lambda }_s`$, weak link capacitance $`C_s`$) $`\mathrm{\Phi }`$ is the total magnetic flux threading the ring and $`Q`$ is the total electric displacement flux (displacement charge) between the electrodes of the weak link; for the field oscillator these flux and charge variables relate, respectively, to the magnetic and electric displacement fluxes threading the equivalent inductor (inductance $`\mathrm{\Lambda }_e`$) and capacitor (capacitance $`C_e`$). In each case these operators satisfy the usual cannonical commutation relation $`[\mathrm{\Phi },Q]=i\mathrm{}`$. We assume that both the SQUID ring and em field oscillators (frequencies $`\omega _s/2\pi \left(=1/2\pi \sqrt{\mathrm{\Lambda }_sC_s}\right)`$ and $`\omega _e/2\pi \left(=1/2\pi \sqrt{\mathrm{\Lambda }_eC_e}\right)`$, respectively) operate in the quantum regime, i.e. such that $`\mathrm{}\omega _s,\mathrm{}\omega _ek_BT`$ for temperature $`T`$.

Given a SQUID ring inductance of $`\mathrm{\Lambda }_s`$ and a weak link capacitance of $`C_s`$, the SQUID ring Hamiltonian is, as usual, translated in external bias flux Clark et al. (1998) and takes the time dependent form

$$H_s=\frac{Q_s^2}{2C_s}+\frac{\mathrm{\Phi }_s^2}{2\mathrm{\Lambda }_s}\mathrm{}\nu \mathrm{cos}\left[\frac{2\pi }{\mathrm{\Phi }_0}\left(\mathrm{\Phi }_s+\mathrm{\Phi }_x(t)\right)\right]Q_s\frac{}{t}\mathrm{\Phi }_x(t).$$
(2)
where $`\mathrm{}\nu /2=0.0215\mathrm{\Phi }_0^2/\mathrm{\Lambda }_s`$ is the Josephson phase coherent coupling energy, $`\mathrm{\Phi }_0\left(=h/2e\right)`$ is the superconducting flux quantum and we note that the last term, $`Q_s\frac{}{t}\mathrm{\Phi }_x(t)`$ in (2) arises from our translating the Hamiltonian by $`\mathrm{\Phi }_x(t)`$ (see Clark et al. (1998) for more details). For this quantum regime we set $`C_s=1\times 10^{16}`$F and $`\mathrm{\Lambda }_s=3\times 10^{10}`$H, parameter values which are appropriate for a mesoscopic SQUID ring operating at a few K Everitt et al. (2001).

The Hamiltonian for the em field mode, which is somewhat simpler but equivalent to that for a simple harmonic oscillator, can be written as

$$H_e=\frac{Q_e^2}{2C_e}+\frac{\mathrm{\Phi }_e^2}{2\mathrm{\Lambda }_e}$$
(3)
while the energy associated with the inductive coupling between the ring and the em field mode is given by

$$H_{es}=\frac{\mu _{es}}{\mathrm{\Lambda }_s}\mathrm{\Phi }_s\mathrm{\Phi }_e$$
(4)
where $`\mu _{es}=0.01`$ is the ring-field magnetic flux linkage factor. For simplicity, and to facilitate the computations, we set $`C_e=C_s`$ and $`\mathrm{\Lambda }_e=\mathrm{\Lambda }_s`$ throughout this work.

Using (1), we showed in our previous work Everitt et al. (2001) that energy could be exchanged between the quantum modes of a SQUID ring-em field system around particular values of the static flux applied to the ring. For the SQUID ring ($`i=s`$) and the em field ($`i=e`$) the time averaged energy expectation values are found from

$$H_i=\underset{\tau \mathrm{}}{lim}\frac{1}{\tau }_0^\tau \psi \left(t\right)\left|H_i\right|\psi \left(t\right)𝑑t$$

Examples of exchange (or transition) regions in these time averaged energies are shown in figure 1 around the arrowed point A (i.e. where the energy difference between two of the ring energies is very close to $`n\mathrm{}\omega _e`$, $`n`$ integer). As in the rest of the paper, the initial state was taken to be $`|\psi \left(0\right)=|1_e0_s`$, i.e. with the em field mode in its first excited Fock state and the SQUID ring in its ground energy state. In computing these averaged energies we clearly could not integrate for the system over an infinite time; instead we calculated the time averaged expectation values for a time $`\tau `$ such that $`H_i`$ did not change noticeably if $`\tau `$ was further increased.

The actual value of $`\mathrm{\Phi }_x`$ around which exchange regions develop depends on circuit (ring and field oscillator) parameters and the initial condition of the system. For example, for the particular case presented in figure 1 the point $`A`$, at $`\mathrm{\Phi }_x=0.42864\mathrm{\Phi }_0`$ (together with its twin at $`\mathrm{\Phi }_x=0.57136\mathrm{\Phi }_0`$, as shown), is where significant energy is exchanged between the SQUID ring and the em field.

The energy exchange, which reaches a maximum at the centre of the exchange region, follows the strength of the coupling between the SQUID ring and the em field mode. Other subtle properties of this two component quantum circuit system also follow this changing coupling strength, in particular the quantum entanglement between the ring and the field mode. Entanglement plays a central role in considerations of quantum technologies Orlando et al. (1999); Makhlin et al. (1999); Averin et al. (1990) where it is often of the first importance that strong entanglement is maintained as the quantum system evolves with time. With this in mind, and as the central result in this paper, we set about trying to create a state of maximal entanglement that remains constant with time. In our chosen system this takes the form

$$|\psi \left(t\right)=\frac{1}{\sqrt{2}}\left(e^{i\phi _1}|1_e0_s+e^{i\phi _2}|0_e1_s\right)$$
(5)
To implement this we now make the external flux time dependent with a functional form (see inset in figure 2 (b)) given by

$$\mathrm{\Phi }_x\left(t\right)=\{\begin{array}{cc}A\hfill & tt_0\hfill \\ \frac{\left(BA\right)}{t_r}\left(tt_0\right)\hfill & t_0<tt_0+t_r\hfill \\ B\hfill & t>t_0+t_r\hfill \end{array}$$
(6)
where $`t_0=326\omega _s^1`$ and the ramp time is $`t_r=16.6\omega _s^1`$. Following this change in bias flux, the operating point has moved from $`A`$ to $`B`$ $`\left(\mathrm{\Phi }_x=0.38\mathrm{\Phi }_0\right)`$ in figure 1, well outside the exchange region. The cut off time $`t_0`$ is the time at which the probability of the system being in state $`|1_e0_s`$ is equal to the probability of it being in state $`|0_e1_s`$. For our two component system we found these probabilities by solving the TDSE with the ring and field mode parameters given above. As there is no exchange of energy between the ring and the field mode at point $`B`$, it seemed reasonable to assume that at this flux bias the system would maintain equal probability of being in each of the states $`|1_e0_s`$ and $`|0_e1_s`$ at times after $`t_0`$. The results of our calculations of these probabilities are presented in figure 2(a), plotted against dimensionless time $`\omega _et`$. We see from these computed solutions that our assumption is essentially correct, at least for the particular system example being studied.

Given this prediction of the behaviour of the circuit system, it is important to establish a measure of the level of entanglement between its two components. In this work we adopt the entropic measure of entanglement first introduced by Adami and Cerf Cerf and Adami (1998), i.e.

$$I_i=S\left(\rho \right)S\left(\rho _i\right)$$

where $`S\left(\rho \right)=Tr\left(\rho \mathrm{ln}\rho \right)`$ and $`\rho `$ is the density operator. If $`I_i`$ is negative this implies that the $`i^{\mathrm{th}}`$ component of the system is entangled with the rest of the system, an inequality that holds even in the presence of dissipation. This computed quantity for our system is shown in figure 2 (b) as a function of $`\omega _et`$. As we might suspect, $`t_0`$ is also the time when the entanglement between the em field and the SQUID ring reaches a maximum. From the figure we see that after approximately a time $`t_0`$ both the ring and field mode probabilities have reached values of 0.5 (figure 2(a)) which then remain constant with time. Similarly, in 2(b) beyond $`t_0`$ the entanglement is strong (at 0.7) and also constant with time.

In the above discussion we have seen from solutions of the TDSE that we can use a time dependent flux (6), applied to the SQUID ring, to make a maximally entangled state of the form (5). Nevertheless, for this entanglement to be of use technologically it must be maintained over some reasonable period of time in the presence of dissipation. In order to estimate the effects of dissipation we adopt the open systems approach usually used in the field of quantum optics. In this approach the evolution of the system is found by solving a master equation for the density operator. The master equation for a system subject to thermal dissipation is given by

$`{\displaystyle \frac{}{t}}\rho `$ $`=[H,\rho ]+{\displaystyle \underset{i=e,s}{}}{\displaystyle \frac{\gamma _i}{2\mathrm{}}}\left(M_i+1\right)\left(2a_i\rho a_i^{}a_i^{}a_i\rho \rho a_i^{}a_i\right)`$
$`+{\displaystyle \frac{\gamma _i}{2\mathrm{}}}M_i\left(2a_i^{}\rho a_ia_ia_i^{}\rho \rho a_ia_i^{}\right)`$
where $`\gamma _i`$ is the damping rate of each system to the thermal bath and the mean photon number $`M_i=2.746\times 10^5`$ in each component is related to the temperature $`T_{b_i}=4.2`$K and frequency $`\omega _{b_i}=\omega _s`$ of each decohering bath via $`M_i=\left[\mathrm{exp}\left(\mathrm{}\omega _{b_i}/k_BT\right)1\right]^1`$.

The effect of dissipation on the probabilities and entanglement entropies of our chosen system are shown in figures 3 and 4 for the cases of weakly and strongly dissipative environments, respectively. Thus, in figure 3 at a $`\gamma `$ value of $`1\times 10^5\omega _s`$, the effect on the probablities is negligible over a few $`t_0`$ while the entanglement can be seen to weaken only gradually over the same time period. For this degree of decoherence the corrsponding decoherence time is of the order of 20ns, which is in line with recent experimental results on superconducting weak link circuits Martinis et al. (2002) (albeit at much larger weak link capacitances ($`10^{12}`$F) and lower temperatures ($`25`$mK)). This suggests that at this level of dissipation the coherence of the coupled system could be used effectively for quantum coherent circuit operations. It is also apparent from the second example of the effect of a dissipative environment that substantially greater values of $`\gamma \left(=1\times 10^4\omega _s\right)`$ i.e. stronger dissipation, still do not entirely suppress entanglement over this same time scale of a few $`t_0`$. As can be seen, at this higher level of dissipation our system starts to become mixed after a sufficient evolution in time. Concomitantly, a time independent, maximally entangled, state is not maintained. Clearly, this result suggests a practical limit to the level of dissipation (decoherence) that can be tolerated if quantum coherent operations are to be performed successfully. With this in mind, we note that these results do not necessarily imply that dissipative effects will be significant, although when designing real circuits environmental effects must be accounted for in its design.

In this paper we have shown that, via a time dependent external bias flux, a SQUID ring can be used to control the quantum state of a coupled em field mode-SQUID ring system, including the degree of entanglement between these two circuit components. As a corollary of this, we showed that the level of this entanglement could be set and maintained, with only a slow degradation in time, when the system was coupled to a weakly dissipative external environment. With the current interest in using superconducting weak link circuits for quantum technologies, this work provides some guide to the time scales over which it may be possible to perform quantum coherent operations in the presence of decohering environments.

###### Acknowledgements.

We would like to thank Dr. T.P. Spiller for interesting discussions on entanglement in quantum circuits.


text:
# 1 Nucleon FFs in the SL region. The magnetic FF for proton (a) and for neutron (c) are normalized to one and divided by 𝐺_𝐷. The proton electric FF is shown as the ratio 𝜇𝐺_{𝐸𝑝}/𝐺_{𝑀𝑝}(b), the dotted line is from (). The neutron electric FF (c) is plotted with the () parametrization (dash-dotted line). The solid line is from (), the dashed line from ().

QCD-view on hadron form factors in space-like and time-like momentum transfer regions

Egle Tomasi-Gustafsson and Michail P. Rekalo 111 Permanent address: NSC-Kharkov Physical Technical Institute, 61108 Kharkov, Ukraine

DAPNIA/SPhN, CEA/Saclay, 91191 Gif-sur-Yvette Cedex, France

Abstract

QCD gives definite predictions for hadron electromagnetic form factors in space-like and time-like momentum transfer regions, such as the quark counting rules, the hypothesis of hadron helicity conservation, and the relations between nucleon and deuteron form factors in the formalism of reduced nuclear matrix elements. Recent precise data about these form factors, obtained in polarization experiments at the Jefferson Laboratory, have essentially changed our view on this subject. QCD-predictions do not apply to these data up to $`Q^2`$=5-6 GeV2 for deuteron and for the electric form factor of proton. An analysis of these data suggests that the asymptotic region will more probably start at $`Q^2`$=20-25 GeV2. We show that the separation of magnetic and electric proton form factors in the time-like region represents the most stringent test of the asymptotic regime and QCD-predictions.

In this talk we will discuss the recent developments in the field of hadron electromagnetic form factors (FFs), due to the very precise and surprising data obtained at the Jefferson Laboratory (JLab), in $`\stackrel{}{e}+pe+\stackrel{}{p}`$ elastic scattering.

The application of the polarization transfer method, proposed about 30 years ago Re68  has been possible only recently, as it needs high intensity polarized beams, large solid angle spectrometers and advanced techniques of polarimetry Mi98 ; Jo00 ; Ga02  in the GeV range.

At these energies, one probes distances of the order of the nucleon size or less, and would expect expect the manifestation of quark degrees of freedom.

FFs, which characterize the internal structure of composite particles, can be experimentally measured and theoretically calculated, thus providing a good playground for the models of nucleon structure.

In a P and T invariant theory, a particle with spin $`S`$ is characterized by $`2S+1`$ electromagnetic form factors, which are complex (real) functions in the time(space)-like region. The nucleon has two FFs, called electric $`G_{EN}`$ and magnetic $`G_{MN}`$, which, a priori, are different.

Before the new, precise data appeared in the space-like (SL) region, the dipole approximation, $`G_D=\left[1+Q^2/a\right]^2`$ with $`a=0.71\text{GeV}^2`$, has been considered a good approximation for three of the four nucleon FFs, $`G_{En}`$being neglected or taken according to Galster .

In the Breit system, FFs are related to the Fourier transform of the charge and magnetic moment distributions and the dipole approximation results from an exponential distribution, the coefficient $`a=0.71`$ GeV2 corresponding to a root mean squared radius $`\sqrt{<r^2>}`$=0.81 fm.

Perturbative QCD gives definite rules about the scaling behavior of the form factors and about helicity conservation, indicating where quark degrees of freedom should be taken explicitely into account. Scaling laws, which give the probability that a hadron remains intact after absorbing a photon of momentum $`Q^2`$, have been formulated in Br76  as $`F=\left[1+Q^2/(n\beta _n^2)\right]^{n1}`$, where $`n`$ is the number of quarks and $`\beta _n`$ the quark average momentum. When applied to proton, they give the same power law as dipole, including the coefficient $`a=0.71`$, derived from a fit to pion form factors.

The traditional way to access FFs in the SL region, is the elastic scattering of electrons on hadrons $`e+he+h`$, and in the time-like (TL) region, the annihilation processes $`e^++e^{}p+\overline{p}`$. The unpolarized elastic $`ep`$ cross section, in one-photon exchange approximation, can be written as a function of the electric $`G_{Ep}`$and magnetic $`G_{Mp}`$proton form factors Ro50 :

$$\frac{d\sigma }{d\mathrm{\Omega }}(epep)=\left(\frac{d\sigma }{d\mathrm{\Omega }}\right)_{Mott}\left[\frac{G_{Ep}^2+\tau G_{Mp}^2}{1+\tau }+2\tau G_{Mp}^2\mathrm{tan}^2\frac{\theta _e}{2}\right],\tau =\frac{Q^2}{4M^2}$$
(1)
with

$$\left(\frac{d\sigma }{d\mathrm{\Omega }}\right)_{Mott}=\frac{\alpha ^2\mathrm{cos}^2(\theta _e/2)E^{}}{4E^3\mathrm{sin}^4(\theta _e/2)}\text{ and }E^{}=\frac{E}{1+2{\displaystyle \frac{E}{M}}\mathrm{sin}^2(\theta _e/2)}$$

where $`M`$ is the proton mass, $`E`$ is the energy of the incident electron, $`E^{}`$ and $`\theta _e`$ are the energy and scattering angle of the outgoing electron, $`\alpha `$ is the fine structure constant. The momentum of the virtual photon is $`Q^2=4EE^{}\mathrm{sin}^2(\theta _e/2)`$ and it is positive in the SL region.

In the TL region the cross section can be expressed as a function of FFs according to the following formula Zi62 :

$$\frac{d\sigma }{d(cos\theta )}(e^+e^{}p\overline{p})=\frac{\pi \alpha ^2}{8M^2\tau \sqrt{\tau (\tau 1)}}\left[\tau |G_{Mp}|^2(1+\mathrm{cos}^2\theta )+|G_{Ep}|^2\mathrm{sin}^2\theta \right],$$
(2)
where $`\theta `$ is the angle between the electron and the antiproton in the center of mass frame.

Eqs. (1) and (2), contain the moduli squared of the FFs, therefore one can not access their sign. Moreover the contribution of $`G_{Mp}`$appears weighted by the factor $`\tau `$: as $`Q^2`$ increases, it becomes the dominant term, making the extraction of $`G_{Ep}`$very imprecise.

In the SL region the collisions of polarized electron and polarized target (or measuring the polarization of the scattered proton), induce an interference term proportional to the product $`G_{Ep}G_{Mp}`$. It is therefore more sensitive to a small contribution of $`G_{Ep}`$.

Data have been obtained at JLab up to $`Q^2=5.6`$ GeV2 and an extension up to 9 GeV2 is in preparation 00111 . The measurement of the polarization of protons at momentum up to $`P_p=5.3`$ GeV/c can be done with a POMME-like polarimeter, as it has been proved at the JINR-LHE, in Dubna dubnacalib . Although the analyzing power for the inclusive reaction $`p+CH_2\text{one  charged particle }+X`$ decreases with increasing incident momentum, it is still sizeable at a proton momentum of 5.3 GeV.

The recent data Jo00 ; Ga02  show a linear deviation from the dipole behavior and can be parametrized by $`\mu `$ $`G_{Ep}`$/ $`G_{Mp}`$$`=1.00.130(Q^2/\text{GeV}^20.04)`$ Ga02  ($`\mu `$ is the proton magnetic moment), definitely proving that the electric and magnetic distributions in the proton are different.

Although different models can reproduce this trend, (as constituent quark models Ca00 , soliton model Ho96 , diquark model Kroll  ..) few of them give a satisfactory description of all four electromagnetic form factors and even fewer of them are applicable in SL and TL regions. We would like to quote here the models based on vector meson dominance (VMD)GK92 . The most recent one Lomon  contains several parameters to fit the world data and includes the asymptotic QCD behavior. Note that ’Il Nuovo Cimento’ reported the first studies, in this field. In Ref. Wataghin  the right trend for the $`G_{Ep}`$was already proposed, and in Zi66  a ’one parameter fit’ gave a good qualitative description of the present data and prescriptions for TL region as well. The best VMD-predictions, concerning the proton FFs, are given by IJL  (Fig. 1). A compilation of the world data reported in the figure can be found in Lomon .

The nucleon FFs are important ingredients for the calculations of the light nuclei structure. One of the consequences of the new data is the revision of the models of the deuteron structure. Following Br76 , let us introduce a generalized deuteron FF, $`F_D(Q^2)`$, $`F_D(Q^2)=\sqrt{A(Q^2)}`$, where $`A(Q^2)`$ is the structure function related to the forward electron deuteron cross section, and a reduced deuteron FF $`f_D(Q^2)`$:

$$f_D(Q^2)=\frac{F_D(Q^2)}{F_N^2(Q^2/4)},$$
(3)
where $`F_N`$ is the nucleon electromagnetic FF. The $`Q^2`$-behavior of $`f_D(Q^2)`$ (at large $`Q^2`$) can be predicted in the framework of pQCD, in the following form:

$$f_D(Q^2)=N\frac{\alpha _s(Q^2)}{Q^2}\left(ln\frac{Q^2}{\mathrm{\Lambda }^2}\right)^\mathrm{\Gamma },$$
(4)
where N is the normalization factor (which is not predicted by QCD), $`\alpha _s`$ is the running QCD strong interaction coupling constant, $`\mathrm{\Lambda }`$ is the scale QCD parameter, and $`\mathrm{\Gamma }`$ is determined by the leading anomalous dimension, here $`\mathrm{\Gamma }={\displaystyle \frac{8}{145}}`$.

In Al99  it was shown that the QCD prediction (4), which can be applied to asymptotic momentum transfer, is working well already for $`Q^22`$ GeV2, with a plausible value of the parameter $`\mathrm{\Lambda }`$ 100 MeV, in agreement with the values determined by many other possible methods pdg .

In Br76  another interesting prediction, concerning the scaling behavior of the reduced deuteron FF, was done:

$$f_R=\left(1+\frac{Q^2}{m_0^2}\right)f_D(Q^2)const,$$
(5)
where $`m_0^2=0.28`$ GeV2 is a parameter related to the $`Q^2`$-behavior of the pion FF. The same data from Al99 , if plotted in the representation of the reduced deuteron FFs, should illustrate the $`Q^2`$-independence of this product.

This result was confirmed by the previous $`A(Q^2)`$ data Ar86 , in the limit of their accuracy. In Fig. (2) we show that the new, more precise data about $`A(Q^2)`$ Al99 , are not consistent with the prediction (5) as they show an evident dependence of the product $`f_R`$ on $`Q^2`$, even for $`Q^2`$2 GeV2. This behavior can not be changed by varying the parameter $`m_0`$.

Although the scaling laws seem to be consistent with cross section measurements, up to 6 GeV2, if one replaces the dipole approximation with other descriptions of the nucleon FFs, taking into account the deviation from dipole for $`G_{Ep}`$, a fit following Eqs. (3) and (4) shows a large instability for $`\mathrm{\Lambda }`$ MPR03 .

One should mention that a calculation based on the impulse approximation, where the nucleon FF are taken from GK92 , satisfactorily reproduces the existing data on the three deuteron FFS ETG01a .

We will go a step further, and look for a definition of the asymptotic region with respect to the analyticity properties of complex functions.

FFs must obey the Phràgmen-Lindelöf theorem Ti39 , which gives a rigorous prescription for the asymptotic behavior of analytical functions: $`lim_{Q^2\mathrm{}}F^{(SL)}(Q^2)=lim_{Q^2\mathrm{}}F^{(TL)}(Q^2)`$. This means that, asymptotically, the TL phase vanishes and the real part of the FFs, $`eF^{(TL)}(Q^2)`$, coincides with the corresponding value $`F^{(SL)}(Q^2)`$.

The Rosenbluth separation has not yet been realized in the TL region. In order to extract the FFs, due to the poor statistics, it is necessary to integrate the differential cross section over a wide angular range. One assumes that the $`G_{Ep}`$-contribution plays a minor role in the cross section and the experimental results are usually given in terms of $`|G_{Mp}|`$, under the hypothesis that $`G_{Ep}`$=0 or $`|G_{Ep}|=|G_{Mp}|`$. The first hypothesis is arbitrary. The second hypothesis is strictly valid at threshold only, but there is no theoretical argument which justifies its validity at any other momentum transfer, where $`s4M^2`$. The prediction for the TL region in Zi62  shows that if $`G_{Ep}`$/$`G_D`$ decreases, $`G_{Mp}`$/$`G_D`$ increases, due to the definitions (dashed and solid line, respectively in Fig. 3).

Therefore a comparison of data in TL and SL region should give an unambiguous indication on the asymptotic region. The experimental data are shown in Fig. 3, normalized to the function $`G_D`$. For a compilation see ETG01 , here updated with recent data An03  (open squares).

The values of $`G_{Mp}`$in the TL region, obtained under the assumption that $`|G_{Ep}|=|G_{Mp}|`$ (open symbols), are larger than the corresponding SL values (solid squares and solid triangles). This has been considered as a proof of the non applicability of the Phràgmen-Lindelöf theorem, or as an evidence that the asymptotic regime is not reached Bi93 .

The magnetic form factor of the proton in the TL region (which is deduced from the hypothesis $`G_{Ep}`$=0 (case 1) or $`G_{Ep}`$=$`G_{Mp}`$(case 2), can be parametrized as: $`G_M^{(TL)}=G_Da/\left(1+s/m_{nd}^2\right)`$, where $`a`$ is a normalization parameter and $`m_{nd}^2=3.6\pm 0.9`$ GeV2 characterizes the deviation from the usual dipole $`s`$-dependence. The extrapolation to higher $`s`$ based on this formula (Fig. 3, dotted line), indicates that the Phràgmen-Lindelöf theorem will be satisfied by this FF, only for $`s(Q^2)20`$ GeV2.

Let us assume now that one of the two proton electromagnetic FFs has reached the asymptotic regime and apply the Phràgmen-Lindelöf theorem to extract the other. This looks as a reasonable hypothesis for $`G_M`$, which shows an early scaling behavior, in accordance with quark counting rules. From Eq. 2 we can deduce $`|G_E|`$, using the existing experimental data about $`\overline{p}+pe^++e^{}`$ (case 3). We report, in Fig. 4, the recent data in TL region, reanalized following the possibilities suggested above. Fig. 4a shows the values of the form factors taking $`G_E=0`$ (circles) and $`|G_E|=|G_M|`$ (squares) respectively). For case 3, where $`G_M=G_D`$, the values of $`|G_E|`$ (triangles) are larger than in cases 1 and 2. This suggests that asymptotics are not reached for $`G_E`$, as the values in the SL and TL regions get more apart. A fourth possibility is taking for $`G_E`$ in the TL region the values from Jo00 ; Ga02  and calculate $`|G_M|`$ (case 4). This affects very little the values of $`G_M`$, due to the kinematical factor $`\tau `$, which weights the magnetic contribution to the differential cross section (stars).

One can express the angular dependence of the differential cross section as a function of the angular asymmetry $`𝒜`$ as:

$$\frac{d\sigma }{d(\mathrm{cos}\theta )}=\sigma _0\left[1+𝒜\mathrm{cos}^2\theta \right],\text{ with }𝒜=\frac{\tau |G_M|^2|G_E|^2}{\tau |G_M|^2+|G_E|^2},$$
(6)
where $`\sigma _0`$ is the value of the differential cross section at $`\theta =\pi /2`$. Fig. 4b shows the angular asymmetry for the different cases. Case 1 and case 2 give, respectively, $`𝒜=1`$ and $`𝒜=(\tau 1)/(\tau +1)`$. The calculated asymmetries are very sensitive to the different underlying assumptions, therefore a precise measurement of this quantity would be very interesting.

Finally, we note that the angular dependence of the cross section, Eq. (2), results directly from the assumption of one-photon exchange, where the spin of the photon is equal 1 and the electromagnetic hadron interaction satisfies the $`C`$invariance. Therefore the measurement of the differential cross section at three angles (or more) would also allow to test the presence of $`2\gamma `$ exchange. The relative role of the $`2\gamma `$ mechanism can increase at relatively large momentum transfer in SL and TL regions, for the same physical reasons, which are related to the steep decreasing of the hadronic electromagnetic FFs, as recently discussed in Re99 .

…instead of Conclusions

It is for me a great pleasure and honour to express my gratitude to Prof. R. A. Ricci, to whom this Conference is dedicated. As supervisor of my ’Tesi di Laurea, Università di Padova, Anno Accademico DCCLVI’, he gave me the basis to progress in the field of Nuclear Physics. Since that time, his enthusiasm, his dynamic personality, his optimism in affording a variety of problems of different nature, his broad interests and the ability to listen, to understand and synthesize information very rapidly, have always been for me a source of motivation and example. I take this opportunity to thank Prof. Ricci, and to wish to him, besides the recognizements due by the italian and world-wide nuclear physics community, to fully enjoy a time (may be) more quiet and to have a lot of satisfaction in his family life.

E. T.-G.


text:
# The Majorana Neutrino Mass Matrix with One Texture Zero and One Vanishing Eigenvalue

## Abstract

Possible patterns of the Majorana neutrino mass matrix $`M`$ with one texture zero and one vanishing eigenvalue are classified and discussed. We find that three one-zero textures of $`M`$ with $`m_1=0`$ and four one-zero textures of $`M`$ with $`m_3=0`$ are compatible with current neutrino oscillation data. The implications of these phenomenological ans$`\ddot{\mathrm{a}}`$tze on the neutrino mass spectrum and the neutrinoless double beta decay are also explored in some detail.

The recent SK , SNO , KamLAND  and K2K  neutrino oscillation experiments have provided us with very convincing evidence that neutrinos are massive and lepton flavors are mixed. At low energy scales, the phenomenology of neutrino masses and lepton flavor mixing can be formulated by the effective neutrino mass matrix $`M`$ in the flavor basis where the charged lepton mass matrix is diagonal. Note that $`M`$ totally involves nine real parameters: three neutrino masses ($`m_1,m_2,m_3`$), three flavor mixing angles ($`\theta _x,\theta _y,\theta _z`$), one Dirac-type CP-violating phase ($`\delta `$) and two Majorana-type CP-violating phases ($`\rho `$ and $`\sigma `$). Although we have achieved some quantitative knowledge on two neutrino mass-squared differences and three flavor mixing angles from current neutrino oscillation experiments, much more effort is needed to determine these parameters to a high degree of accuracy. The most challenging task is to pin down the absolute neutrino mass scale and the CP-violating phases. The latter are entirely unrestricted by today’s experimental data.

However, the presently conceivable sets of feasible neutrino experiments are only possible to determine seven of the nine independent parameters of $`M`$ . This observation implies that some phenomenological hypotheses have to be made for the texture of $`M`$, in order to establish testable relations between the experimental quantities and the parameters of $`M`$. In Ref. , the patterns of $`M`$ with two texture zeros have been classified and discussed. In Ref. , particular attention has been paid to the possibility that $`M`$ may have one vanishing eigenvalue (i.e., $`m_1=0`$ or $`m_3=0`$). Both approaches allow us to determine or constrain the values of neutrino masses and Majorana phases in terms of the relevant experimental observables.

The purpose of this paper is to consider another simple but interesting approach for the study of neutrino masses and lepton flavor mixing at low energy scales. We assume that the Majorana neutrino mass matrix $`M`$ possesses one texture zero and one vanishing eigenvalue. Then we find that three one-zero patterns of $`M`$ with $`m_1=0`$ (normal hierarchy) and four one-zero patterns of $`M`$ with $`m_3=0`$ (inverted hierarchy) are compatible with current neutrino oscillation data. The implications of these phenomenological ans$`\ddot{\mathrm{a}}`$tze on the neutrino mass spectrum and the neutrinoless double beta decay are also discussed in some detail.

In the flavor basis where the charged lepton mass matrix is diagonal, the neutrino mass matrix $`M`$ can be written as

$$M=V\left(\begin{array}{ccc}m_1& 0& 0\\ 0& m_2& 0\\ 0& 0& m_3\end{array}\right)V^T,$$
(1)
where the lepton flavor mixing matrix $`V`$ can be parametrized as

$$V=\left(\begin{array}{ccc}c_xc_z& s_xc_z& s_z\\ c_xs_ys_zs_xc_ye^{i\delta }& s_xs_ys_z+c_xc_ye^{i\delta }& s_yc_z\\ c_xc_ys_z+s_xs_ye^{i\delta }& s_xc_ys_zc_xs_ye^{i\delta }& c_yc_z\end{array}\right)\left(\begin{array}{ccc}e^{i\rho }& 0& \\ 0& e^{i\sigma }& 0\\ 0& 0& 1\end{array}\right)$$
(2)
with $`s_x\mathrm{sin}\theta _x`$, $`c_x\mathrm{cos}\theta _x`$ and so on. It is obvious that $`M`$ depends on nine real parameters. Allowing one entry of $`M`$ to vanish, i.e., $`M_{\alpha \beta }=0`$, we obtain

$$\underset{i=1}{\overset{3}{}}\left(m_iV_{\alpha i}V_{\beta i}\right)=\mathrm{\hspace{0.33em}0},$$
(3)
where $`\alpha `$ and $`\beta `$ run over $`e`$, $`\mu `$ and $`\tau `$. There are totally six possible one-zero textures of $`M`$, illustrated in Table 1 as patterns A, B, C, D, E and F. Pattern A and its consequences have somehow been discussed in Ref. .

Now we consider whether one of three neutrino masses can vanish. Note that the solar and atmospheric neutrino oscillations probe the mass-squared differences

$`\mathrm{\Delta }m_{\mathrm{sun}}^2`$ $``$ $`\left|m_2^2m_1^2\right|,`$ (4)
$`\mathrm{\Delta }m_{\mathrm{atm}}^2`$ $``$ $`\left|m_3^2m_2^2\right|,`$ (5)
respectively. As the present solar neutrino data support $`m_1<m_2`$ , $`m_2`$ must not vanish. Thus we are left with two possibilities : either $`m_1=0`$ (normal hierarchy) or $`m_3=0`$ (inverted hierarchy).

(1) If $`m_1=0`$ holds, Eq. (3) can be simplified to

$$\xi \frac{m_2}{m_3}=\frac{|V_{\alpha 3}V_{\beta 3}|}{|V_{\alpha 2}V_{\beta 2}|},$$
(6)
and

$$\mathrm{arg}\left(\frac{V_{\alpha 3}V_{\beta 3}}{V_{\alpha 2}V_{\beta 2}}\right)=\pm \pi .$$
(7)
Eq. (6) implies a strong constraint on the CP-violating phases hidden in the matrix elements of $`V`$. The ratio of $`\mathrm{\Delta }m_{\mathrm{sun}}^2`$ and $`\mathrm{\Delta }m_{\mathrm{atm}}^2`$ can be calculated by use of Eqs. (4) and (5). We obtain

$$R_\xi \frac{\mathrm{\Delta }m_{\mathrm{sun}}^2}{\mathrm{\Delta }m_{\mathrm{atm}}^2}=\frac{\xi ^2}{\left|1\xi ^2\right|}.$$
(8)
As $`R_\xi 10^2`$ is required by current experimental data, $`\xi ^210^2`$ must hold. Hence one arrives at a normal neutrino mass hierarchy.

(2) If $`m_3=0`$ holds, Eq. (3) can be simplified to

$$\zeta \frac{m_1}{m_2}=\frac{|V_{\alpha 2}V_{\beta 2}|}{|V_{\alpha 1}V_{\beta 1}|},$$
(9)
and

$$\mathrm{arg}\left(\frac{V_{\alpha 2}V_{\beta 2}}{V_{\alpha 1}V_{\beta 1}}\right)=\pm \pi .$$
(10)
Again, Eq. (9) implies a strong constraint on the CP-violating phases hidden in the matrix elements of $`V`$. The ratio of $`\mathrm{\Delta }m_{\mathrm{sun}}^2`$ and $`\mathrm{\Delta }m_{\mathrm{atm}}^2`$ in this case is found to be

$$R_\zeta \frac{\mathrm{\Delta }m_{\mathrm{sun}}^2}{\mathrm{\Delta }m_{\mathrm{atm}}^2}=\mathrm{\hspace{0.33em}1}\zeta ^2.$$
(11)
Because both $`\zeta <1`$ and $`R_\zeta 10^2`$ are required, only a very narrow range of $`\zeta `$ (even smaller than the range $`0.95<\zeta <1`$) can satisfy Eq. (10).

The explicit expressions of $`\xi `$ and $`\zeta `$ for each one-zero texture of the neutrino mass matrix $`M`$ are given in Table 1 in terms of the flavor mixing parameters $`\theta _x`$, $`\theta _y`$, $`\theta _z`$ and $`\delta `$. In view of the recent SK , SNO , KamLAND , K2K  and CHOOZ  data on neutrino oscillations, we have $`\mathrm{\Delta }m_{\mathrm{sun}}^2[5.9,8.8]\times 10^5\mathrm{eV}^2`$, $`\mathrm{sin}^2\theta _x[0.25,0.40]`$ ; $`\mathrm{\Delta }m_{\mathrm{atm}}^2[1.65,3.25]\times 10^3\mathrm{eV}^2`$, $`\mathrm{sin}^22\theta _y[0.88,1.00]`$ ; and $`\mathrm{sin}^22\theta _z<0.2`$ at the $`90\%`$ confidence level. Then we arrive at $`R_\xi =R_\zeta [0.018,0.053]`$ as well as $`\theta _x[30.0^{},39.2^{}]`$, $`\theta _y[34.9^{},55.1^{}]`$ and $`\theta _z[0.0^{},13.3^{})`$. These experimental results allow us to examine which pattern of $`M`$ is phenomenologically favored, although the Dirac phase $`\delta `$ is entirely unrestricted.

Taking arbitrary values of $`\delta [0,2\pi )`$, we find that patterns A, B and C with $`m_1=0`$, in which the magnitudes of $`R_\xi `$ are suppressed by small $`s_z^2`$, can be compatible with current neutrino oscillation data. Patterns D, E and F with $`m_1=0`$ must be discarded, because their predictions for $`R_\xi `$ are too big to agree with the present experimental data. In the $`m_3=0`$ case, patterns B, C, D and F are allowed by current data; and patterns A and E are ruled out. Our numerical analyses show that pattern A with $`m_3=0`$ yields $`R_\zeta 0.5`$ and pattern E with $`m_3=0`$ yields $`\zeta >1`$, thus both patterns must be rejected. We are then left with seven phenomenologically acceptable patterns of the neutrino mass matrix with one texture zero and one vanishing eigenvalue. Subsequently we discuss each of the seven patterns in some detail.

Pattern A with $`m_1=0`$: We have $`\xi s_z^2`$ and $`R_\xi s_z^4`$, whose exact expressions are independent of the Dirac phase $`\delta `$. Considering Eqs. (2) and (6) for $`\alpha =\beta =e`$, we obtain $`\sigma =\pm \pi /2`$. The other Majorana phase $`\rho `$ is completely unrestricted. Because of $`M_{ee}=0`$ in this pattern, the effective mass of the neutrinoless double beta decay $`m_{ee}|M_{ee}|`$ vanishes. This result has no conflict with the present experimental data on $`m_{ee}`$ , namely, $`m_{ee}<0.35`$ eV at the $`90\%`$ confidence level. It is remarkable that three neutrino masses perform a clear hierarchy: $`m_1=0`$, $`m_2=\sqrt{\mathrm{\Delta }m_{\mathrm{sun}}^2}8\times 10^3`$ eV and $`m_3\sqrt{\mathrm{\Delta }m_{\mathrm{atm}}^2}5\times 10^2`$ eV.

Pattern B with $`m_1=0`$: We have $`\xi s_z`$ and $`R_\xi s_z^2`$, whose exact expressions are weakly dependent on the Dirac phase $`\delta `$. Considering Eqs. (2) and (6) for $`\alpha =e`$ and $`\beta =\mu `$, we obtain $`\sigma =(\mathrm{\Delta }_1\pm \pi )/2`$, where

$$\mathrm{tan}\mathrm{\Delta }_1=\frac{c_xc_ys_\delta }{c_xc_yc_\delta s_xs_ys_z},$$
(12)
with $`s_\delta \mathrm{sin}\delta `$ and $`c_\delta \mathrm{cos}\delta `$. The other Majorana phase $`\rho `$ is entirely unrestricted. The neutrino mass spectrum in this pattern is essentially the same as that in pattern A with $`m_1=0`$. Because $`m_{ee}m_3s_z`$ holds, the effective mass of the neutrinoless double beta decay is maximally of $`𝒪(10^2)`$ eV.

Pattern C with $`m_1=0`$: We have $`\xi s_z`$ and $`R_\xi s_z^2`$, whose exact expressions depend weakly on the Dirac phase $`\delta `$. Considering Eqs. (2) and (6) for $`\alpha =e`$ and $`\beta =\tau `$, we obtain $`\sigma =(\mathrm{\Delta }_1^{}\pm \pi )/2`$, where

$$\mathrm{tan}\mathrm{\Delta }_1^{}=\frac{c_xs_ys_\delta }{c_xs_yc_\delta +s_xc_ys_z}.$$
(13)
One can see that $`\mathrm{\Delta }_1^{}=\mathrm{\Delta }_1`$ will hold, if the replacement $`\theta _y\theta _y+\pi /2`$ is made for Eq. (12). Again, the Majorana phase $`\rho `$ is entirely unrestricted. The neutrino mass spectrum in this pattern is essentially the same as that in pattern A or pattern B with $`m_1=0`$. As a result of $`m_{ee}m_3s_z`$, the effective mass of the neutrinoless double beta decay is maximally of $`𝒪(10^2)`$ eV.

Pattern B with $`m_3=0`$: If $`\zeta `$ and $`R_\zeta `$ are expanded in powers of $`s_z`$, one can obtain $`(1\zeta )s_z`$ and $`R_\zeta s_z`$, whose magnitudes depend strongly on the Dirac phase $`\delta `$. The explicit dependence of $`R_\zeta `$ on $`\delta `$ is shown in Fig. 1, where current neutrino oscillation data have been taken into account. Considering Eqs. (2) and (9) for $`\alpha =e`$ and $`\beta =\mu `$, we obtain $`(\rho \sigma )=(\mathrm{\Delta }_2\pm \pi )/2`$, where

$$\mathrm{tan}\mathrm{\Delta }_2=\frac{s_yc_ys_zs_\delta }{s_xc_x(s_y^2s_z^2c_y^2)+(s_x^2c_x^2)s_yc_ys_zc_\delta }.$$
(14)
Because of $`m_3=0`$ and $`m_1m_2=\sqrt{\mathrm{\Delta }m_{\mathrm{atm}}^2}5\times 10^2`$ eV, the mass spectrum of three neutrinos performs an inverted hierarchy. In this case, the effective mass of the neutrinoless double beta decay $`m_{ee}`$ is expected to be of $`𝒪(10^2)`$ eV or smaller.

Pattern C with $`m_3=0`$: If $`\zeta `$ and $`R_\zeta `$ are expanded in powers of $`s_z`$, one can also obtain $`(1\zeta )s_z`$ and $`R_\zeta s_z`$, whose magnitudes depend strongly on the Dirac phase $`\delta `$. The explicit dependence of $`R_\zeta `$ on $`\delta `$ is illustrated in Fig. 1, where the present experimental data have been taken into account. Considering Eqs. (2) and (9) for $`\alpha =e`$ and $`\beta =\tau `$, we obtain $`(\rho \sigma )=(\mathrm{\Delta }_2^{}\pm \pi )/2`$, where

$$\mathrm{tan}\mathrm{\Delta }_2^{}=\frac{s_yc_ys_zs_\delta }{s_xc_x(s_y^2c_y^2s_z^2)+(s_x^2c_x^2)s_yc_ys_zc_\delta }.$$
(15)
It is easy to check that $`\mathrm{\Delta }_2^{}=\mathrm{\Delta }_2`$ will hold, if the replacement $`\theta _y\theta _y+\pi /2`$ is made. The neutrino mass spectrum in this pattern is essentially the same as that in pattern B with $`m_3=0`$. Thus we expect that $`m_{ee}`$ is of $`𝒪(10^2)`$ eV or smaller.

Pattern D with $`m_3=0`$: If $`\zeta `$ is expanded in powers of $`s_z`$, one can obtain $`\zeta c_x^2/s_x^2`$. Then the smallness and positiveness of $`R_\zeta `$ depends on the Dirac phase $`\delta `$ in an indirect way. The numerical dependence of $`R_\zeta `$ on $`\delta `$ is shown in Fig. 2, where current neutrino oscillation data have been taken into account. Considering Eqs. (2) and (9) for $`\alpha =\beta =\mu `$, we get $`(\rho \sigma )=\mathrm{\Delta }_2\pm \pi /2`$, where $`\mathrm{\Delta }_2`$ has been given in Eq. (13). We also have $`m_1m_2=\sqrt{\mathrm{\Delta }m_{\mathrm{atm}}^2}5\times 10^2`$ eV in this pattern. Thus the neutrino mass spectrum is essentially the same as that in pattern B or C with $`m_3=0`$. Again, the size of $`m_{ee}`$ is anticipated to be of $`𝒪(10^2)`$ eV or smaller.

Pattern F with $`m_3=0`$: If $`\zeta `$ is expanded in powers of $`s_z`$, one can obtain $`\zeta c_x^2/s_x^2`$. The smallness and positiveness of $`R_\zeta `$ depends on the Dirac phase $`\delta `$ in an indirect way. We illustrate the numerical dependence of $`R_\zeta `$ on $`\delta `$ in Fig. 2. Considering Eqs. (2) and (9) for $`\alpha =\beta =\tau `$, we get $`(\rho \sigma )=\mathrm{\Delta }_2^{}\pm \pi /2`$, where $`\mathrm{\Delta }_2^{}`$ has been given in Eq. (14). It is found that the neutrino mass spectrum in this pattern is essentially the same as that in pattern B, C or D with $`m_3=0`$. Thus we anticipate $`m_{ee}`$ to be of $`𝒪(10^2)`$ eV or smaller.

In summary, we have examined whether the Majorana neutrino mass matrix $`M`$ with one texture zero and one vanishing eigenvalue can interpret the experimental data on solar, atmospheric and reactor neutrino oscillations. We find that three one-zero patterns of $`M`$ with $`m_1=0`$ (normal hierarchy) and four one-zero patterns of $`M`$ with $`m_3=0`$ are compatible with current data. Thus these simple textures of $`M`$ are interesting and useful for model building, in order to reach some deeper understanding of the origin of lepton masses and flavor mixing.

It is worth mentioning that the one-zero pattern F with $`m_3=0`$ has recently been derived by He  from a modified version of the Zee model . One may examine that other one-zero textures of $`M`$ with $`m_1=0`$ or $`m_3=0`$ are possible to result from similar modifications of the Zee model. If the seesaw mechanism  is taken into account, one may also derive possible one-zero patterns of the heavy right-handed Majorana neutrino mass matrix in a phenomenological way with a few reasonable assumptions (see Ref.  for some detailed discussions).

Of course, a specific texture of the neutrino mass matrix $`M`$ may not be preserved to all orders or at any energy scales in the unspecified interactions from which neutrino masses are generated . Whether $`m_1=0`$ or $`m_3=0`$ is stable against radiative corrections should also be taken into account in a concrete model, into which the simple ans$`\ddot{\mathrm{a}}`$tze of $`M`$ under discussion can be incorporated. Between the electroweak scale ($`m_Z^{}`$) and the mass scale of the lightest right-handed Majorana neutrino ($`m_{\mathrm{R1}}^{}`$) in the seesaw mechanism , however, the running behavior of the effective left-handed neutrino mass matrix $`M`$ is rather simple at the one-loop level (see, e.g., Ref. ):

$$M(m_{\mathrm{R1}}^{})=I_g^1I_t^1\left(\begin{array}{ccc}\sqrt{I_e}& 0& 0\\ 0& \sqrt{I_\mu }& 0\\ 0& 0& \sqrt{I_\tau }\end{array}\right)^1M(m_Z^{})\left(\begin{array}{ccc}\sqrt{I_e}& 0& 0\\ 0& \sqrt{I_\mu }& 0\\ 0& 0& \sqrt{I_\tau }\end{array}\right)^1,$$
(16)
where $`I_g`$, $`I_t`$ and $`I_\alpha `$ (for $`\alpha =e,\mu ,\tau `$) stand respectively for the evolution functions of the gauge couplings, the top-quark Yukawa coupling and the charged-lepton Yukawa coupling eigenvalues. It becomes obvious that the vanishing eigenvalue and texture zeros of $`M`$ at $`m_Z^{}`$ remain to hold at $`m_{\mathrm{R1}}^{}`$, at least at the one-loop level under consideration . This interesting observation implies that our classification of $`M`$ with texture zeros or vanishing eigenvalues at low energy scales does make some sense for model building at a much higher energy scale.

I like to thank X.J. Bi for useful discussions. This work was supported in part by the National Natural Science Foundation of China.


text:
# Radial velocities of “slow movers” – call for observations

## 1 Introduction

During our research on the long-period comets origins we simulate the mechanisms of production of the observable comets directly form the Oort(1950) cometary cloud by means of stellar and galactic perturbations. A natural consequence of that study would be a search for stellar perturbers which passed near the Sun during the last several million years or will have such a chance in the future. Such passages can induce strong asymmetries in the observed population of the long-period comets, as it was showed by Dybczyński (2002a, b, c). Additionally, such passages should be taken into account when studying the past motion of observed long-period comets as it was demonstrated in Dybczyński (2001).

In our first attempt to complete a list of potential stellar perturbers (Dybczyński & Kankiewicz 1999) we found, that the majority of stars with measured parallax (most of them are from the Hipparcos catalogue, ESA 1997) do not have their radial velocities measured what excludes them from further consideration. Similar conclusions were reported by García-Sánchez et al. (1999, 2001): they could attribute radial velocity values to less than 50% of the selected candidate stars. The lack of radial velocities in such cases is most often a consequence of the selection criteria used by observers. When searching for the promising candidate for radial velocity measurements they prefer stars which high proper motion suggests high space velocity. Thus “fast movers” (in the sense of high proper motion) are highly favored. In our research we are interested in the stars moving almost directly towards (or outwards) the Sun what implies very small proper motion. We will call these stars “slow movers”.

## 2 Missing radial velocities

In our study on the stellar perturbations of the Oort cloud comets we decided to use the most precise subset of the Hipparcos catalogue (ESA 1997), the ARIHIP catalogue (Wielen et al. 2001b), containing the results of the Hipparcos mission re-reduced for better proper motion calculations with the aid of the precise ground based astrometric catalogues*.* ARIHIP is available on the Internet at *http://www.ari.uni-heidelberg.de/arihip/*. For our purpose we used the “long term” solutions included in ARIHIP because they seem to be the best available approximations of the long term motion of these stars. For a more detailed explanation see Wielen et al. (2001a).

Our procedure of selecting candidate stars is very similar to that used by García-Sánchez et al. (1999, 2001). For all stars with the reliable parallax in ARIHIP we calculated the miss-distance between a star and the Sun basing on the straight line motion model. As it was shown this simple model is quite sufficient as the first selection tool. In our first list of candidates (hereafter LC) we included all stars with the calculated miss-distance less than 5 pc assuming for all stars $`v_r=100`$ km s-1. This list consisted of 2123 stars. However, there are radial velocities included in ARIHIP for some stars so in the second step we used these values and again adopted $`v_r=100`$ km s-1 for the rest of the stars. We found ARIHIP velocities for 824 stars and as a consequence 738 of them have been rejected because of the recalculated miss-distance greater than 5 pc. Only 86 were adopted as the Oort cloud perturbers (hereafter LP).

Next we searched other sources for the missing radial velocities. The most helpful were the ADS (*http://adswww.harvard.edu*) and CDS (*http://cdsweb.u-strasbg.fr*) services as well as the Bibliographic Catalogue of Stellar Radial Velocities available on the Internet (*http://www.casleo.gov.ar/catalogo/catalogoin.htm*). The latter is the 2002 update of the previously published version (Malaroda et al. 2001). As a result we found additional radial velocity values for 199 stars, 181 of which were subsequently excluded from the list: after recalculations their miss-distance become greater than 5 pc. The remaining 18 stars were added to the list of stellar perturbers (LP).

Finally we obtained two lists: LP, with the confirmed stellar perturbers among the ARIHIP stars with miss-distance less than 5 pc and the LC list of potential perturbers which radial velocities are still unknown. A detailed description of our LC list and its “most wanted” subset are presented in Sect.4 while the influence of the stars from the LP list on the Oort cloud comets motion will be discussed in detail in the separate paper, in preparation.

## 3 Observing plans for slow movers

As the number of slow movers is quite large, reliable measurements of their radial velocities will deserve a significant amount of observing time. The accuracy needed to filter out most of the non-significant objects (with relatively low $`v_r`$ ) is not particularly high and thus it does not put a special demand on the instrumentation to be used. When found, slow movers with a high radial velocity can be studied further with a better spectrograph, yielding absolute accuracy on the level of $`0.1`$ km s-1.

There were several projects aiming at measurements of radial velocities of Hipparcos stars. Up till now they obtained radial velocities of several thousands objects (see for example Grenier et al. (1999); Nidever et al. (2002)), but their selection criteria were different from ours so it may take a lot of time before all slow movers are measured.

To speed up the process of building a database of slow movers radial velocities we would like to get observers interested in this group of star. To avoid unnecessary duplication of efforts we have set up the LC list of candidates for perturbers which is available on the web. The immediate aim of this would to be to mark all stars that have already been measured even if the results are not publicly available yet. Later, as new radial velocities are being published, the list will be gradually changed into a database, which can easily be used for different research projects.

We are also planning our own observations with a newly build fiber-fed echelle spectrograph. While the main aim of this instrument will be radial velocity curves of binary stars, slow movers will also be measured during its free time.

Our new instrument allows us to measure radial velocities of late spectral type stars down to 8-9 mag, with 30-60 min exposures. The location of the telescope ($`\varphi =52^{}`$) sets additional constraints on the available targets. Going through our slow movers list and selecting stars with $`\delta >8^{}`$ and spectral type later then F0 we obtain 18 stars with visual magnitude $`V<7`$ mag, 94 stars with $`7<V<8`$ mag and 158 stars with $`8<V<9`$ mag. We plan to observe them starting from the brightest end and moving on to the dimmer ones. As we progress, our list available at the web-page will be updated.

## 4 Proposed targets

As it was described in Sec.2 after using all radial velocity values available to us we shortened the LC list to 1100 stars. All of them have the Sun miss-distance less than 5 pc assuming $`v_r=100`$ km s-1 . In Figs. 1 and 2 we present distributions of the visual magnitude and spectral type among proposed targets. As it can be easily noted, these stars are relatively bright and the majority of them are of late spectral types. The proposed targets list (LC) is available on the web at the address: *http://www.astro.amu.edu.pl/Science/SlowMovers*/.

The main list is named *targets.lst.* The first column gives the Hipparcos number of the star (strictly ascending), next is the Tycho-2 corresponding number and the third column shows the calculated Sun miss-distance ($`DD`$) in parsecs. The lower is the $`DD`$ value the higher is our interest in this particular object, ie. the probability of confirmation that the star is actually the Oort cloud perturber.

For the user convenience in the next four columns we added the visual magnitude of the star, its position in the sky (J2000 Right Ascension and Declination) and its spectral type copied from the Tycho-2 Spectral Type Catalog (Wright et al. 2003), recently made available at CDS. \>From our 1100 stars 96 are not included in this catalogue but with the help of SIMBAD database we decreased the number of missing spectra down to 48 stars.

The calculated Sun miss-distance of each star on the LC list is computed for the assumed $`v_r=100`$ km s-1 . To obtain the “most wanted” group of targets we repeated all the calculations, now assuming $`v_r=20`$ km s-1 and limiting the Sun miss-distance $`DD`$ to be smaller than 3 pc. This gave us a list of 13 stars of particularly high probability of confirmation as an Oort cloud perturber. They are listed in Tab.1.

## 5 Conclusions

As it was mentioned above, the main purpose of this paper is to raise the interest in measuring radial velocities of the proposed 1100 low proper motion stars. Due to the current location of our instrument we cannot observe stars having declinations below $`8^{}`$. There are 572 such stars in our list. We are also restricted to the brighter end of the list but still the number of stars is quite large: it takes time to determine the absolute radial velocity of the star, especially when several subsequently obtained values are not consistent with each other.

So, here is the “call for cooperation”. We offer some kind of coordination by marking as “done” already measured stars or marking as “planned” other targets upon the information obtained from observers. First such marks were added to our list during the preparation of this manuscript. Following the private communications from Dave Latham we marked stars that have several unpublished radial velocity measurements made with CfA Digital Speedometers (238 such stars) and stars included in the CfA list of 10000 main-sequence stars closer than 100 pc and that have at least one new CfA radial velocity (unpublished) measurement (142 such stars).

Acknowledgements.This research has made an extensive use of the SIMBAD database, operated at CDS, Strasbourg, France and has made use of NASA’s Astrophysics Data System. This work was partially supported by the KBN Grants no. 2P03D01324 (PAD) and 5P03D00220 (TK). This manuscript was prepared with the LyX front-end under the Debian GNU Linux.


text:
# Perturbation of coupling matrices and its effect on the synchronizability in arrays of coupled chaotic systems.

## 1 Introduction

The object of interest in this paper is a coupled array of $`n`$ identical systems with linear coupling:

$$\dot{x}=\left(\begin{array}{c}f(x_1,t)\\ \mathrm{}\\ f(x_n,t)\end{array}\right)(GD)x$$
(1)
where $`x=(x_1,\mathrm{},x_n)`$. The coupling matrix $`G`$ describes the coupling topology of the array, i.e. $`G_{ij}0`$ if there is a coupling between system $`x_i`$ and system $`x_j`$. If $`G_{ij}<0`$ or $`G_{ij}>0`$, we call such a coupling element cooperative or competitive respectively. In this paper we assume that $`G`$ is a symmetric zero row sums matrix whose off-diagonal elements are all nonpositive. This implies that all coupling is cooperative and all eigenvalues of $`G`$ are real. The matrix $`D`$ describes the coupling between a pair of systems. In  it was shown that for a suitable $`D`$ (such as $`D=I`$), the synchronization criterion depends solely on the second smallest eigenvalue of $`G`$, which we denote as $`\lambda _2(G)`$. In particular, the larger $`\lambda _2(G)`$ is, the easier it is to synchronize the array. Furthermore, when the single uncoupled system $`\dot{y}=f(y,t)`$ is a chaotic system, $`\lambda _2`$ must necessarily be positive in order for the array to synchronize. Therefore we will equate $`\lambda _2`$ with synchronizability in the remainder of this paper. In  several coupling topologies were considered. In particular, for nearest neighbor coupling, $`\lambda _20`$ as $`n\mathrm{}`$ and thus such an array of chaotic systems will not synchronize for large $`n`$.

In  it was shown using wavelet analysis how a perturbed coupling matrix $`\stackrel{~}{G}=G+\mathrm{\Delta }G`$ can be constructed such that $`\lambda _2(\stackrel{~}{G})`$ is much larger than $`\lambda _2(G)`$ and thus the synchronizability is greatly improved. The perturbation $`\mathrm{\Delta }G`$ can be thought of as the additional coupling added to the coupling array to affect its synchronization properties. In , $`\stackrel{~}{G}`$ is obtained by modifying a small fraction (the $`LL_1`$ subspace) of the wavelet coefficients of $`G`$. However, all of the matrix elements of $`G`$ are changed in the process, i.e. $`\mathrm{\Delta }G`$ is not sparse. In particular, there are no zero elements in $`\stackrel{~}{G}`$, the perturbed array is fully connected and every chaotic system is connected to every other chaotic system. To illustrate the increase of $`\lambda _2`$, in  the $`LL_1`$ coefficients are multiplied by a factor of $`1101`$. This results in a full matrix $`\stackrel{~}{G}`$ whose off-diagonal elements are up to 8 times larger than the original matrix $`G`$. With such drastic changes to $`G`$, it is not surprising that $`\lambda _2`$ will also change dramatically.

Note that Eq. (1) has a minus sign in front of the term $`(GD)x`$ as opposed to a plus sign in . This means that the coupling matrix here is the negative of the coupling matrix in . In other words, whereas  talks about the second largest eigenvalue, we consider the second smallest eigenvalue. Writing the state equation as in Eq. (1) has the benefit that $`\lambda _2`$ is positive (rather than negative as in ). Furthermore, in this way, $`G`$ and $`\lambda _2`$ correspond to the graph-theoretical notions of Laplacian matrix and algebraic connectivity respectively . This is notationally less cumbersome later on when we use results from the theory of random graphs.

The purpose of this paper is to attempt to answer the question of what type of matrices $`\mathrm{\Delta }G`$ will make $`\lambda _2(\stackrel{~}{G})`$ much larger than $`\lambda _2(G)`$, i.e. greatly improve the synchronizability while keeping $`\mathrm{\Delta }G`$ “small”. This notion of smallness could be the sparseness of $`\mathrm{\Delta }G`$ as this refers to the number of additional coupling added to the original array. It could also be the size of the off-diagonal elements of $`\mathrm{\Delta }G`$ as they relate to the coupling strengths between the chaotic systems in the array. We will show that there exists sparse matrices $`\mathrm{\Delta }G`$ with small elements that will change $`\lambda _2(G)`$ significantly. We also show that uniform coupling will increase the synchronizability the most. On the other hand, we show that local coupling, i.e. $`\stackrel{~}{G}`$ corresponding to a graph with only local connections, will not synchronize when the number of systems is large.

We also show how this analysis allows us to prove and explain the change in synchronizability observed in a coupled array where random coupling elements are added to a nearest neighbor coupled array.

Finally we will classify several classes of coupling matrices considered in the literature according to a synchronizability measure, and show that the more random and less ordered the coupling graph is, the higher the synchronizability.

## 2 Super-additivity of $`\lambda _2`$

###### Definition 1

$`W`$ is the class of symmetric matrices that have zero row sums and nonpositive off-diagonal elements.

It can be shown that matrices in $`W`$ are positive semi-definite and that $`0`$ is a simple eigenvalue if and only if the matrix is irreducible . We assume that the coupling matrices $`G`$ and $`\stackrel{~}{G}`$ are in $`W`$.

###### Lemma 1

If $`A`$ and $`B`$ are in $`W`$, then $`\lambda _2(A+B)\lambda _2(A)+\lambda _2(B)`$.

Proof: see . $`\mathrm{}`$

This lemma shows that $`\lambda _2(\stackrel{~}{G})\lambda _2(G)+\lambda _2(\mathrm{\Delta }G)`$ when $`\mathrm{\Delta }G`$ is also a matrix in $`W`$. Therefore if we can find $`\mathrm{\Delta }G`$ such that $`\lambda _2(\mathrm{\Delta }G)`$ is large, then we can increase $`\lambda _2(\stackrel{~}{G})`$ and thus the synchronizability.

A corollary to Lemma 1 is that if $`A`$, $`B`$ are matrices in $`W`$ and every off-diagonal element of $`A`$ is smaller in absolute value than $`B`$, then $`\lambda _2(A)\lambda _2(B)`$. This means that adding additional cooperative coupling elements ($`G_{ij}<0`$) can only increase $`\lambda _2`$ and improve the synchronizability.

Let us consider three reasonable constraints on the possible coupling matrix $`G`$ that we can choose:

1. $`GW`$.
2. the number of nonzero coupling elements is a percentage of the total number of coupling elements, i.e. only a percentage of the off-diagonal elements of $`G`$ are nonzero.
3. the absolute value of each off-diagonal element is less than some fixed number $`B>0`$.

Under these constraints, the above discussion shows that the matrix $`G`$ with the largest $`\lambda _2`$ must necessarily be matrices in $`W`$ where the off-diagonal elements are either $`0`$ or $`B`$. We call such a coupling configuration a uniform coupling as the coupling between any pair of systems is either zero or the same111This notion of uniform coupling is different than the notion used by some authors to denote fully connected coupling.. Thus given a bound on the size of each coupling element, the uniform coupling is the best coupling configuration in terms of synchronizability.

This motivates us to study the class of matrices in $`W`$ where the off-diagonal elements is either 0 or -1. Such a matrix is the Laplacian matrix of the underlying graph: vertex $`i`$ is connected to vertex $`j`$ if and only if $`G_{ij}0`$ and $`\lambda _2`$ corresponds to the algebraic connectivity of the graph . This allows us to use results in graph theory to study $`\lambda _2(\stackrel{~}{G})`$. In particular, in the rest of the paper, we will mainly consider coupling matrices which is a positive multiple of the Laplacian matrix of the underlying graph and the goal is to find appropriate coupling graphs whose algebraic connectivity are large.

## 3 Local coupling

Consider the following model of local coupling : let the vertices of the coupling graph be the points in the hypercube $`[0,r]^d`$ with integer coordinates, i.e. the vertices are the lattice points inside the hypercube $`[0,r]^d`$ generated by the standard basis vectors. Let the systems in the array be located on these vertices. We say this array is locally connected if there exists a fixed $`\delta `$ such that all the systems which are connected to each system $`x`$ lies in the $`\delta `$-neighborhood of $`x`$. This notion of local coupling includes examples such as nearest neighbor coupling and grid graphs, but also more general local coupling configurations where the vertex degree is bounded but not necessarily constant.

In  it was shown that:

###### Lemma 2

For fixed $`d`$ and $`\delta `$, a locally connected array will have $`\lambda _20`$ as $`n\mathrm{}`$ (i.e. $`r\mathrm{}`$). This is true even if the hypercube is considered as a hypertorus.

This means that locally connected arrays will have low synchronizability for large $`n`$. In particular, locally connected arrays of chaotic systems will not synchronize for large $`n`$. This means that $`\mathrm{\Delta }G`$ should have nonlocal coupling in order to improve the synchronizability of the array. Next we show that random coupling is some sense the best candidate for such nonlocal coupling.

## 4 Random coupling

How much coupling do we need to add in order to make $`\lambda _2`$ large while keeping the size of the off-diagonal element small? In  it was shown that:

###### Lemma 3

For a non-complete graph $`\lambda _2`$ is less than or equal to the minimum vertex degree222For a complete graph, $`\lambda _2(K_n)=n`$..

In particular, for a $`k`$-regular (non-complete) graph, $`\lambda _2k`$. Therefore, if we want $`\lambda _2`$ to be large we need each vertex to connect to a large number of other vertices, i.e. $`\mathrm{\Delta }G`$ cannot be too sparse. Regular random graphs appear to be asymptotically optimal in this respect by Lemma 3 and the following lemma.

###### Lemma 4

There exists a constant $`0<c<1`$ such that almost every $`k`$-regular random graph satisfies $`\lambda _2>ck`$ as $`n\mathrm{}`$ for large enough $`k`$. For even $`k`$, $`c`$ can be chosen close to $`1`$.

Proof: See . $`\mathrm{}`$

For small $`k`$, it can be shown that for $`k3`$, there is a constant $`c_k>0`$ such that almost every $`k`$-regular random graph satisfies $`\lambda _2c_k`$ as $`n\mathrm{}`$. In particular, since $`\lambda _2k\sqrt{k^2i^2}`$  where $`i`$ is the isoperimetric number of the graph, and using the results in , a graph of $`c_k`$ for $`k`$ between $`3`$ and $`15`$ is shown in Figure 1.

This allows us to prove the following theorem which says that there is an $`n`$ by $`n`$ coupling matrix $`\mathrm{\Delta }G`$ with an arbitrarily small percentage of nonzero off-diagonal elements, each of which is arbitrarily small, and which has arbitrarily large $`\lambda _2`$ when $`n`$ is large enough.

###### Theorem 1

Given $`0<p1`$, $`r>0`$, and $`s`$, for all large enough $`n`$, there exists a matrix $`\mathrm{\Delta }G`$ in $`W`$ such that

* $`\lambda _2>s`$,
* the fraction of nonzero elements of $`\mathrm{\Delta }G`$ is $`p`$,
* each off-diagonal element of $`\mathrm{\Delta }G`$ is less than $`r`$.

Proof: Choose $`k`$ an integer larger than $`\frac{s}{cr}`$. By Lemma 4 for large enough $`n`$ and $`k`$, almost every random $`k`$-regular graph will have $`\lambda _2>\frac{s}{r}`$. Choose large $`n`$ such that $`n>\frac{k}{p}`$. Let $`A`$ be the Laplacian matrix of such a graph. Then $`\mathrm{\Delta }G=rA`$ will suffice. $`\mathrm{}`$

So far, the random graph corresponding to $`\mathrm{\Delta }G`$ is required to be a regular graph. Next we show that this requirement can be relaxed. In  a network model is presented where a nearest neighbor configuration is modified by adding a fraction $`0p1`$ of random edges out of all possible edges333The authors of  erroneously associated this model with the small world model in . In  the number of added edges is a fraction $`p`$ of the existing edges in the nearest neighbor coupling. Therefore the number of edges in  is much less. In particular, for $`p=1`$, the graph in  is fully connected, whereas in  the number of edges is only twice the number of edges of the nearest coupling graph. In particular for the one-dimensional nearest neighbor coupling, by Lemma 3 the algebraic connectivity of the model in  is at most 2, whereas it is unbounded in  as $`n\mathrm{}`$.. In  it was shown via numerical experiments that $`\lambda _2(\stackrel{~}{G})`$ increases in a near linear fashion with respect to $`p`$ and to $`n`$. We prove that this is true in general. Let $`G`$ be a locally connected graph and $`\mathrm{\Delta }G`$ is constructed as in  by picking randomly a fraction $`p`$ of all possible edges.

First note that by Lemmas 1 and 3, $`\lambda _2(\stackrel{~}{G})pn+q`$ where $`q`$ is the average number of connections of $`G`$. In  it was shown that $`\lambda _2(\mathrm{\Delta }G)`$ is equal to $`pn+o(n^{\frac{1}{2}+ϵ})`$ in probability for any $`ϵ>0`$. Combining this with Lemmas 1 and 2, it’s clear that $`\lambda _2(\stackrel{~}{G})`$ increases as $`o(pn)`$ in probability as $`n\mathrm{}`$, i.e. linear in $`p`$ and in $`n`$.

Thus in general, given a locally connected array, the addition of a fraction of random edges out of all possible edges will increase the synchronizability by an arbitrarily large amount as the number of vertices increases. Note that the synchronizability of the random edges dominates that of the locally connected graph.

### 4.1 Example

To compare the results in this section with the results in , consider the example in  that was discussed in Section 1 where an array of 64 oscillators with nearest neighbour coupling is perturbed by multiplying the $`LL_1`$ wavelet coefficients of the coupling matrix $`G`$ by a factor of $`1101`$. This increases $`\lambda _2`$ from $`0.0096`$ to $`0.1522`$. However, all entries of the 64 by 64 matrix $`G`$ are changed, with the off-diagonal entries of the perturbation matrix $`\mathrm{\Delta }G`$ ranging from $`7.6`$ to $`9.7`$. On the other hand, suppose we want to change only $`5`$ percent of the off-diagonal entries of $`G`$. By choosing a random graph with $`100`$ edges and corresponding Laplacian matrix $`L`$, we construct a perturbation matrix $`\mathrm{\Delta }G=0.0665L`$ which also increases $`\lambda _2`$ from $`0.0096`$ to $`0.1522`$. However, in this case, only $`5`$ percent of the off-diagonal entries of $`G`$ are perturbed, and the off-diagonal entries of the perturbation matrix $`\mathrm{\Delta }G`$ is either $`0.0665`$ or $`0`$. Thus the improvement in synchronization is the same as in , but the perturbation matrix is much smaller, both in the sense of the number of nonzero elements and in the sense of the magnitude of the nonzero elements.

## 5 Comparison of various classes of coupling graphs

In this section we compare several classes of coupling graphs that have been proposed in the literature for modeling complex networks, both man-made (the internet or social networks) and occuring in nature (neural networks). We can categorize their synchronizability by studying the synchronizability measure $`s(m)=lim_n\mathrm{}\frac{\lambda _2}{m}`$ when it exists, where $`m`$ is the average number of connections defined as 2 times the number of edges divided by the number of vertices, i.e. $`m`$ is the number of nonzero off-diagonal elements of the Laplacian matrix divided by $`n`$. For $`r`$-regular graphs, $`m=r`$. By Lemma 3, $`0s(m)1`$. We will only consider the cases $`m3`$.

### 5.1 Small world networks

In  a small-world network model is introduced where starting from a locally connected graph with $`m`$ edges, $`pm`$ edges ($`0<p1`$) are added randomly. Without the random edges $`\lambda _20`$ as $`n\mathrm{}`$ by Lemma 2. In the cases considered in , $`m`$ grows linearly with $`n`$, i.e. $`m=kn`$. If we modified this model such that the added random edges form a regular random graph444This requires that $`2kp`$ is an integer., then the discussion in Section 4 shows that $`\lambda _2`$ is bounded away from zero as $`n\mathrm{}`$ when $`kp\frac{3}{2}`$. In particular, Figure 1 shows that $`s(m)0.0018`$. Thus by adding the random edges, there is a coupling strength factor $`ϵ>0`$ such that when $`D`$ is replaced by $`ϵD`$, the coupled array will synchronize independent of $`n`$.

### 5.2 Scale-free random networks

In  a scale-free random network model of the internet is presented where the graph is built up as follows: Starting from a initial set of $`k_0`$ (usually set at $`k_0=k`$) vertices, at each iteration a new vertex is added with $`k`$ new edges connecting this vertex to $`k`$ existing vertices. Thus the average connections is $`m=2k`$ as $`n\mathrm{}`$. The probability that the new vertex is connected to vertex $`i`$ is proportional to the vertex degree of vertex $`i`$. This results in a graph whose vertex degree distribution follows a power law. Since there is always at least one vertex with degree $`k`$, $`\lambda _2k`$ and $`s(m)\frac{1}{2}`$ by Lemma 3. Computer simulations555We differ slightly from  in our construction of the scale-free networks. In  the initial configuration is $`k_0`$ vertices with no edges. After many iterations, this could still result in some of these initial $`k_0`$ vertices having degree less than $`k`$, even though every newly added vertex has degree at least $`k`$. Since even one vertex with degree 1 results in $`\lambda _2`$ to drop below $`1`$, we would like to have every vertex of the graph to have degree at least $`k`$. Therefore we will start with the initial $`k_0`$ vertices fully connected to each other when $`k_0=k`$. For the case $`k_0>k`$, we require that the initial $`k_0`$ vertices are connected with each other with vertex degrees at least $`k`$. This will guarantee that the vertex degree of every vertex is at least $`k`$. show that for a fixed $`m`$, $`\lambda _2`$ of such graphs appear to increase monotonically as a function of $`n`$. If this is the case, then since $`\lambda _2`$ is bounded by $`\frac{m}{2}`$, this means that $`\lambda _2`$ converges to a constant value $`\lambda _2^{\mathrm{}}(m)`$ as $`n\mathrm{}`$. We computed via simulations the value of $`s(m)=\frac{\lambda _2^{\mathrm{}}(m)}{m}`$ for various $`m`$ as shown in Figure 2. We see that $`s(m)`$ is increasing and converging to close to the upper bound $`\frac{1}{2}`$. As the scale-free network is not as random as a random network, its synchronizability measure is also less.

For locally connected graphs, $`s(m)=0`$ for all $`m`$ by Lemma 2. For both small-world networks and the model in  $`s(m)1`$ for large $`m`$. The main difference is that in small-world networks, $`m`$ (and thus also $`\lambda _2`$) is small and bounded, whereas in , $`m=pn`$ and therefore $`\lambda _2`$ is unbounded. These results are summarized in Table 1.

Thus we can order these classes of matrices as follows in terms of increasing synchronizability: locally connected graphs, scale-free networks, small-world networks, and the model in . We see that the more disorder and randomness the graph exhibits, the higher the synchronizability.

## 6 Conclusions

We study the conditions under which perturbation of the coupling in a coupled array results in improved synchronizability by reducing it to the analysis of the second smallest eigenvalue of the coupling matrix. We show that local coupling and random coupling form two extremes in their effect on synchronizability. In particular, local coupling will not improve synchronizability in the limit and random coupling can change the synchronizability using arbitrarily small percentage of coupling elements whose strengths are arbitrarily small. As graphs become less ordered and more random, we see that the synchronizability increases. Furthermore, for graphs which is constructed by adding a random component to a structured component, the super-additivity of $`\lambda _2`$ allows us to give a lower bound for the synchronizability of graphs by adding the synchronizability of each of the components. In particular, when the structured component is a locally connected graph, the synchronizability of the random component dominates the synchronizability of the locally connected component.