text:
# On arrangements of real roots of a real polynomial and its derivatives

## Abstract

We prove that all arrangements (consistent with the Rolle theorem and some other natural restrictions) of the real roots of a real polynomial and of its $`s`$-th derivative are realizable by real polynomials.
Key words: arrangement of roots
AMS classification: 12D10

In the present paper we consider a real polynomial of one real variable $`P(x,a)=x^n+a_1x^{n2}+\mathrm{}+a_{n1}`$. We are interested in the question what arrangements between the real roots of $`P`$ and $`P^{(s)}`$ are possible ($`1sn1`$). To define an arrangement means to write down the roots of $`P`$ and $`P^{(s)}`$ in a chain in which every two consecutive roots are connected either by an equality or by an inequality $`<`$. The arrangement $`\alpha `$ is said to belong to the closure of the arrangement $`\beta `$ if it is obtained from $`\beta `$ by replacing some inequalities by equalitites. The results are the first step towards the study of real discriminant sets $`\{a𝐑^{n1}|\mathrm{Res}(P,P^{(s)})=0\}`$.

In an earlier paper \[KoSh\] it is shown that if $`P`$ is hyperbolic, i.e. with $`n`$ real roots, then the standard Rolle restrictions are necessary and sufficient conditions for a root arrangement to be realizable (see Theorems 2 and 4.4 in \[KoSh\]). Namely, denote by $`x_1\mathrm{}x_n`$ the roots of $`P`$ and by $`\xi _1\mathrm{}\xi _{ns}`$ the ones of $`P^{(s)}`$ (which is also hyperbolic). Then one has

$$x_l\xi _lx_{l+s}$$
(1)
for $`l=1,\mathrm{},ns`$ and every arrangement of the roots of $`P`$ and $`P^{(s)}`$ which is consistent with (1) is realizable. One presumes also that the following conditions hold:

A) If a root of $`P`$ of multiplicity $`d>s`$ coincides with a root of $`P^{(s)}`$ of multiplicity $`g`$, then $`g=ds`$ (self-evident).

B) If a root $`\xi `$ of $`P^{(s)}`$ coincides with a root of $`P`$ of multiplicity $`\kappa s`$, then $`\xi `$ is a simple root of $`P^{(s)}`$ (see \[KoSh\], Lemma 4.2) and one has $`\kappa s1`$.

C) If $`x_l=\xi _l`$ or $`x_{l+s}=\xi _l`$, then $`x_l=x_{l+1}=\mathrm{}=x_{l+s}=\xi _l`$ (self-evident for $`s=1`$ and easy to prove by induction on $`s`$ for $`s>1`$).

###### Example 1

If $`n=2`$, $`s=1`$, then there are two possible arrangements (i.e. consistent with (1), A) B) and C)) : $`x_1<\xi _1<x_2`$ and $`x_1=\xi _1=x_2`$. They are both realizable by hyperbolic polynomials.

In the present paper we treat the case when $`P`$ is arbitrary (not necessarily hyperbolic). (Notice that $`P^{(s)}`$ can be hyperbolic even if $`P`$ is not.)

###### Definition 2

Suppose that $`P`$ has $`m`$ conjugate couples of complex roots and $`n2m`$ real roots. Then a priori $`P^{(s)}`$ has at least $`n2ms`$ real roots counted with the multiplicities. Indeed, a real root of $`P^{(i)}`$ of multiplicity $`l1`$ is a root of $`P^{(i+1)}`$ of multiplicity $`l1`$ and between every two real roots of $`P^{(i)}`$ there is a root of $`P^{(i+1)}`$. Iterating this rule $`s`$ times one obtains the existence of $`n2ms`$ real roots of $`P^{(s)}`$ (we call them Rolle roots) which together with the real roots of $`P`$ satisfy conditions (1), A) and B). A Rolle root is multiple only if it coincides with a root of $`P`$ of multiplicity $`>s`$. Eventually, $`P^{(s)}`$ can have $`2m`$ other (non-Rolle) real roots counted with the multiplicities some (or all) of which can coincide with Rolle ones. Which real roots of $`P^{(s)}`$ should be chosen as Rolle and which as non-Rolle ones is not always uniquely defined and when it is not we assume that a choice is made.

###### Example 3

The polynomial $`x^6x^2=x^2(x^21)(x^2+1)`$ has real roots $`x_1=1`$, $`x_2=x_3=0`$, $`x_4=1`$ (and complex roots $`\pm i`$). One has $`P^{}=6x^52x=2x(\sqrt{3}x^21)(\sqrt{3}x^2+1)`$, i.e. $`P^{}`$ has three Rolle roots (and no non-Rolle ones) – $`0`$ and $`\pm 1/3^{1/4}`$ where $`0`$ is a common root for $`P`$ and $`P^{}`$, see A). It has also two complex roots $`\pm i/3^{1/4}`$. One has $`P^{\prime \prime }=30x^42`$, i.e. $`P^{\prime \prime }`$ has two Rolle roots $`\pm 1/15^{1/4}`$, no non-Rolle ones and two complex roots $`\pm i/15^{1/4}`$. One has $`P^{\prime \prime \prime }=120x^3`$, i. e. $`P^{\prime \prime \prime }`$ has a triple real root at $`0`$ and no complex roots. One copy of this real root should be considered as a Rolle one and the other two as non-Rolle ones.

###### Proposition 4

Suppose that a real root of $`P`$ of multiplicity $`d`$ coincides with a real root of $`P^{(s)}`$ of multiplicity $`g`$. Then

1) if $`d>s`$, then one has $`g=ds`$; in this case this is a Rolle root of $`P^{(s)}`$ of multiplicity $`ds`$;

2) if $`0ds`$, then one has $`g2m+1`$ (and if $`g1`$, then $`d<s`$).

Observe that in the above example one has $`m=1`$ and for $`s=3`$ the estimation $`2m+1`$ is attained by the multiplicity of $`0`$ as a root of $`P^{\prime \prime \prime }`$. The proposition generalizes conditions A) and B) in the case of arbitrary $`m`$.

Proof:

Part 1) is self-evident. Prove part 2). If the root is non-Rolle and does not coincide with a Rolle one, then its multiplicity is $`2m`$. If the root is Rolle and does not coincide with a non-Rolle one, then either it coincides with a root of $`P`$ of multiplicity $`>s`$ and we are in case 1) or it is a simple root. Finally, if the root is Rolle and coincides with a non-Rolle one, then the Rolle root must be simple (otherwise there will be a contradiction with part 1)) and the sum of their two multiplicities is $`2m+1`$. $`\mathrm{}`$

###### Definition 5

An arrangement of the real roots of $`P`$ and $`P^{(s)}`$ is called a priori admissible if there exist $`n2ms`$ Rolle roots of $`P^{(s)}`$ in the sense of Definition 2 and if conditions 1) and 2) of Proposition 4 hold.

###### Theorem 6

All a priori admissible root arrangements are realizable by real polynomials of degree $`n`$.

Proof:

$`1^0`$. We explain first in $`1^0`$$`7^0`$ why all a priori admissible arrangements in which the derivative $`P^{(s)}`$ is hyperbolic and which are the least generic are realizable. ”Least generic” means that all non-Rolle roots of $`P^{(s)}`$ coincide with Rolle ones or with roots of $`P`$ . The general case is treated in $`8^0`$$`11^0`$.

To realize an a priori admissible arrangement with $`P^{(s)}`$ hyperbolic and with the necessary multiplicities of the real roots of $`P`$ consider the family of polynomials

$$P(x,w,g,t)=\underset{j=1}{\overset{q}{}}(xw_j)^{m_j}\underset{j=1}{\overset{m}{}}((xg_j)^2+t_j^2)$$
(2)
where $`w_j`$, $`j=1,\mathrm{},q`$, are the real roots of $`P`$, of multiplicities $`m_j`$ ($`w_0=0w_1\mathrm{}w_q1=w_{q+1}`$), and $`g_j\pm it_j`$ are its complex roots (not necessarily distinct), $`t_j0`$, $`0g_j1`$. We allow here equalities between the roots $`w_j`$ for convenience; it will be shown that the necessary arrangement is realizable for roots with strict inequalities between them.

Denote by $`\xi _1\mathrm{}\xi _{ns}`$ the real parts of the roots of $`P^{(s)}`$ ($`n2ms`$ of them are just Rolle roots) and by $`\theta _1\mathrm{}\theta _m`$ the biggest nonnegative imaginary parts of the roots of $`P^{(s)}`$ (recall that for a least generic arrangement one has $`\theta _j=0`$). Set $`\xi _0=0`$, $`\xi _{ns+1}=1`$. (Notice that $`P^{(s)}`$ has not more conjugate couples of complex roots than $`P`$, i.e. not more than $`m`$.) The functions $`\xi _i`$, $`\theta _j`$ are continuous in $`(w,g,t)`$.

$`2^0`$. Suppose that for the desired arrangement of the real roots of $`P`$ and $`P^{(s)}`$ the Rolle and non-Rolle roots of $`P^{(s)}`$ are fixed. Denote the non-Rolle roots by $`u_1\mathrm{}u_{2m}`$. Impose additional requirements upon the numbers $`g_j`$ as follows: if the non-Rolle roots with odd indices $`u_{2p1},u_{2p+1},\mathrm{},u_{2p+2p^{}1}`$ belong to the interval $`[w_j,w_{j+1})`$, $`j<q`$, or to $`[w_q,w_{q+1}]`$, then we require that $`w_jg_p\mathrm{}g_{p+p^{}}w_{j+1}`$. Define the variables $`h_1\mathrm{}h_{q+m}`$ as the union of the variables $`w_j`$ ($`j=1,\mathrm{},q`$) and $`g_i`$ ($`i=1,\mathrm{},m`$) with the order defined above. Hence, they belong to the unit simplex $`\mathrm{\Sigma }_{q+m}`$.

$`3^0`$. In what follows we assume that the variables $`t_j`$ belong to some interval $`[0,N]`$ where $`N>1`$. We define with the help of the variables $`h_j`$, $`t_i`$ continuous functions $`\eta _j`$, $`\zeta _i`$ such that $`(\eta _1,\mathrm{},\eta _{q+m})\mathrm{\Sigma }_{q+m}`$, $`\zeta _i[0,N]`$. The set $`𝒮=\mathrm{\Sigma }_{q+m}\times [0,N]^m`$ is homeomorphic to $`\mathrm{\Sigma }_{q+2m}`$. By the Brouwer fixed point theorem (see \[Do\], p. 57), there exists a fixed point of the mapping $`\tau :𝒮𝒮`$, $`\tau :(h,t)(\eta ,\zeta )`$, i.e. a point where one has $`\eta _j=h_j`$, $`\zeta _i=t_i`$. The functions $`\eta _j`$, $`\zeta _i`$ are defined such that the arrangement of the real roots of $`P`$ and $`P^{(s)}`$ at the fixed point is the required one.

$`4^0`$. Define the functions $`\eta _j`$ by the following rules:

1) We want to achieve the additional conditions (at the fixed point) $`g_p=u_{2p1}`$, $`\mathrm{}`$, $`g_{p+p^{}}=u_{2p+2p^{}1}`$ for all appropriate indices, see $`2^0`$; therefore we set $`\eta _{i_1}=\xi _{i_2}`$ whenever $`h_{i_1}`$ is a variable $`g_{p+l}`$ and $`\xi _{i_2}`$ is the corresponding function $`u_{2p+2l1}`$;

2) If a variable $`h_j`$, which is a root $`w_i`$ of multiplicity $`<s+1`$, must coincide with a simple root $`\xi _k`$ of $`P^{(s)}`$ or, more generally, with the roots $`\xi _k=\xi _{k+1}=\mathrm{}=\xi _{k+l}`$, then we set $`\eta _j=\xi _k`$;

3) If the variables $`h_r<h_{r+1}<\mathrm{}<h_{r+l}`$ (which are all consecutive roots $`w_j`$ and among which there might be roots $`w_j`$ of multiplicity $`s+1`$) lie between the Rolle roots $`\xi _k`$ and $`\xi _{k+v}`$ of $`P^{(s)}`$ and all roots among the roots $`\xi _{k+1}`$, $`\mathrm{}`$, $`\xi _{k+v1}`$ (if $`v>1`$) coincide with roots $`w_j`$ ($`rjr+l`$) of multiplicity $`s+1`$, then we set

$`\eta _{r+j}=\xi _k+(j+1)(\xi _{k+v}\xi _k)/(l+2)`$, $`j=0,1,\mathrm{},l`$.

###### Remark 7

It follows from rules 1) – 3) that there are $`q+m`$ functions $`\eta _j`$ – as many as the variables $`h_j`$.

Recall that the arrangement is least generic, i.e. for every non-Rolle root $`\xi _i`$ of $`P^{(s)}`$ one has either $`\xi _i=\xi _{i_1}`$ where $`\xi _{i_1}`$ is a Rolle one or $`\xi _i=w_{i_2}=h_j`$ for some $`i_2`$, $`j`$. Denote by $`l_1`$, $`\mathrm{}`$, $`l_{2m}`$ the absolute values $`|\xi _i\xi _{i_1}|`$ and $`|\xi _iw_{i_2}|`$ for all $`i`$, $`i_1`$ and $`i_2`$ as above. Set $`\mathrm{\Phi }=l_1+\mathrm{}+l_{2m}`$ and

$$\zeta _i=|t_i\frac{1}{3m}\underset{j=1}{\overset{m}{}}\theta _j\frac{t_i}{3(N+1)^m}|t_1t_2\mathrm{}t_m1|\frac{t_i}{12m}\mathrm{\Phi }|$$
(3)
$`5^0`$. Denote by $`t_{i_0}`$ the greatest variable $`t_i`$ at the fixed point (see $`3^0`$). Observe first that one can assume that $`t_{i_0}>0`$. Indeed, if $`t_{i_0}=0`$, then $`t_i=0`$ for all $`i`$, $`P`$ is hyperbolic and the roots of $`P`$ and $`P^{(s)}`$ define an arrangement $`\alpha `$ from the closure of the desired least generic one $`\beta `$.

###### Lemma 8

If $`t_{i_0}=0`$, then there exists a real-analytic deformation of $`P`$ into a real polynomial which together with its $`s`$-th derivative defines the arrangement $`\beta `$.

The lemma is proved after the theorem. It allows one to consider only the case $`t_{i_0}>0`$.

$`6^0`$. One has

$$\zeta _{i_0}=t_{i_0}\frac{1}{3m}\underset{j=1}{\overset{m}{}}\theta _j\frac{t_{i_0}}{3(N+1)^m}|t_1t_2\mathrm{}t_m1|\frac{t_{i_0}}{12m}\mathrm{\Phi }.$$

Indeed, all roots of $`P^{(s)}`$ lie within the convex hull of all roots of $`P`$ (see \[PoSz\], p. 108). Hence, one has $`\theta _jt_{i_0}`$, $`j=1,\mathrm{},m`$. One has also $`|t_1t_2\mathrm{}t_m1|t_1t_2\mathrm{}t_m+1<(N+1)^m`$ and $`\mathrm{\Phi }4m`$ (because for each term $`l_j`$ one has $`l_j2`$). Thus

$$\frac{1}{3m}\underset{j=1}{\overset{m}{}}\theta _j+\frac{t_{i_0}}{3(N+1)^m}|t_1t_2\mathrm{}t_m1|+\frac{t_{i_0}}{12m}\mathrm{\Phi }<mt_{i_0}/3m+t_{i_0}/3+4mt_{i_0}/12m=t_{i_0}$$
(4)
and for $`i=i_0`$ one can delete the absolute value sign in the right hand-side of (3). But then to have $`\zeta _{i_0}=t_{i_0}`$ one must have $`\theta _j=0`$ for $`j=1,\mathrm{},m`$, $`t_1t_2\mathrm{}t_m1=0`$ and $`l_1=\mathrm{}=l_{2m}=0`$. This means that $`t_j0`$, i.e. no root $`g_j+it_j`$ of $`P`$ will be real, that $`P^{(s)}`$ will indeed be hyperbolic ($`\theta _j=0`$) and that all non-Rolle roots of $`P^{(s)}`$ equal either roots $`w_j`$ of $`P`$ or Rolle roots of $`P^{(s)}`$.

###### Remark 9

The condition $`N>1`$ makes possible the choice of the values of the variables $`t_i`$ so that $`t_1t_2\mathrm{}t_m1=0`$. One can prove by analogy with (4) that $`|\zeta _i|<N`$, i.e. the mapping $`\tau `$ is indeed from $`𝒮`$ into itself.

$`7^0`$. A priori the fixed point assures the existence of an arrangement only from the closure of the necessary one. The fact that at the fixed point no inequality between roots of $`P`$ is replaced by equality is proved by analogy with $`4^0`$$`7^0`$ of the proof of Theorem 4.4 from \[KoSh\] where the case of $`P`$ hyperbolic is considered. The proof there shows that equalities replacing inequalities between roots of $`P`$ imply that a root of $`P`$ of multiplicity $`ms+1`$ is a root of $`P^{(s)}`$ of multiplicity $`ms+1`$ which contradicts part 1) of Proposition 4. In the general case ($`P`$ not necessarily hyperbolic) the proof is essentially the same, the presence of eventual non-Rolle roots can only increase the multiplicity of the root as a root of $`P^{(s)}`$.

Hence, the fixed point provides the necessary arrangement.

$`8^0`$. To obtain (in $`8^0`$$`9^0`$) all arrangements in which $`P^{(s)}`$ is hyperbolic but which are not necessarily least generic we use the same construction but with another function $`\mathrm{\Phi }`$. Namely, consider a family of such functions $`\mathrm{\Phi }`$ depending on a parameter $`b(𝐑_+,0)`$ defined as follows: if instead of $`\xi _i\xi _{i_1}=0`$, see $`4^0`$, one must have $`\xi _i\xi _{i_1}>0`$ or $`\xi _i\xi _{i_1}<0`$ (and no root $`\xi _j`$ or $`w_j`$ lies between $`\xi _i`$ and $`\xi _{i_1}`$), then in $`\mathrm{\Phi }`$ we replace the absolute value $`l_\nu =|\xi _i\xi _{i_1}|`$ by $`|\xi _i\xi _{i_1}b|`$ (resp. by $`|\xi _i\xi _{i_1}+b|`$); in the same way for $`\xi _iw_{i_2}`$, see $`4^0`$. In a sense, we obtain the not least generic arrangements by deforming least generic ones the deformation parameter being $`b`$.

$`9^0`$. Denote by $`F(b)`$ the set of fixed points of the mapping $`\tau `$ from $`3^0`$. For $`b`$ small enough one has $`(\eta ,\zeta )𝒮`$. The set $`F(0)`$ contains all limit points of the family of sets $`F(b)`$ when $`b0`$ and there exists at least one such limit point because all sets $`F(b)`$ (for $`b`$ small enough) are non-empty and belong to $`𝒮`$ which is compact. Hence, one can choose $`b>0`$ small enough and a fixed point of $`F(b)`$ at which there is an inequality between two roots in the arrangement if there is an inequality in the arrangement for $`b=0`$, and the equalities $`\xi _i\xi _{i_1}=0`$ or $`\xi _iw_{i_2}=0`$ where this is necessary are replaced by the desired inequalities.

$`10^0`$. Obtain all arrangements in which $`P^{(s)}`$ is not hyperbolic and which are least generic. Suppose that $`P^{(s)}`$ must have exactly $`m^{}`$ conjugate couples of complex roots. In this case we assume that $`m^{}`$ of the couples of roots $`g_j\pm it_j`$ are replaced by a couple $`\pm iv`$ where $`v>0`$ is “large”, i.e. much bigger than $`N`$. Hence, $`P^{(s)}`$ also has exactly $`m^{}`$ couples of conjugate complex roots with “large” imaginary parts. One has

$$Q:=P/v^{2m^{}}=(1+x^2/v^2)^m^{}\underset{j=1}{\overset{q}{}}(xw_j)^{m_j}\underset{j=1}{\overset{mm^{}}{}}((xg_j)^2+t_j^2),$$

i.e. the family $`Q`$ is a one-parameter deformation of a family of polynomials like (2) (the role of the small parameter is played by $`1/v^2`$) and the existence of the necessary arrangements can be deduced by analogy with $`1^0`$$`7^0`$ (see $`9^0`$ for the role of the small parameter; however, the function $`\mathrm{\Phi }`$ is the one from $`1^0`$$`7^0`$).

$`11^0`$. To obtain the existence of all arrangements (which are not necessarily least generic and with $`P^{(s)}`$ not necessarily hyperbolic) one has to combine $`8^0`$, $`9^0`$ and $`10^0`$. The theorem is proved. $`\mathrm{}`$

Proof of Lemma 8:

$`1^0`$. We assume that $`P`$ has the same number of distinct real roots as in the desired arrangement $`\beta `$. If not, then one can first deform $`P`$ within the class of hyperbolic polynomials (while remaining in the closure of $`\beta `$) to achieve this condition. See \[Ko\] for such deformations.

We begin with two observations:

1) for $`a>0`$, $`\mu 𝐍\{0\}`$ and $`\nu `$ even the polynomial $`Q=x^\mu (x^\nu +a)`$ has a $`\mu `$-fold root for $`x=0`$ and its $`s`$-th derivative for $`s>\mu `$ has a $`(\mu +\nu s)`$-fold one; $`Q`$ has also $`\nu /2`$ couples of conjugate complex roots;

2) with $`a`$, $`\mu `$ and $`\nu `$ as above, the polynomial $`Q_1=x^\mu (x^\nu +a+aQ_2(x,a))`$ where $`Q_2`$ is a polynomial in $`x`$ of degree $`\nu 1`$, $`Q_2(0,a)0`$, has $`\nu `$ complex zeros for $`a`$ small enough and a real $`\mu `$-fold root at $`0`$; to see this set $`a=c^\nu `$, $`x=cy`$; one has $`Q_1(cy,c^\nu )=c^{\mu +\nu }y^\mu (y^\nu +1+Q_2(cy,c^\nu ))`$; the last polynomial has a $`\mu `$-fold root at $`0`$ and $`\nu `$ roots which for $`c`$ small enough are close to the roots of $`y^\nu +1`$, hence, are complex.

$`2^0`$. Suppose that the polynomial $`P`$ of degree $`n`$ realizing with $`P^{(s)}`$ the arrangement $`\alpha `$ has a real root of multiplicity $`\mu +\nu `$ (with $`\nu `$ even) which (in order to obtain the arrangement $`\beta `$) must split into $`\nu /2`$ couples of conjugate complex roots and into a real root of multiplicity $`\mu `$. (If several roots of $`P`$ must split, we make them split one by one.) Suppose in addition that in the deformed polynomial (denoted by $`R`$) the real root of multiplicity $`\mu `$ must coincide with a root of $`R^{(s)}`$ of multiplicity $`\mu +\nu s`$. Assume that the bifurcating root is at $`0`$ and that

$$P=x^{\mu +\nu }(1+h(x)),h(0)=0$$
(5)
($`P`$ is not necessarily monic). Construct the necessary deformation of $`P`$ in the form

$$R(x,a)=x^\mu (x^\nu +a+b_{s\mu }x^{s\mu }+\mathrm{}+b_{\nu 1}x^{\nu 1})(1+g(x,a))$$
(6)
where $`a(𝐑,0)`$ and $`b_i=b_i(a)`$ and $`g(x,a)`$ ($`g(0,a)0`$) are defined such that all equalities of the form $`x_i=\xi _j`$ defining the arrangement $`\beta `$ will be preserved.

$`3^0`$. Suppose first that in (6) one has $`g(x,a)h(x)`$. The condition

$$(A):R^{(s)}\mathrm{has}\mathrm{a}(\mu +\nu s)\mathrm{fold}\mathrm{root}\mathrm{at}0$$

is a triangular linear non-homogeneous system with unknown variables $`b_i`$; the system defines unique functions $`b_i=b_i^{}a`$, $`b_i^{}𝐑`$. This can be checked directly.

Suppose that in (6) one has $`g=h(x)+_{j=1}^ld_jh_j(x,d)`$ where $`d=(d_1,\mathrm{},d_l)(𝐑^l,0)`$ and $`h_j`$ depend smoothly on $`d`$. Then condition (A) defines unique functions $`b_i(a,d)=b_i^{}a+a_{j=1}^ld_j\stackrel{~}{b}_{i,j}(d)`$ where $`b_i^{}𝐑`$ and $`\stackrel{~}{b}_{i,j}`$ are smooth in $`d`$. This can also be checked directly.

$`4^0`$. For each root $`w_j0`$ of $`P`$ of multiplicity $`<s`$ which must be equal to a root $`\xi _i`$ of $`P^{(s)}`$ denote by $`d_j`$ the deviation from its position in a deformation of $`P`$. Admitting such deviations means that in (5) the function $`h`$ should be replaced by $`h(x)+_{j=1}^ld_jh_j(x,d)`$.

Denote by (B) the system of all conditions $`w_j=\xi _i`$ for all such equalities with $`w_j0`$ characterizing the arrangement $`\beta `$.

$`5^0`$. For any deformation $`R^{}(x,a,d)=x^\mu (x^\nu +a+b_{s\mu }x^{s\mu }+\mathrm{}+b_{\nu 1}x^{\nu 1})(1+g(x,d))`$ of $`P`$ (where $`b_k`$ are considered as small parameters) one can find $`d`$ depending smoothly on $`a`$ and $`b_k`$ such that for all $`a`$ small enough all equalities from (B) hold. This follows from Propositions 11 and 13 from \[Ko\] where it is shown that the linearizations of the conditions (B) w.r.t. $`d`$ are linearly independent. (In \[Ko\] their linear independence is proved only when $`P`$ is hyperbolic; this independence is an “open” property, so it holds for all nearby polynomials as well.)

$`6^0`$. The independence of these linearizations implies that for $`a`$ small enough the system of conditions (B) applied to the deformation

$$\stackrel{~}{R}(x,a,d)=x^\mu (x^\nu +a+b_{s\mu }(a,d)x^{s\mu }+\mathrm{}+b_{\nu 1}(a,d)x^{\nu 1})(1+h(x)+\underset{j=1}{\overset{l}{}}d_jh_j(x,d))$$

(with $`b_i(a,d)`$ defined as in $`3^0`$) defines unique $`d_j=d_j(a)`$ smooth in $`a`$. Indeed, the linearizations w.r.t. $`d`$ of the system of conditions (B) from $`6^0`$ and from $`5^0`$ are the same.

On the other hand, $`b_i`$ were defined such that condition (A) holds. Hence, for $`d=d(a)`$ and $`b_i=b_i(a,d(a))`$ (where $`a>0`$ is small enough) the $`(\mu +\nu )`$-fold root of $`P`$ at $`0`$ splits into a real $`\mu `$-fold root at $`0`$ and $`\nu `$ complex roots close to $`0`$ (see observation 2) from $`1^0`$) and $`P^{(s)}`$ has a $`(\mu +\nu s)`$-fold root at $`0`$. The arrangement of the other real roots of $`P`$ and $`P^{(s)}`$ remains the same. $`\mathrm{}`$

Acknowledgement. The present text is written under the impetus of the fruitful discussions with B.Z. Shapiro, the author’s host during a short visit at the University of Stockholm. The author thanks both him and his university for the kind hospitality.

Author’s address: Université de Nice, Laboratoire de Mathématiques, Parc Valrose, 06108 Nice Cedex 2, France. e-mail: kostov@math.unice.fr


format: text

text:
# Large bulk matter search for fractional charge particles

## I INTRODUCTION

We have carried out the largest search for fractional electric charge elementary particles in bulk matter using 70.1 mg of silicone oil. That is, we looked for stable particles whose charge $`Q`$ deviates from $`Ne`$ where $`N`$ is an integer, including zero, and $`e`$ is the magnitude of the charge on the electron. No evidence for such particles was found in this amount of silicone oil. We used our new version loomba  of the Millikan oil drop method containing two innovations compared to the classical method that we used in Halyo et al.halyo . One innovation is that the drop charge is obtained by observing the drop motion in a *horizontal*, alternating electric field compared to the classical use of a vertical electric field halyo ; mar ; savage . The other innovation is the use of an upward flow of air to reduce the vertical terminal velocity of the drop, which enabled us to use larger drops, about 20.6 $`\mu `$m in diameter compared to the 10 $`\mu `$m drops used in our previous experiments.

We define the residual drop charge, $`Q_r=QN_le`$ where $`N_l`$ is the largest integer less than $`Q/e`$. We find the 95% confidence level upper limit on the abundance of fractional charge particles in silicone oil for $`0.18eQ_r0.82e`$ is $`1.17\times 10^{22}`$ particles per nucleon. This experiment was a follow up on our previous search in silicone oil, Halyo et al.halyo , based on 17.4 mg. In that search we found one drop with anomalous charge, but no such charge was found in the present experiment.

In this paper we describe the experimental method and apparatus in Sec. II. In Sec. III we discuss the measurement precision resulting from the various measurement errors and the calibration methods. The data analysis method, including the criteria used to accept drop charge measurements, is discussed in Sec. IV. Here we pay particular attention to the drop spacing criterion necessitated by interactions between adjacent drops. This is the primary limitation on the rate at which drops can be measured and we had to acquire considerable experience to understand this limitation. We conclude with Sec. V, giving our results, comparing our results with other fractional charge searches, and discussing the applicability and extension of this new Millikan oil drop technique to other searches.

## II EXPERIMENTAL METHOD AND APPARATUS

### II.1 Experimental method

The principle of the experimental method is simple. Consider a drop of radius $`r`$, density $`\rho `$, and charge $`Q`$ falling in air through a horizontal electric field of strength $`E`$, as shown in Fig. 1. Applying Stokes’ law the horizontal terminal velocity, $`v_x`$ is

$$v_x=\frac{QE}{6\pi \eta r}$$
(1)
where $`\eta `$ is the viscosity of air. Hence measuring $`v_x`$ gives $`Q`$ providing $`r`$ is known. As explained in Sec. III.1, the drop radius is determined from the $`v_x`$ of integer charge particles. Note that the measurement of $`Q`$ does not depend on the density of the drop and is also independent of the gravitational force on the drop. The electric field alternates in the $`+x`$ and $`x`$ direction so that the drop is moved back and forth along the $`x`$ axis. This cancels some sources of error and allows the drop motion to be viewed in a relatively narrow horizontal area. Previous uses of a static horizontal electric field in the Millikan oil drop method were in 1941 by Hopper and Laby hopper  who measured the electron charge and by Kunkel kunkel  in 1950 who measured the charge on dust particles.

If the drop were falling in still air, the vertical terminal velocity would be given by

$$v_{z,term}=\frac{2r^2\rho g}{9\eta }$$
(2)
where g is the acceleration of gravity. However we use an upward flow of air of velocity $`v_{air}`$ in the $`z`$ direction. Hence the net downward velocity of the drop is

$$v_z=\frac{2r^2\rho g}{9\eta }v_{air}.$$
(3)
As explained in the next section we want $`v_z`$ to be small hence we set $`v_{air}`$ to be close to $`v_{z,term}`$ but slightly smaller.

### II.2 General description of experiment

Figure 2 is a schematic picture of the apparatus. Drops averaging 20.6 $`\mu `$m in diameter are produced at a rate of 1 Hz using a piezoelectrically actuated drop-on-demand microdrop ejector. The drops fall through the upward moving air in the measurement chamber passing through a horizontal, uniform, alternating electric field. In this figure the electric field is perpendicular to the paper. The electric field alternates as a square wave with a frequency of 2.5 Hz and has an amplitude of about $`1.8\times 10^6`$ V/m.

A rectangular measurement region 2.29 mm in the $`x`$ direction by 3.05 mm in the $`z`$ direction is projected by a lens onto the charge-coupled device (CCD) sensor of a monochrome, digital video camera. A light source consisting of a bank of light emitting diodes (LEDs) provides 10 Hz stroboscopic illumination. As the motion of the drop carries it through the measurement region, its image appears on the surface of the CCD. Thus the camera collects 10 frames per second, the drop appearing as a dark image on a bright background.

In addition to the $`v_x`$ motion there is also the $`v_z`$ motion, Eq. (3). Since the camera has a field of view of $`Z=3.05`$ mm in the vertical direction, the 10 Hz stroboscopic illumination leads to acquisition of

$$N_{images}=10Z/v_z=30.5/v_z$$
(4)
images of any given drop, before the drop moves below the viewing area of the camera. Here $`v_z`$ is in mm/s. Hence we get a larger number of images per drop, leading to better charge measurement precision, when $`v_z`$ is small. Of the order of $`N_{images}=15`$ are required.

We give an example of the importance of the upward airflow in obtaining this many images. Consider a typical drop of diameter 20.6 $`\mu `$m with a density of 0.913 g/cm3. From Eq. (2), $`v_{z,term}=11.3`$ mm/s. If there were no upward airflow there would be an average of 2.6 images per drop. To obtain $`N_{images}=15`$, $`v_z`$ must be about 2.0 mm/s. Therefore from Eq. (3), $`v_{air}`$ must be $`11.32.0=9.3`$ mm/s.

Each image from the CCD camera is processed though a framegrabber in a conventional desktop computer, the signal in each pixel being recorded. An analysis program then finds the drop images and calculates the $`x`$ and $`z`$ coordinates of the centroid of the drop image. Using all the images of the drop and knowing the time spacing of the images, $`v_x`$ and $`v_z`$ are then calculated.

### II.3 Drop generator

The drop generator, Fig. 3, is based on the general design principles used in piezoelectrically actuated, drop-on-demand, inkjet print heads. Our generator is designed for flexibility, allowing a variety of liquids to be used and providing ease of control and maintenance. The body of the generator is glass so as to preserve the purity of the liquid, but the ejection aperture at the bottom is micromachined silicon patent . The diameter of the aperture sets approximately the diameter of the drop. Upon application of a short voltage pulse, usually 3 to 20 $`\mu `$s, across the surfaces of the piezoelectric disk, the central hole in the disk contracts in diameter, squeezing the glass tube and sending a pressure pulse down the liquid, ejecting a jet of fluid from the aperture. The forming of a discrete fluid drop from the high speed jet is a complex process with the repeatability of the process and the final diameter of the drop being highly dependent on the properties of the fluid, and on how the fluid is driven. For this reason, the shape and amplitude of the voltage pulse applied to the piezoelectric disk must be specifically tuned for stability and the desired drop size. In addition, it is necessary to experiment with single and double pulsing, varying both the pulse width and the separation between the pulses.

The pressure in the drop generator is maintained slightly below atmospheric pressure by 10 to 30 mm Hg. This helps to retract the excess liquid outside the ejection aperture after the drop has been produced and also prevents leaking of the liquid between pulses.

It is important that the drop generator produce drops of constant radius, the primary reason being that the size of the drops determines the ratio between $`Q`$ and the measured quantity $`v_x`$, from Eq. (1). A secondary reason that we did not initially appreciate has to do with Eqs. 2 and 3. We set $`v_{air}`$ close to $`v_{z,term}`$ so that $`v_z`$ is much smaller than $`v_{z,term}`$. Thus a small change in $`r`$ leads to a relatively large change in $`v_z`$. This dispersion disrupts the consistent spacing between adjacent drops, and a decreased separation between drops is undersirable for reasons discussed in Sec. IV.1.6.

With a clean and newly tuned drop generator we get remarkably uniform drop radii, constant to about $`\pm 0.2\%`$. The drop generator also ejects in a consistent downward vertical direction along the centerline of the airflow tube. At 1 Hz operation, the drop ejector exhibits slow drifts in its characteristics with time scales of the order of a week. These drifts appear as changes in drop size and as destabilization that manifests as the appearance of satellite drops or inconsistent drop production. Typically, these effects can be compensated for by small changes in the drive parameters or adjustments of the air velocity or both. By the end of the first data set, set 1, the drop ejector had destabilized to the point where it had to be removed from the apparatus for cleaning and refilling. Similarly during the taking of data set 2, the drop ejector and air velocity required periodic small adjustments. The end of set 2 was caused by increasing instability in drop production which could not be compensated for. We do not know the reason for this behavior.

In our drop generator the silicone oil drops are produced with a spread of charges, $`|Q|`$ ranging mostly from 0 to about 10 $`e`$. A few percent of the drops have larger $`|Q|`$. As described in Sec. IV.1.1 we used drops with $`|Q|<9.5e`$ to maintain good precision in the charge measurement. We do not know what sets the charge distribution for a particular drop generator. But we have the general observation that silicone oil gives narrow charge distributions, whereas water, mineral oil and most other organic liquids give broad charge distributions, with $`|Q|`$ values as large as several 1000 $`e`$ or even larger.

### II.4 Optical system

Referring to Fig. 2, the stroboscopic light source consists of a rectangular bank of 20 LED’s emitting at 660 nm. The pulse length was about 100 $`\mu `$s. The lens, a 135 mm focal length, $`f`$/11, photographic enlarging lens, images the measurement region onto the face of the CCD camera with a magnification of 2.1.

The rectangular active image area of the CCD camera camera  is 4.8 mm in the horizontal, that is, $`x`$, direction, and 6.4 mm in the vertical, that is, $`z`$, direction. Hence the viewing area in physical space is 2.29 mm horizontally by 3.05 mm vertically. We remind the reader that the electric field is horizontal. The active imaging area is an array of 240 horizontal picture elements (pixels) and 736 vertical picture elements (pixels). We chose this orientation of the array to maximize the vertical distance, maximizing the number of images per drop.

Given the magnification and pixel density of the CCD, one would expect from geometric optics that the shadow of a 20 $`\mu `$m diameter drop would cover 2 pixels horizontally and 5 pixels vertically. The actual observed shadow typically covered 3 pixels horizontally and 7 pixels vertically, and had an intensity variation that was approximately a two dimensional Gaussian. This can be quantitatively described as the convolution of the simple shadow predicted by geometric optics with a point spread function that is a result of the diffractive effects due to the finite aperture of the lens. We do not and should not observe diffractive effects caused by the small size of the drops.

### II.5 Airflow tube and measurement chamber

Figure 4 shows a slightly simplified, dimensioned drawing of the airflow tube and the measurement chamber. A rectangular duct contains the upward flowing air. It is 8.3 mm wide in the direction of the electric field and 31.8 mm wide in the direction perpendicular to the electric field. The field plates that define the measurement chamber are 51 mm high and 28.6 mm wide. The inner surfaces of the plates are machined flat and are in the same plane as the inner surfaces of the walls of the airflow tube. The optic axis of the optical system passes through the transparent side walls of the airflow tube.

The air velocity is sufficiently small, with a Reynolds number on the order of $`R_e=50`$ so that the flow is laminar. The 203 mm length of air flow tube between the measurement region and the air inlet allows the air to settle into its equilibrium flow pattern. At equilibrium, the velocity profile of the air is approximately parabolic across the narrow direction of the channel ($`x`$ axis). Across the long axis, the flow has a roughly constant central region and falls to zero at the boundaries langlois .

### II.6 Electronics

All the electronics of the apparatus, Fig. 5, are hard-wired to give reliable timing, independent of the operation of the computer. A 30 Hz handshaking signal from the CCD is divided down to provide a 10 Hz clock that synchronizes the LED strobe, the electric field switcher and the computer image acquisition. The switching of the electric field, which is driven by the clock signal divided by 4, operates at 2.5 Hz. This results in a cycle where two images are acquired with the electric field in one direction, and then two images with the electric field in the other direction. These relationships between the signals is depicted in the timing diagram of Fig. 5. The drop generator is driven asynchronously at 1 Hz.

### II.7 Data acquisition and storage

Data acquisition was performed by a single desktop computer running Linux. The computer was equipped with two special components: a digital framegrabber that allowed the capture of image data from the camera and a general purpose input/output interface board with digital I/O and A/D conversion capability. The additional inputs allowed the computer to monitor the state of the experiment as well as a variety of environmental variables.

All software was custom written in C. Hardware dependent code was encapsulated into drivers at the kernel level, which allowed a guarantee of synchronization of the software with apparatus by a combination of hardware and software buffering of the data.

The overall strategy was to acquire data from the apparatus and write them to files in *raw* form for later processing off line. Recall that each image frame contains $`736\times 240`$ pixels, each digitized to 8 bit accuracy. Therefore acquiring data at 10 Hz produces a data rate from the CCD camera of about 2 MB/s, much too large to be stored. Since an image contains just a few drops, most of the pixels in an image have just the background signal, which allows the information to be stored using sparse storage. As mentioned in Sec. II.4, the typical image of drop extends over an area of 3 pixels by 7 pixels. In each frame the regions containing drops were isolated using a thresholding operation. The position of each drop was then measured using a simple center of mass algorithm, and for each drop only a surrounding region containing 13 horizontal pixels by 21 vertical pixels of the image is written to the output file.

### II.8 Data collection

The search was carried out in two sets described in Table 1. At the beginning of the Set 1 the drop ejector was operated at 0.5 Hz, then at 1 Hz for the remainder of Set 1 and all of Set 2.

## III CALIBRATION, ERRORS, AND MEASUREMENT PRECISION

### III.1 Electric field and drop radius

Rewrite Eq. (1) in the form

$$v_x=\left(\frac{Q}{6\pi \eta }\right)\left(\frac{E}{r}\right)$$
(5)
Consider nonfractional values of $`Q=ne`$, $`n=0,\pm 1,\pm 2,\mathrm{}`$. Then, as shown in Sec. IV.2, the measured values of $`v_x`$ sharply peak at $`n=0,\pm 1,\pm 2,\mathrm{}`$. From a fit to the center of these peaks, Eq. (5) gives the fitted $`E/r`$ ratio.

The electric field strength, $`E`$, is calculated from the measured voltage across the electric field plates, known to 3%, and the plate separation, 8.25$`\pm `$0.01 mm. The plates are parallel to within 0.1 mrad. Inserting the calculated value of $`E`$ into the fitted $`E/r`$ ratio, we obtain the drop radius $`r`$.

We have two additional checks of the drop radius, one from the measurement of the error caused by the Brownian motion, Sec. III.2, and the other from the measurement of the net downward velocity of the drops, $`v_z`$ in Eq. (3).

The drop radius depends to a moderate extent upon the size and shape of the voltage pulse applied to the drop generator and to a slight extent upon the age and history of the drop generator. However over periods of hours the average drop radius could be taken to be constant, with fluctuations of $`\pm 0.2\%`$ for individual drops. Since the data were analyzed in one-hour-long blocks, the average $`E/r`$ ratio for any given block was known to much better accuracy than this.

### III.2 Brownian motion and drop position measurement errors

The precision of the determination of the drop charge depends upon the precision of the measurement of $`v_x`$. Consider the sequence of position measurements $`x_i`$ of the trajectory of a drop. For two consecutive frames, $`j`$ and $`j1`$, the velocity measurement $`v_{x,j}`$ is given by

$$v_{x,j}\frac{x_jx_{j1}}{\mathrm{\Delta }t}.$$
(6)
Here, $`\mathrm{\Delta }t`$ is the time between successive frames, 0.1 s in our case. Since $`\mathrm{\Delta }t`$ is known with very good precision, the error in measuring $`v_x`$ comes from the error in determining the $`x_i`$ of the drop centers, and from Brownian motion. Take the error in centroiding to be normally distributed with a standard deviation of $`\sigma _c`$.

During the time $`\mathrm{\Delta }t`$ between any two successive measurements of the $`x_i`$ positions, Brownian motion adds a random contribution with standard deviation given by

$$\sigma _b=\sqrt{\frac{kT\mathrm{\Delta }t}{3\pi \eta r}}.$$
(7)
Here $`k`$ is the Boltzmann constant, $`T`$ is the absolute temperature, $`\eta `$ is the viscosity of air and $`r`$ is the drop radius.

The trajectory of uncharged drops, which have no contribution to their trajectories due to the electric field, can thus be written as

$`x_j=x_0+\sigma _{c,j}+{\displaystyle \underset{i=1}{\overset{j}{}}}\sigma _{b,i}`$ (8)
$`j=1,2,\mathrm{},N_{images}`$
with $`x_0`$ set by the initial position of the drop, and where the $`\sigma _{c,i}`$ ($`\sigma _{b,i}`$) are normally distributed with a std. dev. of $`\sigma _c`$ ($`\sigma _b`$). The analysis for charged drops is similar if the effect of the electric field is first subtracted from the observed data points. It then follows that

$$v_{x,j}\mathrm{\Delta }t=\sigma _{b,j}+\sigma _{c,j}\sigma _{c,j1}$$
(9)
and

$$v_{x,j}v_{x,k}\mathrm{\Delta }t^2=\{\begin{array}{cc}2\sigma _c^2+\sigma _b^2,& j=k\hfill \\ \sigma _c^2,& |jk|=1\hfill \\ 0,& \text{otherwise.}\hfill \end{array}$$
(10)
Therefore the total error on any given velocity measurement, $`\sigma _v`$, is given by $`\sigma _v^2=2\sigma _c^2+\sigma _b^2`$, and the centroiding error introduces an anticorrelation with magnitude $`\sigma _c^2`$ in two consecutive velocity measurements, due to the shared position measurement. We use the concept summarized in Eq. (10) and the observed distributions of the $`v_{x,i}`$ (after removal of the contribution due to the alternating electric field) to separate $`\sigma _c`$ from $`\sigma _b`$.

We find that averaged over this experiment

$$\sigma _c=0.31\text{ }\mu \text{m}\text{}\sigma _b=0.47\text{ }\mu \text{m},$$
(11)
in the measurement region. Compared to the size of an individual pixel on the CCD, the centroiding error is small, approximately 1/30 of a pixel. The value of $`\sigma _b`$ obtained provides an independent check on the size of the drops, and is consistent with the size determined from the terminal velocity and the electric field drift velocity.

Equation 11 shows that the Brownian motion error, $`\sigma _b`$, is about the same magnitude as the error involved in finding the drop position, $`\sigma _c`$. Therefore substantially reducing $`\sigma _c`$ through the use of smaller pixels will not by itself substantially reduce the error on the charge measurement, since the Brownian motion error can only be reduced by increasing $`N_{images}`$.

The final charge measurement of a drop is made using a single, detailed best fit to the entire observed trajectory of the drop, and the final error on the charge measurement $`\sigma _q`$ is a result of propagating the errors $`\sigma _c`$ and $`\sigma _b`$ through this calculation.

### III.3 Other sources of errors

We looked for other sources of errors, but all are negligible compared to those in Eq. (10). When we developed the upward air flow method we thought about the possibility that there might be some horizontal air velocity, $`v_{x,air}`$ in the measurement chamber, contributing an error to $`\sigma _v`$ of order $`v_{x,air}\times \mathrm{\Delta }t`$. By studying a large amount of data we found that the distribution of $`v_{x,air}\times \mathrm{\Delta }t`$ had an rms value of 100 nm, and was a fixed property of the measurement region. For comparison $`v_e\times \mathrm{\Delta }t`$ was of the order of 8 $`\mu `$m. These irregularities in $`v_{x,air}`$ are probably due to residual surface imperfections in the electric field plates. Since the irregularities are constant over long periods of time, they can be accurately measured and corrected for. For this analysis, that was not necessary.

Another possible source of error would be a nonuniformity in the electric field in the measurement region giving a horizontal gradient, $`dE/dx`$. This would produce a horizontal force on the drop’s induced electric dipole moment. This dipole force acts in addition to the $`QE`$ force. We found such a dipole force to be negligible compared to the $`QE`$ force.

A small, vertical deceleration of the drops as they fall through the measurement chamber was observed. This amounted to a change of 30 $`\mu `$m/s in the apparent terminal velocity of the drops as they fell through the measurement region, or a systematic uncertainty in the radius of the drop of the order of 0.3%. We believe that the deceleration is due to the evaporation of the drop as it falls. The magnitude of this effect is small enough such that it can be neglected in the calculation of $`v_x`$. As a side note, any systematic uncertainties in the radius of the drops are absorbed by the calibration process described in Sec. III.1, and do not affect the final charge measurement.

Similarly, other possible sources of measurement error such as apparatus vibrations, optical distortions and CCD array distortions, and patch nonuniformities on the electric field plates, were negligible.

## IV DATA ANALYSIS AND RESULTS

### IV.1 Drop selection criteria

In this section we use $`q=Q/e`$, a measure of the drop charge in units of the electron charge. We required that all drops used in the analysis meet the criteria in Table 2. The criteria are designed to maintain a charge measurement accuracy of approximately 0.03 $`e`$ and to reject irregular drops caused by inconsistent operation of the drop generator.

#### IV.1.1 $`q<9.5`$ criterion

For any given drop there is an uncertainty in the radius of approximately 0.2% which contributes to the relative error on $`q`$. The *absolute* error on $`q`$ thus increases linearly with $`q`$. Since the absolute error on $`q`$ must be kept to the order of 0.03 $`e`$, restricting the data sample to drops with $`q<9.5`$ keeps this contribution to less than 0.02 $`e`$. The overall charge distribution is such that only a few percent of the drops have $`q`$ values outside this range.

#### IV.1.2 $`\sigma _q<0.03`$ criterion

Primarily, this criterion is a measure of $`N_{images}`$ of the drop. Brownian motion and centroiding accuracy, characterized by $`\sigma _c`$ and $`\sigma _b`$ as described in Sec. III.2, limit the accuracy of the charge measurement, $`\sigma _q`$. For any given drop, the number of position measurements, $`N_{images}`$, and the state of the electric field during those measurements, in addition to $`\sigma _b`$ and $`\sigma _c`$, determines this accuracy. As noted earlier, $`N_{images}`$ was of the order of 15. If a drop has an exceptionally large radius or is falling too far from the centerline of the airflow tube, $`v_z`$ will be too large and $`N_{images}`$ will be too small.

#### IV.1.3 $`\chi ^2`$ criterion

As mentioned earlier, the final calculation of the charge on a drop is done using a fit to the trajectory of the drop. It was required that the $`\chi ^2`$ probability of the fit to the drop’s trajectory be better than $`10^3`$. This rejects a large class of rare artifacts based on the statistical likelihood that the observed deviations from the fitted trajectory could be attributed to the Brownian motion and centroiding errors. For example, a drop would be rejected if it had an anomalous trajectory due to vibrations in the apparatus or due to its charge having been changed by collision with an ion during measurement.

#### IV.1.4 $`v_z`$ criterion

The net downward velocity of the drop, $`v_z`$, depends upon the drop radius and the upward air velocity, $`v_{air}`$, Eq. (3). This criterion insures consistent drop radii within the hour long data blocks by requiring

$$|v_{z,drop}v_{z,block}|<0.124\text{ mm/s}$$
(12)
where $`v_z`$ is the measured value for one drop and $`v_{z,block}`$ is the average value of $`v_z`$ for all the drops in the one hour data block. Using Eqs. 2 and 3, taking $`v_{air}`$ as fixed and using an average value for $`r`$ of 10 $`\mu `$m, this eliminates any drops with $`r`$ different from the nominal value by more than about $`\pm 0.5\%`$.

Recall that the air velocity is approximately parabolic and that close to the centerline it is given by

$$v_{air}(x)=v_{air,0}\left[1(x/x_w)^2\right]$$
(13)
where $`x`$ is the distance along the $`x`$ axis from the centerline of the airflow tube, $`x_w=4.15`$ mm is the distance to the wall of the tube, and $`v_{air,0}`$ is the air velocity along the centerline. Therefore this criterion indirectly restricts how far the drop can be from the centerline.

#### IV.1.5 $`x`$ deviation criterion

This criterion

$$|xx_{block}|<0.19\text{ mm}$$
(14)
provides a direct constraint on how far a drop may deviate from the centerline in the $`x`$ direction. Here $`x_{block}`$ is the average value of $`x`$ for all the drops in the one hour data block. The purpose of this criterion is to eliminate drops that were produced irregularly. The 0.19 mm upper limit in Eq. (14) was determined by examining the distribution of $`x`$ positions of drops produced during normal operation of the drop generation and setting the upper limit to eliminate the tails.

#### IV.1.6 Minimum distance $`R`$ between any two drops criterion

The drops interact with one another through their induced electric dipole moments and viscous coupling through the air. Consider two drops of radius $`r`$, drop 1 moving with a velocity $`v_x`$ due to the force of the electric field on its charge, as shown in Fig. 6. This motion will move the surrounding air. At the position of drop 2, the velocity of the air in the $`x`$ direction, $`V_{x,\text{disturbed air}}`$, is given by

$$V_{x,\text{disturbed air}}=\frac{3}{4}\frac{v_xr}{R}(1+\text{cos}^2\theta ).$$
(15)
Since drop 2 sits in this disturbed air, its $`v_x`$ due to the force of the electric field on its charge will have superimposed upon it $`V_{x,\text{disturbed air}}`$. This will distort the charge measurement. Therefore $`V_{x,\text{disturbed air}}`$ must be kept small by keeping $`R`$, the distance between the drops, much larger than $`r`$, the radius of the drop. A large separation also serves to minimize the interaction between the induced electric dipole moments of the drops, which increases as the inverse fourth power of the separation. We require

$$R>0.62\text{ mm}$$
(16)
separation between any two drops, which limits these forces to a small fraction of $`QE`$.

#### IV.1.7 Summary and magnitude of drop selection criteria

Table 2 gives the percent of drops removed by each criterion averaged over each of the two data sets. The total percent of drops removed is also given. Since the same drop may be removed by several criteria, the total percent removed is *not* the sum of the percent removed by the individual criteria.

### IV.2 Results

After the application of these criteria we had a final data sample of $`1.7\times 10^7`$ drops of average diameter 20.6 $`\mu `$m. The total mass of the sample was 70.1 mg. Figure 7 shows the charge distribution in units of $`e`$. (The asymmetry of the charge distribution is a property of the drop generator as discussed in Sec. II.3.) We see sharp peaks at integer numbers of charges and no drops further than 0.15 $`e`$ from the nearest integer. We emphasize that there is no background subtraction here, this is all the data after the application of the criteria previously discussed.

To show the shape of the peaks at integer values of $`q`$ we superimpose them in Fig. 8 using the charge distribution, $`q_c`$, defined by $`q_c=qN_c`$ where $`N_c`$ is the signed integer closest to $`q`$. The peaks have a Gaussian distribution with a standard deviation of 0.021 $`e`$. *The absence of non-Gaussian tails is what allows this search method to be so powerful.*

In Fig. 9 we superimpose the valleys between the peaks using the residual charge distribution, $`q_r`$, defined by $`q_r=qN_l`$ where $`N_l`$ is the largest integer less than $`q`$. We did not find any drops with residual charge between 0.15 $`e`$ and 0.85 $`e`$. In this residual charge range there are fewer than $`1.17\times 10^{22}`$ fractional charge particles per nucleon in silicone oil with 95% confidence.

Thus this 70.1 mg search did not confirm the one unusual aspect of our previous 17 mg search, Halyo et al.halyo , where we found 1 drop with a $`q_r`$ of about 0.29 $`e`$. No such charge was found in this search. While it is of course still possible that the fractional charge found in the 17 mg experiment was real, we are inclined to believe that the 17 mg experiment had a very small background that has been eliminated by the improved method of this experiment.

## V CONCLUSIONS

### V.1 Comparison with other fractional charge searches in bulk material

In Table 3 we compare this search with previous, larger sample, searches for fractional charge particles in bulk matter. No evidence for fractional charge particles was found in the searches by Marinelli et al.mari , Smith et al.smi3 , and Jones et al.jon2 , similar to the null result in the present search.

In their superconducting levitometer search in niobium, LaRue et al.lar  claimed to have fractional charge particles with $`e/3`$ and $`2e/3`$. But Smith et al.smi3  who also searched in niobium using a ferromagnetic levitometer method did not find any such evidence in a four times larger sample. At present the results of LaRue et al.lar  are not understood and are generally not accepted.

Our search is by far the largest to date and has the smallest upper limit of any search on the concentration of fractional charge particles in bulk matter. But it is important not to generalize our limit to other kinds of bulk matter for several reasons. First, we do not know what happens to fractional charge particles that are in natural matter when that matter is processed. Note that except for the search in meteoritic material by Jones et al.jon2 , all the material in Table 3 is processed.

Second, if we assume the existence of stable, fractional charge particles, we do not know what natural materials are most likely to have a detectable concentration. Our own thoughts are that the most promising natural material is that found in carbonaceous chondrite asteroids, since they are representative of the primordial composition of the solar system, having not undergone any geochemical or biochemical processes. Hence similar to the motivation of Jones et al.jon2 , our next search will be in meteoritic material from an asteroid.

### V.2 Remarks on further use of this new method

The purpose of the new method loomba  used in this experiment was to allow large drops to be used compared to the classical method, thus increasing the rate at which we could search through a sample and also enabling the use of suspensions of more interesting materials. We have succeeded in doing this, using drops of about 20 $`\mu `$m diameter compared to the approximately 10 $`\mu `$m diameter used in Halyo et al.halyo . In the Appendix we discuss further increasing the search rate by using still larger drops and by using multiple columns of drops to increase the total rate of drop production. We find that with this new method the mass per second search rate can be further increased by a factor of the order of 10, but probably not by a factor of a 100.

*

## Appendix A INCREASING THE SEARCH RATE

Three are several ways in which the mass per second search rate can be increased in this experimental method.

### A.1 Use of larger drops

The first way to increase the search rate is to use larger drops. Maintenance of the precision of the charge measurement requires that $`N_{images}`$ increase in proportion to the drop radius. An increase in $`N_{images}`$ can be accomplished by some combination of a decrease in $`v_z`$ and an increase in the vertical length $`Z`$, Eq. (4). However a significant decrease in $`v_z`$ requires too fine a balance between $`v_{air}`$ and $`r^2`$, Eq. (3). If we keep $`v_z`$ constant, an increase of $`Z`$ can be attained by an increase in the number of vertical direction pixels in the CCD array of the camera. Existing CCD cameras with 10 Hz frame reading rates have twice the 736 vertical pixels used in the present camera and larger arrays will probably be available in the future. Therefore based on this consideration alone, drop diameters of several times 20 $`\mu `$m are feasible.

However, there are two problems that must be considered for drop diameters larger than 30 to 40 $`\mu `$m. The dipole force on a drop in a non-uniform electric field is proportional to the third power of the drop diameter. This force was negligible in this experiment, Sec. III.3, but would have to be considered for much larger drops. The other problem is that the maintenance of a small and constant $`v_z`$, Eq. (3), requires $`v_{air}`$ to increase as the square of the drop diameter, possibly leading to non-laminar flow. Therefore without more design and experimental studies, our conservative conclusion is that the drop diameter is limited to about 30 $`\mu `$m. This would lead to an increase of the mass search rate by a factor of 3 compared to the 20 $`\mu `$m drops used in this experiment.

### A.2 Increase of drop production rate per column of drops

Let the drop production rate for a column of drops be $`n`$ per second. Then the vertical separation between drops in a column is $`R=v_{air}/n`$. The criteria in Sec. IV.1.6 require $`R>0.62`$ mm. Using $`v_{air}=2.0`$ mm/s, this gives an upper limit on $`n`$ of about 3 Hz. However our experience in this experiment, Sec. II.3, strongly suggests that a maximum 1 Hz rate is conservative practice, because of irregularities in drop production.

### A.3 Increase in the number of drop columns

In this experiment we used one column of drops, however the extension to many columns of drops is straightforward. Of course the horizontal separation between adjacent columns must meet the requirements of Eq. (16), a nominal separation of 1 mm is useful for design purposes. The use of multiple columns requires two changes in the experimental design, namely the number of pixels in the horizontal direction in the CCD array must be increased and the space between the electric field plates must be increased. The latter requirement means the alternating potential applied across the plates must also be increased to keep the electric field constant.

Existing 10 Hz frame rate CCD cameras limit the number of columns to three but improvements in these cameras would probably allow five columns. The corresponding increase in the electric plate spacing and the potential difference is straightforward.

### A.4 Correction for drop to drop interactions

It is clear that the primary constraint limiting the density of drops achievable in the measurement chamber is that of Eq. (16). To reiterate, interaction between the drops due to their induced electric dipole moment and viscous coupling requires that there be a minimum separation allowable between drops. In the limit that these interactions are small, both of these effects can be calculated from first principles, for example as in Eq. (15). In principle then, it should be possible to subtract the effect of these perturbing forces from the measured trajectory of each drop. Given this, it would be possible to relax the constraint on $`R`$. Since this possibility requires further study, it is not clear to what extent $`R`$ can be reduced and throughput increased.

### A.5 Summary

Putting these estimates together we can see how to achieve an improvement on the order of 10 times the present mass per second search rate using existing CCD cameras. Future cameras will probably allow a factor of 15 improvement.


text:
# 1 𝑓' vs. 𝛽 at fixed external field strength (𝑛_"ext"=1) for 𝐿_𝑡={4,5,6,7,8} and 𝐿_𝑠=64. Solid lines are the fits Eq. ().

BARI-TH 436/2002

Abelian chromomagnetic fields and confinement

Paolo Cea1,2,***Paolo.Cea@ba.infn.it and Leonardo Cosmai1,Leonardo.Cosmai@ba.infn.it

1INFN - Sezione di Bari, I-70126 Bari, Italy
2Physics Department, Univ. of Bari, I-70126 Bari, Italy

April, 2002

## Abstract

We study vacuum dynamics of SU(3) lattice gauge theory at finite temperature using the lattice Schrödinger functional. The SU(3) vacuum is probed by means of an external constant Abelian chromomagnetic field. We find that by increasing the strength of the applied external field the deconfinement temperature decreases towards zero. This means that strong enough Abelian chromomagnetic fields destroy confinement of color. We discuss some consequences of this phenomenon on confinement and quark stars.

Understanding the mechanism of quark confinement is still a central problem in high energy physics. According to a model conjectured long time ago by G. ’t Hooft  and S. Mandelstam  the confining vacuum behaves as a coherent state of color magnetic monopoles, or, equivalently, the vacuum resembles a magnetic (dual) superconductor. Up to now there is numerical evidence  in favor of chromoelectric flux tubes in pure lattice gauge vacuum. As well there have been extensive numerical studies  of monopole condensation.

To investigate vacuum structure of lattice gauge theories we introduced  a gauge invariant effective action, defined by means of the lattice Schrödinger functional

$$𝒵[U_k^{\mathrm{ext}}]=𝒟Ue^{S_W}.$$
(1)
$`S_W`$ is the Wilson action and the functional integration is extended over links on a lattice $`L_s^3\times L_t`$ with hypertorus geometry and satisfying the constraints ($`x_t`$: temporal coordinate)

$$U_k(x)|_{x_t=0}=U_k^{\mathrm{ext}}(\stackrel{}{x}),(k=1,2,3),$$
(2)
$`U_k^{\mathrm{ext}}(x)`$ being the lattice version of the external continuum gauge field $`\stackrel{}{A}^{\mathrm{ext}}(x)=\stackrel{}{A}_a^{\mathrm{ext}}(x)\lambda _a/2`$. We also impose that links at the spatial boundaries are fixed according to Eq. (2). In the continuum this last condition amounts to the requirement that fluctuations over the background field vanish at infinity.

In terms of the above defined lattice Schrödinger functional the lattice effective action for external static background field $`\stackrel{}{A}^{\mathrm{ext}}(x)`$ is given by:

$$\mathrm{\Gamma }[\stackrel{}{A}^{\mathrm{ext}}]=\frac{1}{L_t}\mathrm{ln}\left\{\frac{𝒵[\stackrel{}{A}^{\mathrm{ext}}]}{𝒵[0]}\right\}.$$
(3)
Indeed, in the continuum limit $`\mathrm{\Gamma }[\stackrel{}{A}^{\mathrm{ext}}]`$ reduces to the vacuum energy in presence of the background field $`\stackrel{}{A}^{\mathrm{ext}}`$. Moreover, the lattice background effective action Eq. (3) is invariant for gauge transformations of the external links $`U_k^{\mathrm{ext}}`$. Different approaches to put background fields on the lattice use external currents  or modified boundary conditions .

At finite temperature we are interested in the thermal partition function  in presence of a static background field. On the lattice we introduced  the background field thermal partition function as

$$𝒵_T\left[\stackrel{}{A}^{\text{ext}}\right]=_{U_k(L_t,\stackrel{}{x})=U_k(0,\stackrel{}{x})=U_k^{\text{ext}}(\stackrel{}{x})}𝒟Ue^{S_W},$$
(4)
where the physical temperature is given by $`T=1/aL_t`$. On a lattice with finite spatial extension we also impose that the links at the spatial boundaries are fixed according to boundary conditions Eq. (2). Thus we see that, sending the physical temperature to zero, the thermal functional Eq. (4) reduces to zero-temperature Schrödinger functional Eq. (1).

At finite temperature the relevant quantity turns out to be the free energy functional which is naturally defined as:

$$F[\stackrel{}{A}^{\mathrm{ext}}]=\frac{1}{L_t}\mathrm{ln}\frac{𝒵_T[\stackrel{}{A}^{\mathrm{ext}}]}{𝒵_T[0]}.$$
(5)
Obviously, if the physical temperature is sent to zero the free energy functional reduces to the vacuum energy functional Eq. (3).

In our previous studies  we found that for U(1) the confining vacuum behaves as a coherent condensate of Dirac magnetic monopoles. On the other hand in SU(2) and SU(3) it seems that there is condensation of Abelian magnetic monopoles and Abelian vortices. So that in SU(2) and SU(3) gauge theories one could look at the confining vacuum as a coherent Abelian magnetic condensate. Moreover, we also found  that a constant Abelian chromomagnetic field at zero temperature is completely screened in the continuum limit, while at finite temperature  it seems that the applied field is restored by increasing the temperature. These results strongly suggest that the confinement dynamics is intimately related to Abelian chromomagnetic gauge configurations. From previous considerations we are led to investigate if deconfinement temperature depends on the strength of an applied external constant Abelian chromomagnetic field. To address this question we performed lattice simulations of pure SU(3) gauge theory at finite temperature and compute the free energy functional Eq. (5) for a static constant Abelian chromomagnetic field. In the continuum a constant Abelian chromomagnetic field is given by:

$$\stackrel{}{A}_a^{\mathrm{ext}}(\stackrel{}{x})=\stackrel{}{A}^{\mathrm{ext}}(\stackrel{}{x})\delta _{a,3},A_k^{\mathrm{ext}}(\stackrel{}{x})=\delta _{k,2}x_1H.$$
(6)
In SU(3) lattice gauge theory we get:

$$\begin{array}{cc}& U_1^{\mathrm{ext}}(\stackrel{}{x})=U_3^{\mathrm{ext}}(\stackrel{}{x})=\mathrm{𝟏},\hfill \\ & U_2^{\mathrm{ext}}(\stackrel{}{x})=\left[\begin{array}{ccc}\mathrm{exp}(i\frac{gHx_1}{2})& 0& 0\\ 0& \mathrm{exp}(i\frac{gHx_1}{2})& 0\\ 0& 0& 1\end{array}\right]\hfill \end{array}$$
(7)
To be consistent with the hypertorus geometry we impose that the magnetic field is quantized as:

$$\frac{gH}{2}=\frac{2\pi }{L_1}n_{\mathrm{ext}},n_{\mathrm{ext}}\text{integer}.$$
(8)
It is easy to verify that the lattice links Eq. (7) give rise to a constant field strength. Therefore, since the free energy functional $`F[\stackrel{}{A}^{\mathrm{ext}}]`$ is invariant for time independent gauge transformations of the background field $`\stackrel{}{A}^{\mathrm{ext}}`$, it follows that $`F[\stackrel{}{A}^{\mathrm{ext}}]`$ is proportional to spatial volume $`V=L_s^3`$ and the relevant quantity is the density $`f[\stackrel{}{A}^{\mathrm{ext}}]`$ of free energy

$$f[\stackrel{}{A}^{\mathrm{ext}}]=\frac{1}{V}F[\stackrel{}{A}^{\mathrm{ext}}].$$
(9)
We want to compute the density of free energy on the lattice. In order to circumvent the problem of computing a partition function which is the exponential of an extensive quantity, we consider the $`\beta `$-derivative of $`f[\stackrel{}{A}^{\mathrm{ext}}]`$ (at fixed external field strength $`gH`$). Indeed this quantity is easy to evaluate numerically since it is related to the plaquette ($`\mathrm{\Omega }=L_s^3\times L_t`$):

$$\begin{array}{cc}\hfill f^{}[\stackrel{}{A}^{\mathrm{ext}}]& =\frac{1}{\mathrm{\Omega }}\underset{x,\mu <\nu }{}\frac{1}{3}\text{Re}\text{Tr}U_{\mu \nu }(x)_0\hfill \\ & \frac{1}{\mathrm{\Omega }}\underset{x,\mu <\nu }{}\frac{1}{3}\text{Re}\text{Tr}U_{\mu \nu }(x)_{\stackrel{}{A}^{\mathrm{ext}}},\hfill \end{array}$$
(10)
where the subscripts on the averages indicate the value of the external field. The generic plaquette $`U_{\mu \nu }(x)=U_\mu (x)U_\nu (x+\widehat{\mu })U_\mu ^{}(x+\widehat{\nu })U_\nu ^{}(x)`$ contributes to the sum in Eq. (10) if the link $`U_\mu (x)`$ is a ”dynamical” one (i.e. it is not constrained in the functional integration Eq. (4)).
Observing that $`f[\stackrel{}{A}^{\mathrm{ext}}]=0`$ at $`\beta =0`$, we may obtain $`f[\stackrel{}{A}^{\mathrm{ext}}]`$ from $`f^{}[\stackrel{}{A}^{\mathrm{ext}}]`$ by numerical integration:

$$f[\stackrel{}{A}^{\mathrm{ext}}]=_0^\beta f^{}[\stackrel{}{A}^{\mathrm{ext}}]𝑑\beta ^{}.$$
(11)
As is well known pure SU(3) gauge system undergoes a deconfinement phase transition by increasing temperature. As a matter of fact, it turns out that the knowledge of $`f^{}[\stackrel{}{A}^{\mathrm{ext}}]`$ at finite temperature can be used to estimate the critical temperature $`T_c`$. To this end it suffices to evaluate $`f^{}[\stackrel{}{A}^{\mathrm{ext}}]`$ as a function of $`\beta `$ for different lattice temporal sizes $`L_t`$. Indeed we found that $`f^{}[\stackrel{}{A}^{\mathrm{ext}}]`$ displays a peak in the critical region. In Figure 1 we display the peak regions for different values of $`L_t`$. We see clearly that the pseudocritical coupling, which is the value of $`\beta `$ at the peak, depends on $`L_t`$. In order to determine the pseudocritical couplings $`\beta ^{}(L_t)`$ we parameterize $`f^{}(\beta ,L_t)`$ near the peak as:

$$\frac{f^{}(\beta ,L_t)}{\epsilon _{\mathrm{ext}}^{}}=\frac{a_1(L_t)}{a_2(L_t)[\beta \beta ^{}(L_t)]^2+1}.$$
(12)
In Eq. (12) we normalize $`f^{}`$ to $`\epsilon _{\mathrm{ext}}^{}`$, the derivative of the classical energy due to the external applied field

$$\epsilon _{\mathrm{ext}}^{}=\frac{2}{3}[1\mathrm{cos}(\frac{gH}{2})]=\frac{2}{3}[1\mathrm{cos}(\frac{2\pi }{L_1}n_{\mathrm{ext}})].$$
(13)
We restrict the region near $`\beta ^{}(L_t)`$ until the fits Eq. (12) give a reduced $`\chi ^2/\text{ d.o.f.}1`$. Once determined $`\beta ^{}(L_t)`$ we estimate the deconfinement temperature as

$$\frac{T_c}{\mathrm{\Lambda }_{\mathrm{latt}}}=\frac{1}{L_t}\frac{1}{f_{SU(3)}(\beta ^{}(L_t))},$$
(14)
where

$$f_{SU(N)}(\beta )=\left(\frac{\beta }{2Nb_0}\right)^{b_1/2b_0^2}\mathrm{exp}\left(\beta \frac{1}{4Nb_0}\right)$$
(15)
$`N`$ being the color number, $`b_0=(11N)/(48\pi ^2)`$, and $`b_1=(34N^2)/(3(16\pi ^2)^2)`$.

In order to obtain the continuum limit critical temperature we have to extrapolate $`T_c/\mathrm{\Lambda }_{\mathrm{latt}}`$, given by Eq. (14), to the continuum.

In Figure 2 we display $`T_c/\mathrm{\Lambda }_{\mathrm{latt}}`$ as a function of $`aT_c`$ for two different values of the chromomagnetic background field. To extract the continuum limit of the critical temperature we follow Ref.  and perform a linear extrapolation to the continuum of our data for $`T_c/\mathrm{\Lambda }_{\mathrm{latt}}`$. For comparison we also display in Fig. 2 the continuum limit of the critical temperature without background field

$$\frac{T_c}{\mathrm{\Lambda }_{\mathrm{latt}}}=\mathrm{\hspace{0.17em}29.67}\pm \mathrm{\hspace{0.17em}5.47}.$$
(16)
We see clearly that the continuum limit deconfinement critical temperature does depend on the applied Abelian chromomagnetic field. Therefore we decided to vary the strength of the applied external Abelian chromomagnetic background field to study quantitatively the dependence of $`T_c`$ on $`gH`$ .

To this aim we performed numerical simulations on $`64^3\times L_t`$ lattices with $`n_{\text{ext}}=1,2,3,5,10`$. Since we are measuring a local quantity, such as average plaquette, a low statistics (from 1000 up to 5000 configurations) is required in order to get a good estimate of $`f^{}[\stackrel{}{A}^{\mathrm{ext}}]`$. Statistical errors have been estimated using the jackknife resampling method modified to take into account the statistical correlations between lattice configurations. As a check of possible finite volume effects we also simulated our gauge system on $`128\times 64^2\times L_t`$ lattice with $`n_{\text{ext}}=3`$. Indeed, we found that our numerical data do not display appreciable finite volume effects.

Following previous steps we are able to determine the critical temperature as a function of the external chromomagnetic field. In Figure 3 we display our determination of $`T_c`$ for different values of the applied field strength. We see that the critical temperature decreases by increasing the external Abelian chromomagnetic field. If the magnetic length $`a_H1/\sqrt{gH}`$ is the only relevant scale of the problem, then for dimensional reasons one expects that

$$T_c^2gH.$$
(17)
Indeed we try to fit our data with

$$\frac{T_c(gH)}{\mathrm{\Lambda }_{\text{latt}}}=\frac{T_c(0)}{\mathrm{\Lambda }_{\text{latt}}}+\alpha \sqrt{gH}.$$
(18)
We find a satisfying fit (see Fig. 3) with $`\alpha =42.4\pm 7.4`$ and $`T_c(0)/\mathrm{\Lambda }_{\text{latt}}`$ that agrees within errors with Eq. (16). Remarkably, we see that there exists a critical field

$$gH_c0.68$$
(19)
such that $`T_c=0`$ for $`gH>gH_c`$.

Let us conclude by stressing the main results of this paper. We found that there is a critical field $`gH_c`$ such that, for $`gH>gH_c`$, the gauge system is in the deconfined phase. As a consequence, we see that there is an intimate connection between Abelian chromomagnetic fields and color confinement. The existence of a critical chromomagnetic field is not easily understandable within the coherent magnetic monopole condensate picture of the confining vacuum. On the other hand, it could be explained if the vacuum behaves as an ordinary relativistic color superconductor. Thus we have to reconcile two apparently different aspects. From one hand, the confining vacuum does display condensation of both Abelian magnetic monopoles and vortices, on the other hand the relation between the deconfinement temperature and the applied Abelian chromomagnetic field consistent with Eq. (17) suggests that the vacuum behaves as a condensate of a color charged scalar field whose mass is proportional to the inverse of the magnetic length. A natural candidate for such tachyonic color charged scalar field is the famous Nielsen-Olesen unstable mode . The eventual condensation of the Nielsen-Olesen modes makes the vacuum a dynamical color superconductor. However, it should be stressed that, unlike the case of an elementary color charged scalar field, the chromomagnetic condensate cannot be uniform due to gauge invariance of the vacuum. Indeed, as showed  by R. P. Feynman (in 2+1 dimensions) gauge invariance of the ground state disorders the system in such a way that there are not long range color correlations. In this way the disordered chromomagnetic condensate should lead to the screening of any external Abelian chromomagnetic charge, thus explaining the condensation of both Abelian chromomagnetic monopoles and vortices.

Another important aspect of this work might be related to astrophysics applications. Indeed, aside from the above considerations, our result implies that it is possible to deconfine quarks and gluons even at low temperature and low density. As a consequence follows the exciting possibility of low density compact quark stars. Obviously, further progress on this matter needs the equation of state appropriate for the description of deconfined quarks and gluons in an Abelian chromomagnetic field.


text:
# Polarons with a twist

## Abstract

We consider a polaron model where molecular *rotations* are important. Here, the usual hopping between neighboring sites is affected directly by the electron-phonon interaction via a twist-dependent hopping amplitude. This model may be of relevance for electronic transport in complex molecules and polymers with torsional degrees of freedom, such as DNA, as well as in molecular electronics experiments where molecular twist motion is significant. We use a tight-binding representation and find that very different polaronic properties are already exhibited by a two-site model – these are due to the nonlinearity of the restoring force of the twist excitations, and of the electron-phonon interaction in the model. In the adiabatic regime, where electrons move in a low-frequency field of twisting-phonons, the effective splitting of the energy levels increases with coupling strength. The bandwidth in a long chain shows a power-law suppression with coupling, unlike the typical exponential dependence due to linear phonons.

polarons, phonons, molecular electronics

There has been considerable interest in the study of electronic transport in organic materials. They provide technologically relevant alternatives for new electronic and optical materials, and offer an opportunity to gain deeper understanding of quantum transport processes. Moreover, transport in molecules also constitutes the basis for important biological phenomena, such as photosynthesis and other charge transfer events. book1  Recent developments in the fabrication of high electron mobility systems using organic semiconductors, have further emphasized the need to understand mobile polaronic behavior in complex molecular crystals. schon  In addition, a series of beautiful experiments in the area of molecular electronics have underscored the importance of atomic vibrations in the measured current-voltage characteristics of these systems. Reed ; Reifenberger ; Zhitenev ; McEuen

In typical inorganic materials the atomic degrees of freedom (DOF) associated with translation are the most important in electronic transport. In contrast, the twisting DOF (or “librons”) of molecular complexes may also play important roles in organic and molecular solids. In a wide variety of organic linear-chain crystals of stacked molecules, the librons may in fact be essential. For instance, it was found that libron scattering was important in the understanding of conductivity mechanisms in the TTF-TCNQ system, ttf1 ; ttf2 ; ttf3 ; ttf4  as well as in phthalocyanine-based “molecular metals”. hale  More recently, the conductivity of DNA and constituent base crystals has been studied experimentally, expa ; expb ; expc ; expd ; expe  and theoretically. th1 ; yu ; bru ; la ; C00  Since DNA can be viewed as a “molecular wire”, a one-dimensional twisted chain of stacked base pairs with a somewhat flexible structure, B00  it is likely that parallel rotation of adjacent base pairs will be important in electronic transfer mechanisms. Notice that these modes would have a stronger impact on the electronic overlaps along the chain than other flexural modes. In fact, a variable range hopping theory where the twist angle fluctuations between base pairs are taken into account was used recently to explain the temperature dependence of electrical conductivity along the DNA double helix. yu  Bruinsma et al. have also studied charge migration along DNA based on a model in which large structure fluctuations, including relative rotation between adjacent bases, are considered. bru

There exist a few models of electron-libron interactions in the literature. Many of them have considered the situation of small rotation angles in a harmonic or nearly harmonic approximation. ttf3 ; ttf4 ; hale ; bru ; B00  Although other authors have considered large angle rotations, yu ; bar  the molecular kinetic energy term was neglected, and thus the inertial backreaction effects in molecules were not taken into account. Our model considers the possibility of large rotation effects, which are not negligible in systems such as (TTF)7I5, and others with similarly flexible structure. We also include the full kinetics of the molecule in order to include backreaction effects. Our description is fully quantum mechanical and includes strong anharmonic interactions, shown to be essential for polaron self-trapping. RBK  From the full solution of a two-site model, we find that in the adiabatic regime, where electrons move in a low frequency field of twisting-phonons, the splitting of the energy levels (and therefore the effective intersite hopping) increases with phonon coupling strength. In the opposite, non-adiabatic regime, the effective hopping constant increases with phonon frequency, but saturates at high values. For a long chain model with local-twists, we find a power-law suppression of the bandwidth, unlike the exponential result for translational polarons. These results arise from the different coupling scheme and the nonlinear twist-polaron restoring force, as we now describe.

The twist polaron model. The on-site interaction of the electron with phonon modes tends to localize the electron. On the other hand, the molecular transfer (hopping) integral in a given atomic or molecular configuration results in electron delocalization. In an extended system, the hopping integral competes with the phonon interaction, and depending on the relative strengths of the different terms, one obtains a mobile polaron (with a heavier effective mass), or a self-localized polaron. fro ; holst  Our model allows for such physical behavior as well. Most importantly, however, the phonon enters here modifying also the hopping integral and is then substantially different to typical models for translation-phonon polarons. pol-book

For simplicity, let us first consider a two-site model. The dynamics of the system is described by the Hamiltonian

$`H`$ $`=`$ $`(c_1^+c_1+c_2^+c_2)\epsilon ^{}\mathrm{cos}\theta `$
$`(t_0+t_1\mathrm{cos}\theta )(c_1^+c_2+c_2^+c_1)`$
$`{\displaystyle \frac{\mathrm{}^2}{2I}}{\displaystyle \frac{^2}{\theta ^2}}K\mathrm{cos}\theta ,`$
where the first term represents the possible on-site energy modulation due to the molecular rotations with angular DOF $`\theta `$, and the $`c_i^+`$ operators describe the two electronic wave functions (“sites”). The second term is the rotation-dependent hopping term. Notice that the molecular rotation would either enhance or suppress electron hopping in this model, depending on whether $`t_1`$ and $`t_0`$ have the same or opposite sign. A microscopic calculation of stacked molecules, such as the base pairs in DNA, or the phthalocyanine stacks, allows one to estimate the effect of molecular twists on the electron hopping or overlap. B00  Recent results for guanosine crystals, moli  for example, show clearly that eclipsing of molecules in successive “layers” is essential for substantial electron hopping between them. This twist dependence is introduced in our model via the hopping amplitude $`t_1`$, and then plays the role of a coupling constant. The last two terms in (1) are the relative kinetic energy for the molecules and the interaction potential characterized by the libron restoring ‘force constant’ $`K`$. Here, $`\theta =\theta _1\theta _2`$, where $`\theta _i`$ is the rotation angle of each molecule, and $`I`$ is the reduced moment of inertia for the relative rotation of the two adjacent molecules (see inset Fig. 1). An overall rotation of the pair is gapless and neglected. Notice also that the first term in (1) leads to a constant shift of $`K`$, which is then implicit in the following discussion. Notice that this model considers nonlinear coupling to the electrons via $`t_1`$, RBK  as well as a non-linear restoring force for the relative rotation of this two-molecule stack. On the other hand, the model does not consider the possible mixing of rotation DOF with translation displacements, although this coupling is for the most part small. bru  We believe that our model contains the main physics of polaronic excitations corresponding to the rotation DOF.

In the following, we use the energy unit $`\mathrm{}^2/2I`$ for convenience, so that all energy quantities are dimensionless from now on. In order to solve the eigenvalue problem $`H\psi =E\psi `$, we use

$$\begin{array}{c}|\psi =\delta (\theta )c_1^+|0+\beta (\theta )c_2^+|0,\hfill \end{array}$$
(2)
which yields, from $`(H\epsilon 1)\psi =0`$, a set of secular equations for $`\delta `$ and $`\beta `$. Defining $`A=\delta +\beta `$, and $`B=\delta \beta `$, we have

$$\begin{array}{c}(\frac{^2}{\theta ^2}+\epsilon +t_0+p\mathrm{cos}\theta )A=0\hfill \\ (\frac{^2}{\theta ^2}+\epsilon t_0+q\mathrm{cos}\theta )B=0,\hfill \end{array}$$
(3)
with $`q=Kt_1`$, and $`p=K+t_1`$. The solutions are periodic Mathieu functions, math  with properties which depend on the values of the coefficients $`p`$ and $`q`$.

As we will see below, the rotation DOF in the model produce modifications of the electronic hopping. In order to reveal the physics more clearly, we analyze the two-site problem in more detail in the adiabatic and non-adiabatic regimes.

Non-adiabatic regime. In the non-adiabatic regime $`Kt_0,t_1`$, the electrons move in a field of high frequency twist-phonons. We first neglect fully the hopping terms. We find two degenerate states

$$\begin{array}{c}\psi _1=\left(\begin{array}{c}\varphi _0\\ 0\end{array}\right),\psi _2=\left(\begin{array}{c}0\\ \varphi _0\end{array}\right),\hfill \end{array}$$
(4)
where $`\varphi _0`$ is the normalized ground state (with eigenvalue $`\epsilon _0`$), solution of the equation

$$\begin{array}{c}\frac{d^2\varphi }{d\theta ^2}+[\epsilon +K\mathrm{cos}\theta ]\varphi (\theta )=0.\hfill \end{array}$$
(5)
To calculate the first order correction, we take the wave function $`\psi =\alpha \psi _1+\beta \psi _2`$, so that the secular equation becomes

$$\begin{array}{c}\left(\begin{array}{cc}\epsilon _0\epsilon & (t_0+t_1^{})\\ (t_0+t_1^{})& \epsilon _0\epsilon \end{array}\right)\left(\begin{array}{c}\alpha \\ \beta \end{array}\right)=0,\hfill \end{array}$$
(6)
where the renormalized twist-dependent hopping constant is

$$t_1^{}=t_1\varphi _0^{}\mathrm{cos}\theta \varphi _0d\theta .$$
(7)
The eigenvalues are simply $`\epsilon _{1,2}=\epsilon _0\pm |t_0+t_1^{}|`$, so that the splitting of the two originally degenerate eigenvalues is $`\mathrm{\Delta }\epsilon =2|t_0+t_1^{}|`$. It is easy to see from Eq. (7) that $`t_1^{}t_1`$. Therefore, $`\mathrm{\Delta }\epsilon T`$, where $`T=2t_0+2t_1`$ is the hopping constant “without rotation”, i.e. $`\theta 0`$.

In the regime of large $`K`$ ($`K1`$), the local twist potential is effectively harmonic and one can write the solution as a superposition of local harmonic oscillators about $`2m\pi `$,

$$\varphi _0=\underset{m}{}\frac{\alpha ^{1/2}}{\pi ^{1/4}}e^{\alpha ^2(\theta 2m\pi )^2/2},$$
(8)
where $`\alpha ^4=K/2`$, and the eigenvalue is $`\epsilon _0=K+\alpha ^2`$. Correspondingly,

$$t_1^{}=t_1e^{1/4\alpha ^2}t_1\left(1\frac{1}{4\alpha ^2}\right),$$
(9)
and $`\mathrm{\Delta }\epsilon =2t_0+2t_1(11/4\alpha ^2)T`$, with $`\mathrm{\Delta }\epsilon `$ approaching $`T`$ in the large $`K`$ limit. In this case, the molecule is weakly oscillating about $`\theta =0`$ (mod $`2\pi `$), with a relatively minor impact on the electron hopping. Notice that the resulting “bandwidth narrowing” is similar in form to the result of a two-site vibrational polaron. pol-book  There, the effective bandwidth in the non-adiabatic regime is $`\stackrel{~}{J}=Je^{\gamma ^2/2}`$, where $`J`$ is the original bandwidth and $`\gamma `$ is the electron-phonon coupling constant. Notice however the different role that $`K`$ plays here, and that the bandwidth is only weakly suppressed in this case.

When $`K`$ is small ($`1Kt_0,t_1`$), the ground state function can be written as, math

$$\varphi _0(\theta )\frac{1}{\sqrt{2\pi (1+K^2/2)}}(1+K\mathrm{cos}\theta ),$$
(10)
and the eigenvalue is $`\epsilon _0=\frac{1}{2}K^2`$. We find that $`t_1^{}=t_1\frac{K}{1+K^2/2}Kt_1`$, and $`\mathrm{\Delta }\epsilon =2t_0+2Kt_1T`$. Notice that in this limit of soft restoring force, the molecules are nearly uniformly rotating at all angles, and strongly suppress the electron intersite hopping.

We see that with increasing $`K`$, the renormalized twist-dependent hopping constant $`t_1^{}`$ increases and saturates to $`t_1`$, as shown schematically in Fig. 1. We should emphasize that one can view the resulting $`\mathrm{\Delta }\epsilon `$ as akin to a bandwidth in the long chain limit. These results then indicate that the twisting DOF strongly prevent the electron from propagating along the chain. The “stiffer” the rotation mode or libron, the easier for the electrons to move from site to site.

Adiabatic regime. In this limit, $`Kt_0,t_1`$, the electron moves ‘rapidly’ in a field of low frequency ‘slow’ phonons. We assume a wave function of the form

$$\begin{array}{c}\left(\begin{array}{c}\alpha (\theta )\\ \beta (\theta )\end{array}\right)=S(\theta )\left(\begin{array}{c}s_1(\theta )\\ s_2(\theta )\end{array}\right),\hfill \end{array}$$
(11)
and first consider the electron moving in a frozen background,

$$\begin{array}{c}\left(\begin{array}{cc}E_0(\theta )& (t_0+t_1\mathrm{cos}\theta )\\ (t_0+t_1\mathrm{cos}\theta )& E_0(\theta )\end{array}\right)\left(\begin{array}{c}s_1(\theta )\\ s_2(\theta )\end{array}\right)=0.\hfill \end{array}$$
(12)
The corresponding eigenvalues $`E_0^\pm =\pm |t_0+t_1\mathrm{cos}\theta |`$ play the role of a potential energy term in the equations for $`S(\theta )`$,

$$\left\{\frac{^2}{\theta ^2}+K\mathrm{cos}\theta E_0^\pm (\theta )+\epsilon \right\}S(\theta )=0.$$
(13)
These equations can also be solved in terms of Mathieu functions. It is easy to find that when $`Kt_0,t_1`$, and $`t_11`$, the eigenvalues of the problem are given by

$`\epsilon _1`$ $`=`$ $`\epsilon _{op}t_0K+{\displaystyle \frac{K}{\sqrt{8t_1}}}`$ (14)
$`\epsilon _2`$ $`=`$ $`\epsilon _{op}+t_0+K{\displaystyle \frac{K}{\sqrt{8t_1}}},`$
where $`\epsilon _{op}=t_1+\sqrt{t_1/2}`$ is the polaronic energy shift due to the electron-libron off-site interaction. In this case, $`\mathrm{\Delta }\epsilon =2t_0+2KK/\sqrt{2t_1}`$. Unlike the non-adiabatic case, $`\mathrm{\Delta }\epsilon `$ saturates to $`2t_0+2K`$ for increasing $`t_1`$, and it therefore shows a strong bandwidth narrowing. Notice, however, that the two-site adiabatic problem for the small Holstein polaron yields an even stronger suppression of tunneling. In that case, pol-book  one finds an exponential narrowing with coupling, $`e^{\gamma ^2/2}`$, similar to the non-adiabatic case mentioned above, but with large $`\gamma `$ values.

Long chain rotational polaron. The long-chain generalization of the two-site model is an interesting physical problem with additional mathematical complexity which will be reported elsewhere. toberep  Here, we focus on a simpler local-twist model involving local angular coordinates. The local-twist hamiltonian is given by

$$H=\underset{i}{}\left\{\frac{\mathrm{}^2}{2I}\frac{^2}{\theta _i^2}K\mathrm{cos}\theta _i\epsilon _1c_i^+c_i\mathrm{cos}\theta _it_0(c_i^+c_{i+1}+c_{i+1}^+c_i)\right\},$$
(15)
which includes anharmonic restoring force and local non-linear electron-phonon coupling via $`\epsilon _1`$ (the more general $`t_1`$-dependence is also discussed in Ref. toberep, ).

In the non-adiabatic regime, $`K1`$, one can generalize the two-site approach and identify two limiting cases. For strong restoring force, $`K\epsilon _1`$, it is possible to show that the polaron energy shift is $`\epsilon _1`$, while the bandwidth is slightly reduced to $`\mathrm{\Delta }=2t_0(1\epsilon _1/K)`$. In the case of softer twist frequency, $`\epsilon _1K1`$, the bandwidth is given by $`2t_0(K/\epsilon _1)^{1/8}`$. This power law suppression of hopping is much weaker than the usual Holstein ‘translational’ polaron in long chains. pol-book  The different dependence arises mostly from the local anharmonic coupling in this case, and would result in the polaron state having a higher mobility along the chain, despite the strong coupling constant.

Conclusion. We have introduced the idea of librons interacting with electrons in a complex molecular system. The strong coupling of rotational DOF yields a unique class of twisting polarons. Inclusion of anharmonic coupling to the electron via hopping and a nonlinear restoring potential, were shown to result in rather different polaronic behavior. We find electronic bandwidth suppression thanks to the excitation of twisting phonons, just as in the case of the usual vibrational DOF, but with a comparatively weaker effect in our case. This difference is due to both of the nonlinear characteristics mentioned, and provides a possible route for differentiation in experiments. We anticipate that the twist polaronic features could be tested by performing experiments in a variety of different molecules, where the restoring force is changed systematically, for example. We expect that our model will provide some insights on the low-temperature transport properties of interesting systems, such as recent molecular electronics experiments, and helical polymer “soft” crystals. soft  Notice however that a detailed comparison with experiments will require one to consider the coupling of these polarons to charge or current leads. Such effects in vibrational polarons have only recently been studied, N01  and one would anticipate qualitative differences, which will be discussed elsewhere.

We acknowledge support from US DOE grant no. DE-FG02-91ER45334, NSF NIRT grant no. 0103034, and from the Condensed Matter and Surface Sciences Program at Ohio University. Helpful discussions with the NIRT group at OU are also appreciated.