text:
# Discrete dynamical systems with 𝑊(𝐴⁽¹⁾_{𝑚-1}×𝐴⁽¹⁾_{𝑛-1}) symmetry

## 1 Introduction

It is known that continuous and discrete Painlevé equations admit affine Weyl group symmetry. These symmetries are realized as birational Bäcklund transformations and can be described in universal way for general root systems. On the other hand, in recent studies on (ultra-)discrete integrable systems, certain affine Weyl group representation has been obtained. This representation is “tropical” (i.e. given in terms of subtraction-free birational mapping) and has a combinatorial counterpart through the ultra-discretization procedure. The aim of this paper is to reveal the relation between these two affine Weyl group representations. In fact, we show that they are essentially the same.

In section 2, we recall the definition of birational transformations $`\pi ,r_0,\mathrm{},r_m`$ and $`\rho ,s_0,\mathrm{},s_n`$ acting on $`mn`$ variables $`x_{ij}`$ ($`i=1,\mathrm{},m`$, $`j=1,\mathrm{},n`$). The main result (Theorem 1.1) states that these transformations give a representation of affine Weyl group of type $`A_{m1}^{(1)}\times A_{n1}^{(1)}`$. We give the proof of Theorem 1.1 in section 3. In section 4, we describe the discrete dynamical systems arising from our affine Weyl group representation. The systems can be viewed as a version of discrete Toda equation and its generalizations. We also discuss the relation with $`q`$-Painlevé equation studied in .

## 2 Affine Weyl group actions

Let $`m`$, $`n`$ be positive integers and $`K=𝐂(x)`$ be the field of rational functions on $`mn`$ variables $`x_{ij}`$, $`i=1,2,\mathrm{},m`$, $`j=1,2,\mathrm{},n`$. We extend the indices $`i,j`$ of $`x_{ij}`$ for $`i,j𝐙`$ by the condition $`x_{i+m,j}=x_{i,j+n}=x_{ij}`$.

Define algebra automorphisms $`\pi `$, $`\rho `$, $`r_i`$ and $`s_j`$ ($`i𝐙/m𝐙`$, $`j𝐙/n𝐙`$) on the field $`K`$ as follows:

$`\pi (x_{ij})=x_{i+1,j},\rho (x_{i,j})=x_{i,j+1},`$ (1)
$`r_i(x_{ij})=x_{i+1,j}{\displaystyle \frac{P_{i,j1}}{P_{ij}}},r_i(x_{i+1,j})=x_{ij}{\displaystyle \frac{P_{ij}}{P_{i,j1}}},`$ (2)
$`r_k(x_{ij})=x_{ij},(ki,i+1),`$ (3)
$`s_j(x_{ij})=x_{i,j+1}{\displaystyle \frac{Q_{i1,j}}{Q_{ij}}},s_j(x_{i,j+1})=x_{ij}{\displaystyle \frac{Q_{ij}}{Q_{i1,j}}},`$ (4)
$`s_k(x_{ij})=x_{ij},(kj,j+1),`$ (5)
where

$`P_{ij}={\displaystyle \underset{a=1}{\overset{n}{}}}\left({\displaystyle \underset{k=1}{\overset{a1}{}}}x_{i,j+k}{\displaystyle \underset{k=a+1}{\overset{n}{}}}x_{i+1,j+k}\right),`$ (6)
$`Q_{ij}={\displaystyle \underset{a=1}{\overset{m}{}}}\left({\displaystyle \underset{k=1}{\overset{a1}{}}}x_{i+k,j}{\displaystyle \underset{k=a+1}{\overset{m}{}}}x_{i+k,j+1}\right).`$ (7)
The following is the main result of this note.

###### Theorem 2.1

$`\pi ,r_0,r_1,\mathrm{},r_{m1}`$ and $`\rho ,s_0,s_1,\mathrm{},s_{n1}`$ generate the extended affine Weyl group $`\stackrel{~}{W}(A_{m1}^{(1)})`$ and $`\stackrel{~}{W}(A_{n1}^{(1)})`$. Moreover these two actions $`\stackrel{~}{W}(A_{m1}^{(1)})`$ and $`\stackrel{~}{W}(A_{n1}^{(1)})`$ mutually commute.

In our previous work , the first part of this theorem has been proved and the second part was conjectured. We will prove Theorem 2.1 in the next section.

###### Example 2.2

For the case of $`(m,n)=(2,3)`$. Let us write the $`x_{ij}`$ variables as

$$X=\left[\begin{array}{ccc}x_{11}& x_{12}& x_{13}\\ x_{21}& x_{22}& x_{23}\end{array}\right]=\left[\begin{array}{ccc}x_1& x_2& x_3\\ y_1& y_2& y_3\end{array}\right]$$
(8)
Then the actions $`r_1,r_0`$, $`\pi `$, $`s_1,s_2,s_0`$ and $`\rho `$ are given as follows.

$$\begin{array}{c}r_1(x_1)=y_1\frac{x_1x_2+x_1y_3+y_2y_3}{x_2x_3+x_2y_1+y_3y_1},r_1(y_1)=x_1\frac{x_2x_3+x_2y_1+y_3y_1}{x_1x_2+x_1y_3+y_2y_3},\hfill \\ r_1(x_2)=y_2\frac{x_2x_3+x_2y_1+y_3y_1}{x_3x_1+x_2y_2+y_1y_2},r_1(y_2)=x_2\frac{x_3x_1+x_2y_2+y_1y_2}{x_2x_3+x_2y_1+y_3y_1},\hfill \\ r_1(x_3)=y_3\frac{x_3x_1+x_2y_2+y_1y_2}{x_1x_2+x_1y_3+y_2y_3},r_1(y_3)=x_3\frac{x_1x_2+x_1y_3+y_2y_3}{x_3x_1+x_2y_2+y_1y_2},\hfill \end{array}$$
(9)
$`\pi (x_i)=y_i`$, $`\pi (y_i)=x_i`$, $`r_0=\pi r_1\pi `$ and

$$\begin{array}{ccc}s_1(x_1)=x_2\frac{x_1+y_2}{x_2+y_1},\hfill & s_1(x_2)=x_1\frac{x_2+y_1}{x_1+y_2},\hfill & s_1(x_3)=x_3,\hfill \\ s_1(y_1)=y_2\frac{x_2+y_1}{x_1+y_2},\hfill & s_1(y_2)=y_1\frac{x_1+y_2}{x_2+y_1},\hfill & s_1(y_3)=y_3,\hfill \end{array}$$
(10)
$`\rho (x_i)=x_{i+1}`$, $`\rho (y_i)=y_{i+1}`$ ,$`i𝐙/(3𝐙)`$, $`s_2=\rho s_1\rho ^1`$ and $`s_0=\rho ^2s_1\rho ^2`$.

## 3 Proof of Theorem 2.1

The first part (affine Weyl group relations): Though this part of the theorem has already been proved in  by direct computation, we give more conceptual proof here. We prove that $`s_i`$ ($`i=1,\mathrm{},n1`$) satisfy the same relations as those for the adjacent permutations $`\sigma _i=(i,i+1)S_n`$. The relations containing $`s_0`$ follow from the cyclic symmetry of the definition (4,5,7).

Using the identities,

$$x_{i,l+1}Q_{i1,l}x_{il}Q_{il}=k_{l+1}k_l,k_l=\underset{i=1}{\overset{m}{}}x_{il}.$$
(11)
we see that the automorphisms $`s_j`$ satisfy the relations,

$$x_iy_i=x_i^{}y_i^{},x_i+y_{i+1}=x_i^{}+y_{i+1}^{},$$
(12)
where $`x_i=x_{ij}`$, $`y_i=x_{i,j+1}`$, $`x_i^{}=s_j(x_{ij})`$ and $`y_i^{}=s_j(x_{i,j+1})`$. The relations (12) can be written in matrix form as

$$A(x)A(y)=A(x^{})A(y^{}),$$
(13)
where

$$A(x)=\left(\begin{array}{ccccc}x_1& 1& & & \\ & x_2& 1& & \\ & & \mathrm{}& \mathrm{}& \\ & & & x_{m1}& 1\\ z& & & & x_m\end{array}\right).$$
(14)
The relations (12) essentially characterize the transformation $`s_j`$. More generally we have the following,

###### Lemma 3.1

For given $`mn`$ variables $`x_j=(x_{1j},x_{2j},\mathrm{},x_{mj})`$ and $`mn`$ unknowns $`x_j^{}=(x_{1j}^{},x_{2j}^{},\mathrm{},x_{mj}^{})`$, a system of algebraic equations

$$A(x_1)A(x_2)\mathrm{}A(x_n)=A(x_1^{})A(x_2^{})\mathrm{}A(x_n^{}),$$
(15)
has $`n!`$ solutions. Each of the solution corresponds to a permutation $`\sigma S_n`$ and characterized by additional condition

$$x_{1j}^{}x_{2j}^{}\mathrm{}x_{mj}^{}=x_{1\sigma (j)}x_{2\sigma (j)}\mathrm{}x_{m\sigma (j)}.$$
(16)
Proof of Lemma 3.1. ¿From the definition of the transformation $`s_j`$ (4,5,7), it follows that

$$A(s_j(x_1))\mathrm{}A(s_j(x_n))=A(x_1)\mathrm{}A(x_n),$$
(17)
and

$$s_j(x_{1k}\mathrm{}x_{mk})=x_{1\sigma _j(k)}\mathrm{}x_{m\sigma _j(k)},k=1,\mathrm{},n,$$
(18)
where $`\sigma _j=(j,j+1)`$. This implies that for any permutation $`\sigma `$ given as a product $`\sigma =\sigma _{i_1}\mathrm{}\sigma _{i_k}`$,

$$x_{ij}^{}=s_{i_1}\mathrm{}s_{i_k}(x_{ij}),$$
(19)
is a solution of (15) and (16). Hence we see that there exist at least $`n!`$ solutions for (15).

Let us next show that such solution is unique, namely we prove that the number of solutions for (15) is at most $`n!`$. We discuss the case $`n=3`$ for simplicity. Consider a system of algebraic equations for $`3m`$ unknown variables $`(x_i^{},y_i^{},z_i^{})`$ given by

$$A(x)A(y)A(z)=A(x^{})A(y^{})A(z^{}),$$
(20)
or

$$a_i=a_i^{},b_i=b_i^{},c_i=c_i^{},$$
(21)
where

$$a_i=x_i+y_{i+1}+z_{i+2},b_i=x_iy_i+x_iz_{i+1}+y_{i+1}z_{i+1},c_i=x_iy_iz_i;$$
(22)
$`a_i^{},b_i^{}`$ and $`c_i^{}`$ are defined similarly. By simple elimination of variables $`y_i^{}`$ and $`z_i^{}`$, one has algebraic equations for $`x_i^{}`$

$$x_{i1}^{}x_i^{}x_{i+1}^{}a_{i1}x_i^{}x_{i+1}^{}+b_ix_{i+1}^{}c_{i+1}=0,$$
(23)
or equivalently

$$x_{i+1}^{}\psi _i=R_{i+1}\psi _{i+1},R_{i+1}=\left[\begin{array}{ccc}& 1& \\ & & 1\\ c_{i+1}& b_i& a_{i1}\end{array}\right],\psi _i=\left[\begin{array}{c}1\\ x_i^{}\\ x_{i1}^{}x_i^{}\end{array}\right].$$
(24)
Due to the periodicity $`x_{i+m}^{}=x_i^{}`$, we have

$$(x_1^{}x_2^{}\mathrm{}x_m^{})\psi _1=(R_m\mathrm{}R_2R_1)\psi _1.$$
(25)
Hence, $`\psi _1`$ is one of the three eigenvectors of $`R_m\mathrm{}R_2R_1`$. Once $`\psi _1`$ is determined other $`\psi _i`$’s and $`x_i^{}`$’s are determined uniquely by using (24) repeatedly. Then we have the equation for $`y_i^{}`$’s and $`z_i^{}`$’s

$$A(x^{})^1A(x)A(y)A(z)=A(y^{})A(z^{}).$$
(26)
By similar argument we see that this equation has two solutions. Hence we have at most $`3!`$ solutions for (20).

If $`\sigma =\sigma _{i_1}\mathrm{}\sigma _{i_k}=\sigma _{j_1}\mathrm{}\sigma _{j_l}`$ then

$$s_{i_1}\mathrm{}s_{i_k}(x_{1r}\mathrm{}x_{mr})=x_{1\sigma (r)}\mathrm{}x_{m\sigma (r)}=s_{j_1}\mathrm{}s_{j_l}(x_{1r}\mathrm{}x_{mr}),$$
(27)
and we have $`s_{i_1}\mathrm{}s_{i_k}=s_{j_1}\mathrm{}s_{j_l}`$ due to Lemma 3.1. This means that all the relations satisfied by $`\sigma _i`$’s also hold for $`s_i`$’s.

The second part (commutativity): Let us put $`G_i(u)=1+uE_{i+1,i}`$ for $`i[1,m1]`$ and $`G_0(u)=1+{\displaystyle \frac{u}{z}}E_{1,m}`$. We also put

$$M=A(x_1)A(x_2)\mathrm{}A(x_n),$$
(28)
where $`x_j=(x_{1j},\mathrm{}x_{mj})`$. We remark that Lemma 3.1 also implies the invariance of $`M`$ under the actions of $`s_j`$ ($`j=1,2,\mathrm{},n1`$). The action of $`r_k`$ on $`M`$ is described as follows.

###### Lemma 3.2

$$G_k(u)M=r_k(M)G_k(u),u=\frac{h_kh_{k+1}}{P_{k,0}}.$$
(29)
Proof. Putting

$$u_j=\frac{h_kh_{k+1}}{P_{kj}}(j=0,1,\mathrm{},n),$$
(30)
the action of $`r_k`$ can be rewritten in the form

$$\begin{array}{c}r_k(x_{kj})=x_{kj}u_j,r_k(x_{k+1,j})=x_{k+1,j}+u_{j1},\\ r_k(x_{ij})=x_{ij}(ik,k+1),\end{array}$$
(31)
where we have used the identity

$$x_{k+1,j}P_{k,j1}x_{kj}P_{kj}=h_{k+1}h_k,h_k=\underset{j=1}{\overset{n}{}}x_{kj}.$$
(32)
The relations above are equivalent to the matrix equation

$$G_k(u_{j1})A(x_j)=A(r_k(x_j))G_k(u_j)(j=1,\mathrm{},n).$$
(33)
Hence we have

$$G_k(u_0)A(x_1)\mathrm{}A(x_n)=A(r_k(x_1))\mathrm{}A(r_k(x_n))G_k(u_n).$$
(34)
Due to the periodicity $`u_0=u_n(=u)`$, we have

$$G_k(u)M=r_k(M)G_k(u),u=\frac{h_kh_{k+1}}{P_{k,0}},$$
(35)
where $`M=A(x_1)\mathrm{}A(x_n)`$.

We remark that the parameter $`u`$ can also be determined from $`M`$ by the formula

$$u=\frac{M_{kk}M_{k+1,k+1}}{M_{k,k+1}}|_{z=0}.$$
(36)
The commutativity of the actions of the two Weyl groups in Theorem 2.1 follows from Lemma 3.2, because the matrix $`M`$ is invariant under the actions of $`s_j`$ ($`j=1,2,\mathrm{},n1`$). The commutativity $`r_is_0=s_0r_i`$ follows from $`r_0=\rho ^1s_1\rho `$ and $`r_i\rho =\rho r_i`$. Thus, the proof of Theorem 2.1 is completed.

In this proof, we have treated the actions $`s_i`$ and $`r_j`$ asymmetrically. However the argument can be applicable for both of them, because the roles of $`s_i`$ and $`r_j`$ can be exchanged with each other by the transposition of the matrix $`X=(x_{ij})`$.

## 4 Discrete dynamical systems

The affine Weyl group $`\stackrel{~}{W}(A_{m1}^{(1)})\times \stackrel{~}{W}(A_{n1}^{(1)})`$ has a translation subgroup $`𝐙^{m1}\times 𝐙^{n1}`$. If we consider a part of these translations as providing a discrete dynamical system, the commutant of them can be viewed as its Bäcklund transformations.

###### Example 4.1

When $`m=2`$, the translation subgroup of $`\stackrel{~}{W}(A_1^{(1)})`$ is generated by $`T=\pi r_1`$. Its explicit actions $`T(x_i)=\overline{x}_i`$ etc. are determined as the nontrivial solution of

$$\overline{x}_j\overline{y}_j=x_jy_j,\overline{x}_j+\overline{y}_{j+1}=x_{j+1}+y_j,$$
(37)
where $`x_j=x_{1j}`$ and $`y_j=x_{2j}`$ ($`j𝐙/n𝐙`$). This system is a version of the discrete Toda equation. The explicit time evolution is given by

$$\overline{x}_j=x_j\frac{P_{j1}}{P_j},\overline{y}_j=y_j\frac{P_j}{P_{j1}},$$
(38)
where

$$P_j=\underset{a=1}{\overset{n}{}}\stackrel{a1}{\stackrel{}{y_{j+1}\mathrm{}y_{j+a1}}}\stackrel{na}{\stackrel{}{x_{j+a+1}\mathrm{}x_{j+n}}}.$$
(39)
The Weyl group $`\stackrel{~}{W}(A_{n1}^{(1)})`$ commutes with this discrete evolution $`T`$. The explicit actions are

$$\begin{array}{cc}s_j(x_j)=x_{j+1}\frac{x_j+y_{j+1}}{x_{j+1}+y_j},\hfill & s_j(y_j)=y_{j+1}\frac{x_{j+1}+y_j}{x_j+y_{j+1}},\hfill \\ s_j(x_{j+1})=x_j\frac{x_{j+1}+y_j}{x_j+y_{j+1}},\hfill & s_j(y_{j+1})=y_j\frac{x_j+y_{j+1}}{x_{j+1}+y_j},\hfill \end{array}$$
(40)
$`s_j(x_k)=x_k`$, $`s_j(y_k)=y_k`$, ($`kj,j+1`$).

###### Example 4.2

When $`m=3`$, we use the variables $`x_j=x_{1j}`$, $`y_j=x_{2j}`$ and $`z_j=x_{3j}`$. Then $`\pi (x_j)=y_j`$, $`\pi (y_j)=z_j`$, $`\pi (z_j)=x_j`$, $`\rho (x_j)=x_{j+1}`$, $`\rho (y_j)=y_{j+1}`$ and $`\rho (z_j)=z_{j+1}`$ ($`j𝐙/n𝐙`$). The translation subgroup of $`\stackrel{~}{W}(A_2^{(1)})`$ is generated by

$$T_1=\pi r_2r_1,T_2=r_1\pi r_2.$$
(41)
We describe the action of $`T_1`$, then the action of $`T_2=\pi T_1\pi ^1`$ is given by the rotation. The actions $`T(x_i)=\overline{x}_i`$ etc. are determined from the system of algebraic equations

$$\begin{array}{c}\overline{x}_j\overline{y}_j\overline{z}_j=x_jy_jz_j,\hfill \\ \overline{x}_j\overline{y}_j+\overline{x}_j\overline{z}_{j+1}+\overline{y}_{j+1}\overline{z}_{j+1}=y_jz_j+y_jx_{j+1}+z_{j+1}x_{j+1},\hfill \\ \overline{x}_j+\overline{y}_{j+1}+\overline{z}_{j+2}=y_j+z_{j+1}+x_{j+2}.\hfill \end{array}$$
(42)
as a unique rational solution such that

$$\overline{x}_1\mathrm{}\overline{x}_n=y_1\mathrm{}y_n,\overline{y}_1\mathrm{}\overline{y}_n=z_1\mathrm{}z_n,\overline{z}_1\mathrm{}\overline{z}_n=x_1\mathrm{}x_n.$$
(43)
The system (42) can be regarded as a generalization of the discrete Toda equation. The explicit formulas for the time evolution $`T_1`$ are given in the form

$$\overline{x}_j=x_j\frac{F_{j1}}{F_j},\overline{y}_j=y_j\frac{G_{j1}F_j}{G_jF_{j1}},\overline{z}_j=z_j\frac{G_j}{G_{j1}}.$$
(44)
Here $`F_j`$ and $`G_j`$ are polynomials in $`x`$, $`y`$, $`z`$ such that $`F_j=\rho ^j(F_0)`$, $`G_j=\rho ^j(G_0)`$. When $`n=3`$ for example, $`F_0`$ and $`G_0`$ are given explicitly as

$$\begin{array}{c}F_0=x_1x_2x_{3}^{}{}_{}{}^{2}+x_1x_2x_3y_1+x_2x_3y_1y_2+x_3y_1y_2z_1+x_1x_2x_3z_2\hfill \\ +y_1y_2z_1z_2+x_2x_3y_1z_3+x_3y_1z_1z_3+x_2x_3z_2z_3,\hfill \\ G_0=x_2x_3+x_3z_1+z_1z_2.\hfill \end{array}$$
(45)
We remark that the affine Weyl group $`\stackrel{~}{W}(A_{n1}^{(1)})`$ provides Bäcklund transformations for the time evolutions $`T_1`$, $`T_2`$.

Our dynamical system has the following conserved quantities.

###### Proposition 4.3

The characteristic polynomial of the matrix

$$M=A(x_1)\mathrm{}A(x_n)$$
(46)
is invariant under the actions of the affine Weyl group $`\stackrel{~}{W}(A_{m1}^{(1)})\times \stackrel{~}{W}(A_{n1}^{(1)})`$.

In fact, by Lemma 3.1, the matrix $`M`$ itself is invariant under the action of $`s_i`$. On the other hand, the action of each $`r_j(M)`$ is conjugate to $`M`$ as we have seen in Lemma 3.2.

###### Example 4.4

For $`(m,n)=(2,3)`$, we have the following invariants.111 By the $`mn`$ duality, the same invariants can also be obtained as the characteristic polynomial of a $`2\times 2`$ matrix.

$$\begin{array}{c}det[A(x)A(y)+wI]=w^3+z^2+w^2(x_1y_1+x_2y_2+x_3y_3)+\hfill \\ z(x_1x_2x_3+y_1y_2y_3)+w(x_1x_2y_1y_2+x_2x_3y_2y_3+x_3x_1y_3y_1)+\hfill \\ zw(x_1+x_2+x_3+y_1+y_2+y_3)+x_1x_2x_3y_1y_2y_3.\hfill \end{array}$$
(47)
Since the Weyl group actions are tropical  (i.e. subtraction-free birational mappings), there exists a combinatorial counterpart obtained by the ultra-discretization

$$a\times ba+b,a+b\mathrm{min}(a,b).$$
(48)
Corresponding combinatorial dynamical system can be viewed as a periodic version of the Box-Ball systems. Similarly to the original Box-Ball system  and their generalizations , this dynamical system also has soliton type solutions.

The transformations as in Example 4.1 are quite similar to those for the $`q`$-Painlevé equation $`q`$-$`P_{IV}`$ and their Bäcklund transformations. In fact, the realization of the Weyl group as in (29), (36) is defined exactly in the same way as that of the Bäcklund transformations of the Painlevé equations  and their Lie-theoretic generalization studied in , Section 4. The $`q`$-$`P_{IV}`$ can be viewed as a non-autonomous deformation of our integrable system for $`(m,n)=(2,3)`$ of this paper.

In order to see this surprising coincidence more precisely, let us reformulate Example 4.1. Putting $`x_jy_j=b_j`$, $`a_j=b_j/b_{j+1}`$ and $`x_jx_{j+1}=\phi _jb_j`$, we have

$$\overline{b}_j=b_j,\frac{\phi _j+1}{\overline{\phi _j}+a_j}=\frac{x_{j+1}}{\overline{x_j}}.$$
(49)
Eliminating the $`x_j`$ variables, one obtain

$$(1+\phi _j)(1+\frac{1}{\phi _{j+1}})=(1+\frac{\overline{\phi _{j+1}}}{a_{j+1}})(1+\frac{a_j}{\overline{\phi _j}}).$$
(50)
In terms of the variables $`a_j`$ and $`\phi _j`$, the time evolution (38) is rewritten in the form

$$\begin{array}{c}\overline{a}_j=a_j,\overline{\phi }_j=a_j\phi _{j+1}\frac{g_{j+2}}{g_j},\hfill \\ g_j=1+\phi _j+\phi _j\phi _{j+1}+\mathrm{}+\phi _j\mathrm{}\phi _{j+n2}.\hfill \end{array}$$
(51)
The Weyl group $`\stackrel{~}{W}(A_{n1}^{(1)})`$ commutes with this discrete evolution. In terms of variables $`a_j`$ and $`\phi _j`$, the Bäcklund transformations (40) are given by

$$\begin{array}{ccc}s_j(a_j)=\frac{1}{a_j},\hfill & s_j(a_{j+1})=a_{j+1}a_j,\hfill & s_j(a_{i1})=a_{i1}a_j,\hfill \\ s_j(\phi _j)=\frac{\phi _j}{a_j},\hfill & s_j(\phi _{j+1})=\phi _{j+1}\frac{a_j+\phi _j}{1+\phi _j},\hfill & s_j(\phi _{i1})=\phi _{i1}\frac{1+\phi _j}{1+\phi _j/a_j},\hfill \end{array}$$
(52)
and $`s_j(a_k)=a_k`$, $`s_j(\phi _k)=\phi _k`$ for $`kj,j\pm 1`$ (mod $`n`$).

In this setting, the parameters $`a_j`$ are subject to the constraint

$$qa_0a_1\mathrm{}a_{n1}=1.$$
(53)
However, the Weyl group representation in (52) and its commutativity with the evolution (51) survive also for $`q1`$. In this non-autonomous situation, the equation (51) gives a generalized $`q`$-Painlevé equation which has the $`\stackrel{~}{W}(A_{n1}^{(1)})`$ Bäcklund transformations (52); the case $`n=3`$ recovers the $`q`$-$`P_{IV}`$ equation . The non-autonomous versions of the cases with $`m3`$ could be considered as a $`q`$-Painlevé hierarchy with multi-discrete times.


format: text

text:
# Distribution of Interacting Ionic Particles in Disordered Media

## Abstract

Equilibrium distribution of interacting ionic particles in a charged disordered background is studied using the nonlinear Poisson-Boltzmann equation. For an arbitrarily given realization of the disorder, an exact solution of the equation is obtained in one dimension using a mapping of the nonlinear Poisson-Boltzmann equation to a self-consistent Schrödinger equation. The resulting density profile shows that the ions are delocalized, despite what the equivalent Schrödinger formulation in one dimension would suggest. It is shown that the ions are not distributed so as to locally neutralize the background, presumably due to their mutual interactions.

Dynamics of interacting charged particles in a heterogeneous Coulombic background, which is a medium with a quenched distribution of charge centers, is a complex phenomenon involving long equilibration times, a multitude of time scales and various intermediate regimes, and competing mechanisms . All these complications, which are primarily due to the long-range nature of Coulomb interaction, have so far made it extremely difficult for even rather sophisticated simulation techniques to provide us with a conclusive picture.

Perhaps a most preliminary step is to determine whether the equilibrium distribution of the mobile ions corresponds to a localized or an extended state, in analogy to the quantum mechanical (QM) problem of disordered electronic systems . The answer to this question is determined by a delicate balance between the following competing effects: (i) attraction of the mobile ions to the opposite fixed charges and the tendency towards local neutralization, (ii) repulsion between similarly charged mobile ions, and (iii) entropy of the mobile ions. (Note that for simplicity of the argument we have assumed that all the mobile ions are similarly charged whereas all the fixed charges are of the opposite type.) While attraction tends to localize the mobile ions in the neighborhood of the charge centers, repulsion and entropy both favor delocalization. One can imagine that determining the outcome of this balance in simulations could be a very difficult task due to the progressively long escape times for the deep Coulomb potential wells. On the other hand, it is also not clear if one can trust approximate theoretical schemes, which may artificially weaken or strengthen any of the above tendencies, in making feasible predictions to this end.

For example, it has been shown that if we neglect the mutual repulsion between the mobile ions, the answer to the above question depends on the dimensionality . One can show that the disorder-averaged density-density correlation function for the mobile ions behaves as $`\overline{c(\mathrm{𝟎})c(𝐱)}\mathrm{exp}\left[\left(x/\xi \right)^{2d}\right]`$ in $`d`$ dimensions , and therefore the ions are localized for $`d<2`$ (with $`\xi `$ being the localization length). On general grounds one might then expect that the interaction between the mobile ions could be effectively taken into account by introducing the notion of screening, which would change the density-density correlation function to $`\overline{c(\mathrm{𝟎})c(𝐱)}\mathrm{exp}\left[\left(x/\xi \right)^{2d}g(\kappa x)\right]`$, where $`g`$ is an exponentially decaying function and $`\kappa ^1`$ is the corresponding Debye screening length that is related to the density of the mobile ions. The appearance of the exponential damping now controls the divergence for $`d<2`$, and consequently the mobiles ions will be delocalized in any dimension. The simplistic implementation of screening using Debye-Hückel theory, however, is based on a linearized approximation of the nonlinear Poisson-Boltzmann (PB) theory, which is known to break down in the vicinity of the charge centers where the mobile ions could presumably be trapped for relatively long times. Therefore, it seems that a reliable answer to the localization question can only be obtained through an exact treatment of the full nonlinear PB formulation. This is what I attempt to do in this article.

I consider a one dimensional medium with a quenched distribution of charge centers and a finite density of oppositely charged mobile ions, and show that one can map the corresponding nonlinear PB equation to a self-consistent Schrödinger equation. An exact solution of the resulting equation is then presented for an arbitrary realization of the disorder. The resulting density profile shows that the ions are in an extended state. It is shown that the mutual interactions of the mobile ions could lead to a gross violation of the anticipated local neutrality.

Consider a set of fixed negatively charged ions described by the given charge density $`\rho _\mathrm{f}(𝐱)`$, embedded in a cloud of positively charged mobile ions with concentration $`c(𝐱,t)`$. The electric potential $`\varphi `$(x,t) is then given by the solution of the Poisson equation:

$$^2\varphi =\frac{4\pi }{ϵ}(ec+\rho _\mathrm{f}),$$
(1)
in which $`ϵ`$ is the dielectric constant of the background, and $`e`$ is the electron charge. The dynamics of $`c(𝐱,t)`$ is governed by a conservation law of the form $`_tc+𝐉=0`$, in which $`𝐉=Dc+c𝐯`$ is the current density, where $`𝐯(𝐱,t)`$ is the velocity field corresponding to the mobile ions, and $`D`$ is their diffusion constant. In a mean-field approximation, the velocity field of the ions at each point is determined by the local value of the electric field $`𝐄(𝐱,t)`$ as $`𝐯=\mu 𝐄=\mu \varphi `$, where $`\mu `$ is the electric mobility of the ions. Using the self-consistency relation, the continuity equation can be written as

$$_tc=D^2c+\mu (c\varphi ).$$
(2)
The above two equations describe the nonlinear dynamics of mobile ions in the presence of a given distribution of fixed ions . At equilibrium we have $`_tc=0`$, and Eq.(2) yields $`c\mathrm{exp}(\mu \varphi /D)`$, where equilibration at a temperature $`T`$ implies an Einstein relation $`D=(\mu /e)k_BT`$. Inserting the Boltzmann form for $`c`$ back into the Poisson equation then yields the PB equation.

Let us now restrict ourselves to the one dimensional case where the equations prove to be more tractable. Eliminating the concentration field $`c`$ from the two equations, one obtains an equation for the electrostatic potential, which can be integrated twice to yield

$`_t\varphi (x,t)`$ $``$ $`_t\varphi (\mathrm{},t)=D_x^2\varphi +{\displaystyle \frac{\mu }{2}}(_x\varphi )^2`$ (4)
$`+{\displaystyle \frac{4\pi }{ϵ}}\left[\mu {\displaystyle _{\mathrm{}}^x}𝑑x^{}\rho _\mathrm{f}(x^{})_x\varphi (x^{})+D\rho _\mathrm{f}(x)\right],`$
where the boundary condition $`E(\pm \mathrm{},t)=0`$ has been implemented. The above equation can be called a dynamical PB equation .

We focus, for simplicity, on the specific example of a series of $`N+1`$ negatively charged permeable plates with charge (number) density $`\sigma `$ and surface area $`A`$, which are placed at random locations $`x_n`$ on the $`x`$-axis, i.e., $`\rho _\mathrm{f}(x)=e\sigma _{n=0}^N\delta (xx_n)`$.

Using a so-called Cole-Hopf transformation $`\varphi (x,t)\varphi (\mathrm{},t)=(2D/\mu )\mathrm{ln}W(x,t)`$, Eq.(4) can be written as a Schrödinger-like equation for $`W(x,t)`$:

$`_tW(x,t)`$ $`=`$ $`D_x^2W(x,t)`$ (5)
$`+`$ $`{\displaystyle \frac{2D}{\lambda }}{\displaystyle \underset{n=0}{\overset{N}{}}}\left[\epsilon _n(t)\mathrm{\Theta }(xx_n)\delta (xx_n)\right]W(x,t),`$ (6)
where

$$\epsilon _n(t)=\left(\frac{_xW(x_n^+,t)+_xW(x_n^{},t)}{W(x_n,t)}\right),$$
(7)
and $`\lambda ϵD/(\pi e\sigma \mu )=ϵk_BT/(\pi e^2\sigma )`$ is a characteristic length scale corresponding to the relative strength of the electrostatic and thermal interactions. Note that Eq.(5) is still nonlinear due to the self-consistent coupling of the so-called asymmetry fields $`\epsilon _n(t)=\mu E(x_n,t)/D`$.

The value of the asymmetry field at each point can be calculated from Gauss’ law as

$$\epsilon _n(t)=\frac{2}{\lambda }\left(\frac{2Q_n}{e\sigma A}\right),$$
(8)
in which $`Q_n`$ is the net amount of (uncompensated) charge on the left hand side of point $`x_n`$. (Overall neutrality then implies that $`Q_n`$ should be the net amount of charge on the right hand side of point $`x_n`$.) In other words, $`\epsilon _n(t)`$ measures the relative deviation from local neutrality as seen at $`x_n`$.

The dynamics described by Eq.(5) is difficult due to the time dependence of the $`\epsilon _n(t)`$ fields. We thus focus on the long time limit of the dynamics, when the ions are expected to equilibrate in the disordered medium. In this limit, the asymmetry fields are time independent, and the Schrödinger equation reads

$$\frac{\lambda }{2}_x^2W(x)+V(x)W(x)=0,$$
(9)
where the effective potential is given as

$$V(x)=\underset{n=0}{\overset{N}{}}\left[\epsilon _n\mathrm{\Theta }(xx_n)+\delta (xx_n)\right].$$
(10)
Note the peculiar feature that the delta function peaks in the potential, which correspond to the attractive sites for the mobile ions, appear as effective repulsive potential wells in the equivalent QM description!

To complete the QM formulation, we need a normalization for $`W(x)`$. Using the equilibrium relation $`c\mathrm{exp}(\mu \varphi /D)`$, and the definition of $`W`$, we can relate $`c`$ to $`W`$ as $`c(x)1/W^2(x)`$. Conservation of the overall number of mobile ions can thus be written as

$$\frac{1}{L}_{L/2}^{L/2}𝑑x\left[\frac{1}{W(x)}\right]^2=1,$$
(11)
where $`L`$ is the system size.

We can now solve the nonlinear self-consistent QM problem described by Eqs.(7), (9), and (11), because the potential is a constant in each interval. The solution of Eq.(9) is a series of plane waves:

$$W(x)=\{\begin{array}{cc}A_0+B_0x,\hfill & \mathrm{for}x<x_0,\hfill \\ & \\ A_je^{ik_j(xx_j)}\hfill & \\ +B_je^{ik_j(xx_j)},\hfill & \mathrm{for}x_{j1}<x<x_j,\hfill \\ & \\ A_{N+1}+B_{N+1}x,\hfill & \mathrm{for}x_N<x,\hfill \end{array}$$
(12)
where

$$k_j=\left[\frac{2}{\lambda }\underset{n=0}{\overset{j1}{}}\epsilon _n\right]^{1/2},$$
(13)
for $`j=1,\mathrm{},N`$. We should then match the neighboring solutions across the boundary points based on the following three requirements: (i) $`W(x)`$ should be continuous, (ii) the jump in $`_xW/W`$ should be set to $`2/\lambda `$ as implied by the delta function peaks in $`V(x)`$, and (iii) the asymmetry fields should satisfy the self-consistency relation of Eq.(7). These criteria translate into a series of self-consistency equations for the wave vectors:

$`e^{2ik_j\mathrm{\Delta }_j}`$ $`=`$ $`\left[{\displaystyle \frac{1+\frac{\lambda ^2}{4}(k_{j+1}^2k_j^2)ik_j\lambda }{1+\frac{\lambda ^2}{4}(k_{j+1}^2k_j^2)+ik_j\lambda }}\right]`$ (14)
$`\times `$ $`\left[{\displaystyle \frac{1+\frac{\lambda ^2}{4}(k_{j1}^2k_j^2)ik_j\lambda }{1+\frac{\lambda ^2}{4}(k_{j1}^2k_j^2)+ik_j\lambda }}\right],`$ (15)
for $`j=1,\mathrm{},N`$, where $`\mathrm{\Delta }_j=x_jx_{j1}`$ and $`k_0=k_{N+1}=0`$. For any given set of $`\mathrm{\Delta }_j`$s, one should solve these coupled nonlinear equations to determine the wave vectors. The matching requirements together with the normalization condition of Eq.(11), can then give us the full solution, and thus the density profile of the ions.

In order to determine whether the ions are in an extended state, or localized near the adsorbing sites, we need to find out whether the $`k_j`$s that satisfy Eq.(14) are real or imaginary. Instead of trying to solve for $`k_j`$s exactly, I shall first resort to approximate schemes. Physically we expect to have a tendency towards localization at low temperatures or for strong disorder, and nearly free ions with uniform distribution at high temperatures or for weak disorder.

Let us first focus on the case of strong disorder where we have $`\lambda \mathrm{\Delta }_j`$, and try to solve Eq.(14) perturbatively. I find

$`k_j`$ $`=`$ $`{\displaystyle \frac{\pi }{\mathrm{\Delta }_j}}{\displaystyle \frac{2\pi \lambda }{\mathrm{\Delta }_j^2}}+{\displaystyle \frac{4\pi \lambda ^2}{\mathrm{\Delta }_j^3}}\left(8\pi {\displaystyle \frac{\pi ^3}{6}}\right){\displaystyle \frac{\lambda ^3}{\mathrm{\Delta }_j^4}}`$ (17)
$`+{\displaystyle \frac{\pi ^3\lambda }{4\mathrm{\Delta }_j^2}}\left({\displaystyle \frac{\lambda ^2}{\mathrm{\Delta }_{j+1}^2}}+{\displaystyle \frac{\lambda ^2}{\mathrm{\Delta }_{j1}^2}}\right)+O\left({\displaystyle \frac{\lambda ^4}{\mathrm{\Delta }^5}}\right),`$
for the first few terms, which is real. This proves, strictly speaking, that if the ions are localized, the localization length should be greater than $`\mathrm{\Delta }^5/\lambda ^4`$. However, it looks rather unlikely that higher order corrections start to develop imaginary components, which implies that the ions should be delocalized even in the strong disorder limit.

To check this, I have solved Eq.(14) numerically for various choices of $`\{\mathrm{\Delta }_j\}`$, and the results for the values of $`\{k_j\}`$ were always real. For example, for the case of 11 charged centers with $`\mathrm{\Delta }_1/\lambda =0.1`$, $`\mathrm{\Delta }_2/\lambda =0.5`$, $`\mathrm{\Delta }_3/\lambda =2`$, $`\mathrm{\Delta }_4/\lambda =0.8`$, $`\mathrm{\Delta }_5/\lambda =0.01`$, $`\mathrm{\Delta }_6/\lambda =7`$, $`\mathrm{\Delta }_7/\lambda =3.2`$, $`\mathrm{\Delta }_8/\lambda =0.3`$, $`\mathrm{\Delta }_9/\lambda =1.6`$, and $`\mathrm{\Delta }_{10}/\lambda =0.0002`$, I find $`k_1\lambda =33.84`$, $`k_2\lambda =7.33`$, $`k_3\lambda =1.42`$, $`k_4\lambda =4.95`$, $`k_5\lambda =316.69`$, $`k_6\lambda =0.399`$, $`k_7\lambda =0.757`$, $`k_8\lambda =12.599`$, $`k_9\lambda =1.933`$, and $`k_{10}\lambda =15710`$, which are all real. This confirms that the ions are in an extended state.

We can also calculate the asymmetry fields from the above results as

$$\epsilon _j=\frac{\lambda }{2}\left(k_{j+1}^2k_j^2\right)\frac{\lambda }{2}\left(\frac{\pi ^2}{\mathrm{\Delta }_{j+1}^2}\frac{\pi ^2}{\mathrm{\Delta }_j^2}\right),$$
(18)
to the leading order. This shows that unlike what one might expect, the distribution of ions is not predominantly determined by the criterion of local neutralization of the background, which would yield $`\epsilon _j=0`$.

Now let us address the peculiar features of the formulation, and try to understand them better: (i) Why do the delta function peaks, corresponding to adsorbing sites in the real problem, appear as repulsive in the effective potential $`V(x)`$? (ii) While we expect a Schrödinger equation with random potential in one dimension to lead to localization , why do we find that the ions are delocalized? To answer these questions, we should note the unusual normalization condition of Eq.(11) above. A true wave function for this problem should be defined as $`\mathrm{\Psi }(x)=1/W(x)`$ according to Eq.(11), leading to a corresponding Schrödinger equation as

$$\frac{\lambda }{2}_x^2\mathrm{\Psi }(x)+\left[\lambda \left(\frac{_x\mathrm{\Psi }(x)}{\mathrm{\Psi }(x)}\right)^2V(x)\right]\mathrm{\Psi }(x)=0,$$
(19)
which is highly nonlinear.

While the sign of the effective potential is now reversed corresponding to attractive delta function peaks, an additional effective repulsive self-consistent potential has also appeared. This potential, which is proportional to the the electrostatic pressure $`(\varphi )^2`$, disfavors pile-up of the ions in the vicinity of the attractive sites and promotes delocalization. This extra self-consistent term, which is a direct consequence of the repulsive interactions between the mobile ions, explains the discrepancy between the predictions of this theory and conventional Schrödinger equation in a random potential.

In conclusion, I have presented an exact analysis of the problem of interacting mobile ions in a random background of fixed ions in one dimension. The density profile of the ions at equilibrium is shown to be extended (as opposed to localized) as a consequence of the interactions between the mobile ions. Since the one dimensional problem in the strong disorder limit that has the best chance of leading to localization is proven to result in an extended state, it is plausible to assert that the same result should hold true in general, i.e., in higher dimensionalities and for weaker disorder. Although I have only considered positively charged mobile ions in a negative random background, the formulation can be readily generalized to contain both species of mobile ions and fixed ions.

I am grateful to M.P.A. Fisher, M. Kardar, M. Mezard, and M. Sahimi for invaluable discussions and comments. This research was supported in part by the National Science Foundation under Grant No. DMR-98-05833.


text:
# Variable UV Absorption in the Spectrum of MRC 2251–178 1footnote 11footnote 1Based on observations with the NASA/ESA Hubble Space Telescope, obtained at the Space Telescope Science Institute, which is operated by AURA, Inc., under NASA contract NAS 5–26555.

## 1 Introduction

The identification of QSO-intrinsic narrow absorption lines222Hereafter, we abbreviate narrow absorption line with the acronym NAL and refer to a QSO that hosts a truly intrinsic NAL as a NALQSO. has exploded over the past four years with the advent of high resolution spectroscopy with large ground-based telescopes. With this explosion has come new insights on how to identify these absorbing systems - that is, how to separate them from randomly distributed intervening material. Historically, large ensembles of intrinsic absorbers have been identified using statistical arguments, either demonstrating an excess of absorbers in some redshift path (Foltz et al., 1986; Anderson et al., 1987) or by correlating the velocity distribution (with respect to the QSO emission redshift) of absorbers with QSO properties (Richards et al., 1999). Such surveys have demonstrated that, in high redshift QSOs, there is an excess of absorbers that lie within 5000 km s-1 of the QSO emission redshift (termed “associated” absorbers). It is now generally accepted that a large fraction of associated absorbers are indeed intrinsic. However, it is not truly known how large this fraction is. Moreover, Richards et al. (1999) demonstrated that a significant fraction of non-associated (that is, high ejection velocity) absorbers may also be intrinsic.

The two smoking guns, according to Barlow & Sargent (1997), that identify an intrinsic system are: (1) variability of the spectral profiles and/or equivalent widths over time \[e.g., Hamann et al. (1995); Barlow et al. (1997); Hamann et al. (1997b, c); Aldcroft et al. (1997)\]; and (2) the signature of partial coverage \[e.g., Barlow & Sargent (1997); Barlow et al. (1997); Hamann et al. (1997a); Ganguly et al. (1999)\]. Thusfar, all cases of confirmed QSO-intrinsic NALs are at high redshift, where the resonant UV doublets of metals are shifted into the optical. At low redshift, Ganguly et al. (2001) estimate from statistical excesses that $`25`$% of QSOs have intrinsic NALs. As a first step in determining that true fraction of intrinsic NALs, and deciphering the dominant cause of variability (if one exists) we have undertaken an HST/STIS snap shot program to identify intrinsic NALs in a sample of QSOs where associated NALs have been detected in archived HST/FOS spectra. In this effort, we also hope to look for high ejection velocity absorbers at low redshift and decipher the nature of the variability. We note that, up until now, only six QSOs have been shown to host intrinsic absorption through time variability - all at $`z_{\mathrm{em}}>1.5`$. Furthermore all cases of associated, time variable absorption have been been low luminosity AGN (Walter et al., 1990; Shull & Sachs, 1993; Maran et al., 1996; Weymann et al., 1997; Crenshaw et al., 2000).

In this letter, we report the first confirmation at low redshift of an associated NAL that is truly intrinsic to a QSO, MRC 2251–178 \[$`z_{\mathrm{em}}=0.066092`$; Bergeron et al. (1983)\]. The radio-quiet MRC 2251–178 is well-known for having a highly variable X-ray warm absorber, the first one known (Halpern, 1984). Monier et al. (2001) first reported an associated UV absorber and, assuming that the same gas produced the X-ray warm absorption, inferred that the mass loss rate from the disk and the accretion rate were comparable.

In §2, we present our HST/STIS snapshot of MRC 2251–178. In §3, we show that the associated UV absorption has varied over the past four years (thus, verifying its intrinsic nature). In §4, we discuss the implications of variability on analyses of intrinsic absorption lines and make suggestions for future work on MRC 2251–178

## 2 Data

The sample for our snap shot program is drawn from QSOs in the HST/FOS archive. The sample is restricted to only those QSOs which have been observed with one of the “high” resolution ($`230`$ km s-1) FOS gratings (G130H, G190H, G270H) and for which there is evidence of associated absorption in the form of a C IV $`\lambda \lambda 1548,1550`$, N V $`\lambda \lambda 1393,1403`$, O VI $`\lambda \lambda 1032,1038`$, or Ly $`\alpha `$ feature. The sample, comprised of 37 QSOs, includes 13 NALQSO candidates from Ganguly et al. (2001), which derived from a subsample of the HST Quasar Absorption Line Key Project QSOs, as well as 24 QSOs which were not part of the Key Project. The QSOs from the Key Project (Bahcall et al., 1993; Schneider et al., 1993; Bahcall et al., 1996; Jannuzi et al., 1998) were obtained from D. Schneider, B. Jannuzi, and S. Kirhakos. The non-Key Project QSOs were reduced and supplied by S. Kirhakos. All QSOs were obtained fully reduced with effective continuum333The effective continuum is the sum of the continuum and broad emission lines from the QSO. It is the spectrum that is presumably incident on the absorbing gas. fits according to Schneider et al. (1993).

The primary goal of the snap shot program is to detect variability of associated NALs. To that end, the QSOs are observed using HST/STIS with the G230L grating and the $`52\times 0.2`$” slit. The dispersion of the grating is 1.58 Å per pixel, with two pixels per resolution element. This corresponds to 300–600 km s-1 across the spectrum, compared to the roughly constant 230 km s-1 resolution of the FOS spectra. In Fig 1, we show the full spectrum of MRC 2251–178 (darker histogram) resulting from a 500 s exposure taken on 5 November 2000. The signal-to-noise per pixel at 2376 Å (the central wavelength) is $`25`$. The comparison HST/FOS spectrum, also shown in Fig 1 (lighter histogram), was taken on 2 August 1996 with all three 230 km s-1 resolution gratings (Monier et al., 2001). The relevant portion of the spectrum (taken with the G190H grating) has a signal-to-noise similar to the STIS spectrum. The flux calibration from the reduction pipelines for both the STIS-G230L and FOS-G190H spectra are accurate to within 3%. The wavelength calibration for the FOS-G190H spectrum is accurate to 0.11 diodes (0.2 Å) while that of the STIS-G230L spectrum is accurate to 0.2–0.5 pixels (0.3–0.8 Å). A STIS-G140M spectrum, obtained 3 February 1998 by J. Stocke, spans the wavelength range 1194–1300 Å (in two separate grating tilts). Although this spectrum covers the associated absorption from Ly $`\alpha `$ and N V transitions, there is, unfortunately, no overlap with our STIS-G230L. As a result, no direct comparison can be made to the STIS-G140M spectrm.

## 3 Analysis & Results

It is apparent by visual inspection of the STIS (Fig. 1) and FOS (Monier et al. 2001: Fig. 1) spectra that there has been a change in the state of both the emission and associated absorption lines of MRC 2251–178. The STIS spectrum has a generally higher flux and the associated C IV absorption has apparently disappeared. To examine the magnitude of the variability, we convolved the FOS spectrum with a Gaussian (575 km s-1 FWHM) and resampled the result with a 1.58 Å/pixel dispersion to mimic the expected STIS-G230L spectrum around the C IV emission line.

In Fig. 2, we show the C IV emission and associated absorption lines from the observed FOS spectrum (top panel; solid histogram), and the observed STIS spectrum (bottom panel; solid histogram). Overplotted on each are the effective continuum fits (smooth curve) and the expected STIS spectrum (dotted histograms). In the top panel, the expected STIS spectrum reflects what should have been observed if there were no variability in either the C IV emission or associated absorption lines. In the bottom panel, the expected STIS spectrum shows what should have been observed if there were no change in the properties of the associated absorber (that is, if only the shape of the C IV emission line and the continuum flux changed). With STIS-G230L, the associated C IV absorption is sampled by seven pixels (3.5 resolution elements) and should have been detected in the observed spectrum. (The Galactic Al II line at 1670 Å, which is sampled by less than a pixel in the STIS spectrum, is completely washed out by the instrumental profile, making it indistinguishable from the C IV and He II emission lines.) Furthermore, Monier et al. (2001) measured a C IV$`\lambda 1548`$ equivalent width of $`1.09\pm 0.09`$ Å in the FOS spectrum. Using the unresolved feature detection method from Schneider et al. (1993), the absorber is not formally detected in the STIS spectrum down to a $`3\sigma `$ equivalent width limit of $`0.19`$ Å. This is more than a $`10\sigma `$ change. It is, therefore, clear and robust that there has been variability in the properties of the associated absorber.

## 4 Discussion

### 4.1 Variability in the UV spectrum

The change in the state of associated C IV absorption happened sometime in the last four years. Although, we cannot say accurately on what timescale that variability occurred, it is now safe to say that this is the first clear demonstration of truly intrinsic narrow absorption in a low redshift QSO. Variability in absorption lines is generally attributed to one of two scenarios. In the first scenario, bulk motion, the absorption gas has an appreciable motion transverse to the line of sight and the total column density through the gas changes with time. The observational signature is that the equivalent width of all lines due to this gas should change in the same way (that is, all get stronger or all get weaker) with time. In the second scenario, ionization/recombination, the physical state of the absorber changes such that the absorbing gas becomes more or less ionized with time. Observationally, this would mean that the equivalent widths of higher-ionization lines and lower-ionization lines should change in opposite directions. This assumes, however, that the kinematics of the gas does not change so that the column densities of the ions drive the changes in the equivalent widths. Unfortunately, the four year separation of the two UV observations and the lack of coverage of relevant species by the STIS G230L observation (e.g., Ly $`\alpha `$, N V $`\lambda \lambda 1393,1403`$, Si IV $`\lambda \lambda 1394,1403`$, Si II$`\lambda 1260`$, C II$`\lambda 1335`$) do not allow one to distinguish between these two causes of the variability in the spectrum of MRC 2251–178.

Nevertheless, since we know that the variability timescale is less than 4 years, we can provide an upper limit on the distance between the absorbing material and the ionizing sources. If we assume that the cause of variability is recombination of C IV to C III and that the temperature of the absorber is about $`T_\mathrm{e}=2\times 10^5`$ K as in UM 675 (Hamann et al., 1995), then the minimum density is $`n_\mathrm{e}=1/t_{\mathrm{rec}}\alpha >3000\mathrm{cm}^3`$, where the recombination coefficient is $`\alpha =2.8\times 10^{12}\mathrm{cm}^3\mathrm{s}^1`$ (Arnaud & Rothenglug, 1985). This implies a maximum distance of 2.4 kpc, if C IV is the dominant species. This is consistent with the estimate from Monier et al. (2001). While not very restrictive, the distance estimate does rule a number of candidates for the absorbers, including the intergalactic medium of the host cluster and most of the interstellar medium of the host galaxy. The remaining possibilities are the material mixed with broad emission line region (BELR), material outside the BLR, and outflowing/ejected material.

### 4.2 Suggestions for Future Work

The next step in understanding the intrinsic absorber in MRC 2251–178 is to decipher the cause of variability - bulk motion or ionization/recombination - and then determine if there is a direct physical relationship between the gas and the X-ray warm absorber. In principle, these issues can be tackled, but they require a regimented observing schedule. To efficiently determine both the underlying cause of variability and the geometry of absorbing gas relative to the central engine, periodic observations at medium to high resolution should be obtained of the rest-frame wavelength range 1200–1600 Å (e.g., with the E140M echelle on STIS). This covers transitions of Si II, Ly $`\alpha `$, N V, C II, Si IV, and C IV. The frequency of the observations should be less than (or on the order of) once in two months. If the QSO and/or absorber properties vary within the course of observations, this schedule should provide the variability timescale. In addition, at high resolution, photoionization modeling of the individual components should constrain the physical conditions of the gas. Furthermore, measurements of the partial coverage fraction are possible (and necessary) with the resolved components. The combination of the partial coverage fraction, the physical conditions of the gas (e.g. the ionization parameter, density), and the variability timescale can constrain the distance from the central engine, the size, and the transverse velocity of the absorbers.

Once the range of conditions in the UV absorbing gas are known, one can then consider its relationship to the X-ray warm absorber. This, however, is difficult since the X-ray warm absorber is known to vary on short timescales ($`<1`$ year; Halpern 1984). Thus, simultaneous observations of MRC 2251–178 in both the UV and the X-ray are required so that inferred relationships are not affected by the changing conditions. In the X-ray, observations can be carried out with either the Chandra/HETG or XMM-Newton/RGS to cover O VII and O VIII. Ultraviolet observations should be carried out with either HST/STIS-E140M, which covers a range of ionization states, and/or with FUSE, which covers the O VI $`\lambda \lambda 1032,1038`$ doublet and the Lyman series.

Support for this work was provided under grants HST–GO–08681.01–A and STSI AR–08763.01–A. We are grateful to Sofia Kirhakos for providing the uniformly and fully reduced FOS archive with continuum and emission line fits.


text:
# Breakdown of Migdal–Eliashberg theory via catastrophic vertex divergence at low phonon frequency

## I Introduction

The use of Migdal–Eliashberg (ME) theory for the study and analysis of electron-phonon systems is widespread. For conventional superconductors such as lead, where the electron-phonon coupling is higher than can be treated with BCS theory, ME theory has been extremely successful for understanding the superconducting properties. Recently, a lot of researchers have been interested in an apparent kink in the electronic dispersion of cuprate supercondutors as determined from angle resolved photo-emission spectroscopy (ARPES), and their analysis / interpretation typically uses the related Engelsberg–Schrieffer (ES) result (i.e. just the lowest order Fock diagram) Engelsberg and Schrieffer (1963); Lanzara et al. (2001); Cuk et al. (2004). The ES result divides excitations into long-lived (coherent) low-energy excitations, and rapidly decaying high-energy excitations, with a kink at the phonon energy. ES based analysis of ARPES results suggests a large electron-phonon coupling in the cuprates, with estimates of the dimensionless coupling constant lying between $`0.3`$ and $`1.5`$ (depending on doping) and very large phonon frequencies of $`50`$ meV Lanzara et al. (2001). There has also been a development of a maximum entropy technique for analysing cuprate superconductors which makes use of Eliashberg theory Hwang et al. (2007). The analysis in Ref. Hwang et al., 2007 determines a dimensionless electron-phonon coupling of $`\lambda >1`$, but this is well above the $`\lambda 1`$ value where one would normally expect perturbation theory to fail. The region of intermediate electron-phonon coupling is also relevant to a number of other materials. For example, electron paramagnetic resonance measurements on the manganites support a strong “electron-phonon” coupling leading to Jahn–Teller polarons Zhao et al. (1996). In light of the current experimental situation, and the importance of the conclusions of Refs. Cuk et al., 2004; Hwang et al., 2007, a full re-analysis of the perturbation theory is of high significance.

The neglect of vertex corrections suggested by Migdal leads to a theory where an infinite set of Feynman diagrams may be summed. According to an analysis carried out by Migdal, this theory should be valid in the physical regime of electron-phonon problems, where the phonon energy is significantly smaller than the intersite hopping (“Migdal’s theorem”) Migdal (1958); Eliashberg (1960) specifically that the vertex corrections are small when $`\lambda \omega _0/ϵ_f1`$ ($`\omega _0`$ is the phonon frequency and $`ϵ_f`$ the Fermi energy) 111Migdal’s theorem it is not a theorem in the mathematical sense, so we prefer to call it Migdal’s analysis in this article.. On the basis of Migdal’s analysis, it is often believed that Migdal–Eliashberg theory is applicable above $`\lambda 1`$ in the adiabatic limit because of the apparently small size of the vertex corrections, even though in general perturbative approaches break down (i.e. the functional form of the self-energy becomes incorrect) when the coupling constant becomes large. In this article, we revisit Migdal’s analysis by examining the large-$`d`$ limit (local approximation).

The large-$`d`$ limit has been an effective workshop for determining the validity of approximate schemes. The quantum Monte-Carlo (QMC) solution in the large-$`d`$ limit Freericks et al. (1994); Hirsch and Fradkin (1983) has been compared with a number of different diagrammatic approaches Freericks and Jarrell (1994); Freericks et al. (1998, 1993) demonstrating that self-consistent second-order perturbation theory worked better than Migdal-Eliashberg theory at intermediate phonon frequencies and strong coupling ($`\omega _00.2t`$, $`U>t`$). The QMC solution is difficult for very low phonon frequency and low temperature. QMC self-energies and Green functions are generated along the Matsubara axis. As such, it is not easy to make quantitative conclusions about the point of breakdown of ME theory as small differences at individual Matsubara frequencies may result in relatively large differences in the spectral functions. Using an alternative approach, Benedetti et al. Benedetti and Zeyher (1998) found the formation of more than one extremum in the path-integral formulation at very low phonon frequencies and intermediate coupling, leading to a breakdown of ME theory. Work on the finite-dimensional Holstein model also shows deficiencies in ME theory. On the basis of a comparison with exact diagonalization results, Alexandrov et al. Alexandrov et al. (1994) have shown that ME theory may break down at intermediate couplings even in the adiabatic limit.

In this article, we investigate the Holstein model of electron-phonon interactions Holstein (1959), which employs a number of approximations. The phonon dispersion is flat corresponding to independently moving ions and phonon anharmonicity is neglected, resulting in a Hamiltonian,

$$H=t\underset{<ij>\sigma }{}c_{i\sigma }^{}c_{j\sigma }+\underset{i\sigma }{}(gx_i\mu )n_{i\sigma }+\frac{M\omega _0^2}{2}\underset{i}{}x_i^2+\frac{1}{2M}\underset{i}{}p_i^2,$$
(1)
where $`c_{i\sigma }^{}`$ $`(c_{i\sigma })`$ create (annihilate) electrons at site $`i`$ with spin $`\sigma `$, $`x_i`$ is the local ion displacement, $`p_i`$ the ion momentum, $`M`$ the ion mass, $`t`$ the electron hopping parameter and $`\mu `$ the chemical potential. A mean-field Coulomb pseudopotential may also be included, but in the normal state, it is just absorbed into the chemical potential. In the following, $`t=0.5`$ and all energies (including temperature) are measured in units of $`2t`$.

An expression for the effective interaction between electrons can be obtained by performing a linear transformation in the phonon variable to remove the electron-phonon term Bickers and Scalapino (1989). The resulting interaction is retarded, with Fourier components:

$$U(i\omega _s)=\left(\frac{g^2}{M\omega _0^2}\right)\frac{\omega _0^2}{(\omega _s^2+\omega _0^2)}=U\frac{\omega _0^2}{(\omega _s^2+\omega _0^2)},$$
(2)
where $`\omega _s=2\pi sT`$ are the Matsubara frequencies for Bosons. Taking the limit $`\omega _0\mathrm{}`$, $`g\mathrm{}`$, while keeping the ratio $`g/\omega _0`$ finite, leads to an attractive Hubbard model Hubbard (1963) with instantaneous attraction of magnitude $`U=g^2/M\omega _0^2W\lambda `$, where $`W`$ is the half band-width (note that $`U`$ is related to the bipolaron binding energy for the strong coupling two-electron problem). In the opposite limit ($`\omega _00`$, $`M\mathrm{}`$, keeping $`M\omega _0^2\kappa `$ finite) the phonon kinetic-energy vanishes, leaving only a static variable $`x_i`$ representing the phonon subsystem in the Hamiltonian, $`H=t_{<ij>\sigma }c_{i\sigma }^{}c_{j\sigma }+_{i\sigma }(gx_i\mu )n_{i\sigma }+\frac{\kappa }{2}_ix_i^2`$. Thus the problem looks like that of a single electron in a disordered potential. The large-$`d`$ limit of this model was exactly solved by Millis et al. Millis et al. (1996) and extended to deal with long range order by Ciuchi et al Ciuchi and de Pasquale (1999). Thus, the phonon frequency may be thought of as a parameter for tuning the level of electronic correlation Hague and d’Ambrumenil (2005).

In this paper, we investigate the validity of the Migdal–Eliashberg approach within the dynamical mean-field theory formalism. As we carry out the self-consistency on the real-axis self-energies, we can investigate the behaviour of ME theory over a full range of temperatures and phonon frequencies. We begin by calculating self-energies, spectral functions, the quasi-particle weight and renormalized phonon frequency. To determine the validity of the theory, we calculate an expression for the lowest-order correction to the vertex function, and evaluate its magnitude. We find that the theory can also break down in the low frequency regime, contrary to the standard interpretation of Migdal’s analysis. Finally, we calculate resistivity curves in the regime where the lowest-order corrections are small.

## II Migdal–Eliashberg theory in the local approximation

In this article, we use the local approximation or dynamical mean-field theory (DMFT) as a way of analysing the Migdal-Eliashberg theory. The self-energy of correlated-electron systems is momentum independent in limit of large dimensions Metzner and Vollhardt (1989), and approximately momentum independent in 3D. In large-$`d`$, lattice models map onto Anderson impurity models with a self-consistent hybridisation Georges et al. (1992).

The self-consistent DMFT equations may be obtained by rewriting the action for the model to be considered in terms of an effective single-site action Georges et al. (1996)

$$S_{\mathrm{eff}}=_0^\beta 𝑑\tau _0^\beta 𝑑\tau ^{}\underset{\sigma }{}c_{i\sigma }^{}(\tau )𝒢_0^1(\tau \tau ^{})c_{i\sigma }(\tau ^{})+S_{\mathrm{int}}.$$
(3)
Here $`𝒢_0(i\omega _n)`$ plays the rôle of the host Green function in the equivalent impurity model. If one assumes that correlations carried via the bath between electrons entering and exiting a single site can be neglected (true in the case of large coordination number), the degrees of freedom associated with all but one site can be integrated out. This leads to the following self-consistent equation Georges et al. (1996),

$$𝒢_0^1(i\omega _n)=i\omega _n+\mu +G^1(i\omega _n)R[G(i\omega _n)],$$
(4)
where $`G(i\omega _n)`$ is the site-local (impurity) Green function, and is itself a functional of $`𝒢_0(i\omega _n)`$. $`\omega _n=2\pi T(n+1/2)`$ are the Fermionic Matsubara frequencies. $`R[x]`$ is the reciprocal function of the Hilbert transform, defined as $`R[\stackrel{~}{D}(\xi )]=\xi `$. The Hilbert transform of the non-interacting density of states, $`𝒟(\epsilon )`$, is defined as $`\stackrel{~}{D}(\xi )_{\mathrm{}}^+\mathrm{}𝑑\epsilon 𝒟(\epsilon )/(\xi \epsilon )`$.

When electrons move on a tight-binding hypercubic lattice, the bare DOS takes the form of a Gaussian Metzner and Vollhardt (1989), $`𝒟(\epsilon )=\mathrm{exp}(\epsilon ^2/2t^2)/t\sqrt{2\pi }`$ and there is no simple expression for the reciprocal function. Introducing the modified Dyson equation,

$$𝒢_0^1(i\omega _n)=G^1(i\omega _n)+\mathrm{\Sigma }(i\omega _n),$$
(5)
(where $`\mathrm{\Sigma }(i\omega _n)`$ is the electron self-energy), equation (4) may be rewritten as,

$$R[G(i\omega _n)]=i\omega _n+\mu \mathrm{\Sigma }(i\omega _n),$$
(6)
and inverted to give an expression for the Green function in terms of the self-energy,

$$G(i\omega _n)=\frac{d\epsilon 𝒟(\epsilon )}{i\omega _n+\mu \mathrm{\Sigma }(i\omega _n)\epsilon }.$$
(7)
This also allows the approximation to be interpreted as a course graining of the momentum space Hettler et al. (2000). To complete the scheme, a form for the electronic self-energy must be calculated in terms of the non-interacting Weiss field. This is normally approximate. Then a self-consistent procedure is followed: Compute the Green function from equation (7), the Weiss field from equation (5) and then re-calculate the self-energy until convergence is reached.

The application of ME theory within DMFT corresponds to computing the self-energy and phonon propagator from the diagrammatic equations in figure 1(a and b). The self-consistent solution of such a self-energy corresponds to the summation of all Feynman diagrams which contain no vertex corrections. In the low (non-zero) phonon-frequency limit, Migdal’s analysis indicates a condition, $`U\omega _0t^2`$, for the neglect of corrections to the vertex function Migdal (1958).

The solution of the Dyson equation seen in figure 1(a) results in the full phonon propagator,

$$D(\omega )=\frac{\omega _0^2}{(\omega +i\eta )^2\omega _0^2[1+\mathrm{\Pi }_0(\omega )]},$$
(8)
which is then used for the calculation of the electron self-energy. The phonon polarisation bubble (dimensionless self-energy) $`\mathrm{\Pi }_0(i\omega _s)`$ is given from perturbation theory,

$$\mathrm{\Pi }_0(i\omega _s)=2UT\underset{m}{}G(i\omega _m)G(i\omega _s+i\omega _m),$$
(9)
and may be analytically continued by introducing the spectral representation, $`G(i\omega _n)=𝑑x\rho (x)/(i\omega _nx)`$ where $`\rho (\omega )=\mathrm{Im}[G(\omega +i\eta )]/\pi `$ and performing the sum over Matsubara frequencies to give,

$$\mathrm{Im}[\mathrm{\Pi }_0(\omega )]=2U_0^\omega 𝑑x\rho (x)\rho (\omega x),$$
(10)
at absolute zero, and

$$\mathrm{Im}[\mathrm{\Pi }_0(\omega ,T)]=U\mathrm{sinh}(\frac{\omega }{2T})_{\mathrm{}}^{\mathrm{}}𝑑x\frac{\rho (x)}{\mathrm{cosh}(\frac{x}{2T})}\frac{\rho (\omega x)}{\mathrm{cosh}(\frac{\omega x}{2T})},$$
(11)
at finite $`T`$. The full spectral function is used so that all diagrams with no vertex corrections are included, consistent with the proper interpretation of Migdal’s analysis.

The lowest-order skeleton diagram shown in figure 1(b),

$$\mathrm{\Sigma }(i\omega _n)=UT\underset{s}{}G(i\omega _ni\omega _s)D(i\omega _s),$$
(12)
may be analytically continued in the same way as the polarisation bubble to give,

$$\mathrm{Im}[\mathrm{\Sigma }(\omega )]=U_0^\omega 𝑑x\rho (x)\sigma (\omega x),$$
(13)
with $`\sigma (\omega )=\mathrm{Im}[D(\omega +i\eta )]/\pi `$. A more complicated expression applies at finite temperature,

$$\mathrm{Im}[\mathrm{\Sigma }(\omega )]=U\mathrm{cosh}(\frac{\omega }{2T})_{\mathrm{}}^{\mathrm{}}\frac{dx\rho (\omega x)\sigma (x)}{2\mathrm{cosh}(\frac{\omega x}{2T})\mathrm{sinh}(\frac{x}{2T})}.$$
(14)
At each stage, the Kramers–Kronig relation is used to compute the real parts of the electron and phonon self-energies, for instance $`\mathrm{Re}[\mathrm{\Pi }_0(\omega )]=\mathrm{P}𝑑x\mathrm{Im}[\mathrm{\Pi }_0(x)]/\pi (\omega x)`$, where P denotes the principal integral.

After the self-consistent procedure has converged, physical properties are calculated, including the quasi-particle weight (inverse effective mass),

$$Z=\frac{m_0}{m^{}}=1\frac{\mathrm{\Sigma }}{\omega }|_{\omega =0},$$
(15)
the effective phonon frequency $`\mathrm{\Omega }`$ from,

$$\mathrm{\Omega }^2\omega _0^2(1+\mathrm{Re}[\mathrm{\Pi }_0(\mathrm{\Omega })])=0,$$
(16)
and the optical conductivity Pruschke et al. (1993),

$$\mathrm{Re}[\sigma (\omega )]=\frac{\pi }{\omega }_{\mathrm{}}^{\mathrm{}}𝑑\epsilon 𝒟(\epsilon )_{\mathrm{}}^{\mathrm{}}𝑑\nu \rho (\epsilon ,\nu )\rho (\epsilon ,\nu +\omega )[f(\nu )f(\nu +\omega )],$$
(17)
where $`f(x)`$ is the Fermi-Dirac distribution and $`\rho (ϵ,\nu )=\mathrm{Im}[1/(\nu +i\eta ϵ\mathrm{\Sigma }(\nu ))]/\pi `$ (taking the limit, $`\omega d\omega `$, the DC conductivity is recovered.)

Although DMFT is approximate, the formalism can be expected to work in the non-interacting limit, where course graining will give the exact non-interacting DOS for the tight-binding model (regardless of dimension). In the opposite limit ($`t0`$), the neglect of loops through the host is justified, and the formalism should be exact. Such propagation of correlation is also small when the coordination number is high (e.g. FCC lattices).

## III Spectral properties

We have solved the DMFT equations using the self-energies in equations (13) and (14) to find electron spectral functions. We show the evolution of these functions with coupling for $`\omega _0=0.125`$ in figure 2(a) and $`\omega _0=0.5`$ in figure 2(b). The spectral functions show two features. There is a low energy peak with a width defined by the phonon frequency, and a high energy shoulder. This behaviour can be related to the two limits of the Holstein model. At low energy scales ($`\omega <\omega _0`$), electrons interact via virtual processes and the behaviour is essentially Hubbard-like (correlated). At frequencies greater than $`\omega _0`$, the spectral weight is reduced, phonons may be created, and static behaviour emerges (in the static limit, phonons may always be created as no energy is required to deform the lattice). The central peak narrows with increased coupling, until a critical value is reached. Here, the theory has a Brinkman–Rice-like transition with a diverging effective mass Brinkman and Rice (1970).

In figure 3, we plot the imaginary part of the self-energy. For low coupling strengths (small $`U`$ and energy scales ($`|\omega |<\omega _0`$), $`\mathrm{Im}[\mathrm{\Sigma }(\omega )]`$ is small because electrons cannot create phonons, consistent with the Englesberg–Schrieffer analysis). However, as the coupling increases, electrons can be scattered by bipolaron resonances and $`\mathrm{Im}[\mathrm{\Sigma }(\omega )]`$ rises sharply on either side of the Fermi-energy. The gap in the weight of the self-energy is frequently assumed in the analysis of ARPES data, but as we can see here, when the electron-phonon coupling approaches the band width (i.e. $`\lambda 1`$), the self-energy has a more complicated form, without the simple picture of coherent and incoherent quasi-particles. Indeed, it has been demonstrated that in low dimensions, vertex corrections are required to achieve a sharp discrimination between coherent and incoherent dressed electrons Hague (2003).

The onset of this “transition” can be studied by examining the inverse quasi-particle mass (quasi-particle weight), shown in figure 4 for two values of phonon frequency. As the coupling increases, $`Z`$ becomes smaller and eventually vanishes at a critical value of coupling (the effective mass diverges).

We also examine the phonon spectral function in figure 5 for $`\omega _0=0.125`$ and various couplings. As the coupling is increased, the phonon modes soften (figure 6). The effective frequency does not tend to zero as quickly as the quasi-particle weight. We will revisit this point later in this article.

It is clearly the case that ME theory does not correctly describe the strong coupling limit. It is known from the exact solution of the static limit Millis et al. (1996) and from approximate “iterated-perturbation theory” Georges and Kotliar (1992) and QMC Jarrell (1992) solutions of the Hubbard model, that sub-bands should form at strong coupling. As our calculations demonstrate, this sub-band formation is not properly reproduced within ME theory and this makes it likely that effects of vertex corrections become significantly more important at strong coupling. Since this is at variance with the traditional interpretation of Migdal’s analysis, we are motivated to re-examine the vertex function in the next section.

## IV Breakdown of Migdal–Eliashberg theory

The neglect of the lowest-order vertex correction shown in figure 1(c) and all higher-order corrections is central to Migdal–Eliashberg theory. In this section, we compute the lowest-order correction and use this to define the region of validity for ME theory. In this sense, we are revisiting Migdal’s analysis to understand why there is a contradiction between low frequency results computed with ME theory and advanced numerical methods.

The ratio of first to zeroth order vertices at finite temperature may be written as

$$\frac{\mathrm{\Gamma }_1(i\omega _n,i\omega _x)}{\mathrm{\Gamma }_0}=TU\underset{s}{}𝒢_0(i\omega _ni\omega _s)𝒢_0(i\omega _ni\omega _si\omega _x)D_0(i\omega _s).$$
(18)
As temperature tends to zero, this sum becomes an integral which may easily be evaluated ($`x`$ and $`\omega `$ are defined as in figure 1),

$$\frac{\mathrm{\Gamma }_1(i\omega ,ix)}{\mathrm{\Gamma }_0}=\frac{U}{2\pi }𝑑s𝒢_0(i\omega is)𝒢_0(i\omega ixis)D_0(is).$$
(19)
In order for vertex corrections to be unimportant, the ratio $`\mathrm{\Gamma }_1/\mathrm{\Gamma }_0`$ must be small. As an example, we show the lowest order correction to the vertex function for (a) $`\omega _0=0.0556`$ and $`U=0.214`$ (b) $`\omega _0=0.75`$ and $`U=0.214`$ in figure 7. For the larger phonon frequencies, the ratio is clearly greater than 10% and the approximation is expected to be significantly changed by the inclusion of vertex corrections. For the weak coupling $`U`$ shown here, the corrections are much less pronounced for small phonon frequencies.

In Figure 8 we show the magnitude of the central peak of the vertex correction (always the largest part of the function) for a range of couplings at $`\omega _0=0.0556`$. Between $`U=0.54`$ and $`U=0.55`$ the vertex correction diverges, and ME theory clearly breaks down. As it is necessary to pass through the divergence to reach $`U>0.545`$, the theory is not applicable for strong couplings.

Figure 9 shows the lowest value of $`U^{}=U/(t+U)`$ for given $`\omega _0^{}=\omega _0/(t+\omega _0)`$ at which the ratio $`\mathrm{\Gamma }_1/\mathrm{\Gamma }_0`$ exceeds 0.1 (the primed quantities are chosen so that $`U^{}=1`$ represents an unprimed value $`U=\mathrm{}`$, and a primed value of $`U^{}=0.5`$ represents an unprimed value $`U=t`$). Below this value ME theory can be expected to give accurate results, while the theory begins to break down above this line. For $`\omega _0^{}0.2`$ ($`\omega _00.125`$) the standard Migdal criterion works well and correctly predicts the region in which ME theory is applicable. However for frequencies $`\omega _0^{}0.2`$ ($`\omega _00.125`$), the Migdal criterion misses the divergence in the vertex correction associated with the divergence in the effective mass. This breakdown happens at $`Ut`$ or $`\lambda 1`$. The line of breakdown due to this vertex divergence levels off as the phonon frequency approaches zero and does not scale as $`1/U`$ (Migdal’s criterion). These results imply that there is a breakdown in ME theory in the adiabatic limit.

We note that at high phonon frequencies ($`U,t\omega _0`$), our line does not tend to zero. This is because the contribution $`\mathrm{\Gamma }_1/\mathrm{\Gamma }_0(x0,\omega ^{}0)`$ does not tend to infinity. Instead, the function gets wider as the phonon frequency increases (see figure 7, panel b), and some combination of the magnitude and width of the vertex is probably a better measure of the importance of vertex corrections. At large $`\omega _0`$, the breakdown is quite different to that at low phonon frequencies as the second order diagrams (one of which contains a vertex correction) develop the same functional form and magnitude. However the breakdown at $`U=0.1`$ ($`U^{}=0.17`$) is not too bad an estimate. At this value of $`U`$, there is only a small mass renormalisation due to vertex corrections.

Our results are consistent with that of Benedetti and Zeyher Benedetti and Zeyher (1998), which was obtained with a different technique using a semi-circular density of states ($`𝒟_{SC}(ϵ)=\sqrt{4t^2ϵ^2}/2\pi t^2`$), and predicted breakdown for $`\lambda _{SC}1.25`$ in the extreme adiabatic limit. We have included their result in the figure, choosing the bandwidth parameter by matching the bare DOS at the Fermi energy for the two cases (the conversion factor is $`\lambda =𝒟_{SC}(0)\lambda _{SC}/𝒟(0)=\lambda _{SC}\sqrt{2/\pi }`$) which is reasonable when considering low energy excitations at half-filling, and corresponds to a critical $`\lambda 0.997`$ for the Gaussian DOS used here (see the diamond in figure 9). Millis et al. Millis et al. (1996) also predicted a breakdown in ME theory at $`Ut`$ from their solution of the static limit of the infinite-dimensional Holstein model, consistent with the breakdown that we have found. The solution of the static limit predicts the formation of sub-bands, which is something not reproduced in the ME solution. The lack of subbands in ME theory demonstrates that higher order diagrams are essential for the description of the sub-band formation, and therefore for the description of the strong coupling limit.

We suggest that the breakdown of ME theory should be understood in the following heuristic manner. As the coupling increases, the host spectral function ($`𝒢_0`$) narrows and spectral weight is shifted towards lower energy scales. There is then an effective bandwidth, $`t_{\mathrm{eff}}<t`$, for quasi-particles close to the Fermi surface. The condition for the applicability of ME theory then becomes $`U\mathrm{\Omega }t_{\mathrm{eff}}^2Z^2`$. As seen in figure 5, there is also a reduction in the effective phonon frequency with increasing interaction strength. Although the renormalisation of the phonon frequency helps to drive against the transition, we note that the band narrowing effect is much stronger, as can be seen by comparing figures 4 and 6.

## V Resistivity Curves

Using the finite temperature form of the self-energy, it is possible to calculate resistivity curves. These are shown in figure 10 at various electron-phonon coupling and phonon frequency. At high temperatures, the tendency is to linear behavior, with a gradient depending on $`U`$, but independent of $`\omega _0`$. The effect of phonon frequency is most dramatically seen at $`T\omega _0/2`$, where a point of inflection is seen, before the curve tends to low temperature $`T^2`$ behavior, consistent with a weakly renormalized electron gas. We note that long range order has not been considered in this study, so no superconducting transition is seen. The negative curvature of the resistivity curve at intermediate temperatures is of interest, since that curvature has recently been interpreted as an onset of the Mott limit, but may have a different interpretation as an intermediate coupling phenomenon Lortz et al. (2005).

## VI Conclusions

We have computed spectral functions, resistivity curves, self-energies, quasi-particle weight and effective phonon frequency for the Holstein model at half-filling using Migdal–Eliashberg theory, within the DMFT framework, and revisited Migdal’s analysis of the strength of the vertex correction. By analysing the first order vertex correction, we have defined the region of applicability of ME theory. We find that ME theory breaks down at intermediate coupling in the adiabatic limit, showing that Migdal–Eliashberg theory should only be trusted at weak coupling within this framework. The coupling at which breakdown occurs corresponds to a divergence in the effective mass, indicating that both long range order and vertex corrections should be included to correctly describe the strong coupling regime.

The question remains: What was lacking from the analysis of Migdal, which estimated the magnitude of the vertex corrections. This can be understood from the magnitude of the kinetic energy. In the calculations, the kinetic energy can be seen to be decreasing due to the localisation at the non-analyticities. Thus in Migdal’s estimate, the renormalised kinetic energy, and not the bare kinetic energy should have featured. The renormalised phonon frequency should also have been used, but it tends to zero more slowly than the kinetic energy. At the non-analyticity, where the kinetic energy tends to zero, it can be seen that a vertex divergence is expected from the modified version of Migdal’s estimate.

Our results are timely because of a recent upsurge of interest in the role of electron-phonon coupling in the cuprate superconductors. There are many experimental results which have been analysed and interpreted using ES/ME theory, in spite of the large coupling constant that has been predicted from that analysis. Many groups are using Migdal’s analysis to justify their claims, without a careful consideration of the self-consistent effect of coupling on the kinetic energy of the electrons (polarons). Exact numerical analysis shows that large electron-phonon couplings significantly reduce the kinetic energy of electrons Hague et al. (2006) and thus affect the internal consistency of Migdal’s analysis. Therefore, the any very large $`\lambda `$ determined from ES style analysis of experimental results is not consistent with a theory that neglects vertex corrections (i.e. a more detailed analysis with higher order effects included is necessary). Without that analysis, the only conclusion that can be reached is that the coupling is large, but no reliable value for that coupling can be determined. The main conclusion of a breakdown is not expected to change in 1D, 2D or 3D. In three dimensional systems, the role of spatial fluctuations should only be significant very close to the bi-polaron instability. The role of spatial fluctuations on the spectral properties in 2D has been analysed at weak coupling, and is found to lead to quantitative differences in the coherent excitations Hague (2003). Spatial fluctuations have a major role in the superconducting state, and it has also been shown that the Eliashberg approach to superconductivity is inadequate for optical phonon mediated $`d`$-wave superconductivity Hague (2006), and incomplete for $`s`$-wave superconductivity Hague (2005). We therefore urge researchers to carefully analyse the internal consistency of their theories when dealing with strong coupling materials.

## VII Acknowledgements

The authors would like to thank F.Gebhard and F.Essler for useful discussions. JPH thanks the EPSRC for partial funding of this work.


text:
# Quantum dots in magnetic fields: thermal response of broken symmetry phases

## Abstract

We investigate the thermal properties of circular semiconductor quantum dots in high magnetic fields using finite temperature Hartree-Fock techniques. We demonstrate that for a given magnetic field strength quantum dots undergo various shape phase transitions as a function of temperature, and we outline possible observable consequences.

The quantum mechanical characteristics (i.e. ground- and excited-state energy levels) of semiconductor quantum dots were recently investigated as a function of increasing magnetic field and electron number . Unique structures that break rotational symmetry within a quantum dot and develop as a function of increasing magnetic field were also reported in theoretical investigations  and associated with experimental observations of charge-density redistribution . In this paper, we explore the thermal characteristics of the broken-symmetry phase of the charge-density distribution, showing that temperature can be used to induce transitions between phases of differing spatial symmetry.

Thermal properties of, and transitions in, strongly correlated quantum mechanical systems have been studied in various contexts. Perhaps the most obvious are the normal-conducting to super-conducting transitions in condensed matter systems. Many nuclear systems exhibit an intrinsic deformation in their ground states, and deformed-to-spherical transitions have been studied in many models (see for example Refs. ). As a result, observable consequences, such as enhanced moments of inertia, decreasing quadrupole transition moments, and decreasing correlated pair-transfer amplitudes may be investigated.

In studies of many-body phenomena in quantum dots, experimental efforts have focused on mapping the magnetic field dependence of their ground-state structure by measuring the chemical potential via capacitance spectroscopy . Cusps and steps in the chemical potential were found to clearly separate different ranges of magnetic fields . These features were identified with phase transitions in the charge density of the quantum dot. At magnetic field strengths on the order of a few tesla, all electrons become spin-polarized initiating the maximum density droplet (MDD) phase , in which the density is constant and homogeneous at the maximum value that can be reached in the lowest Landau level. The stability of the MDD is determined by a competition among the kinetic energy, external confinement, the Coulomb repulsion between electrons, and the attraction created by the Coulomb exchange term. For increasing magnetic field, the charge-density distribution of the droplet reconstructs  with a ring of electrons breaking off from the MDD phase. This edge reconstruction has been shown via mean-field  and density functional theory  calculations to result from a rotational symmetry-breaking phase transition from the MDD to a Wigner molecule or Wigner crystal phase. These calculations are in good qualitative agreement with recent experimental results where instabilities of the MDD state and other transitions in the high magnetic field region were accompanied by a redistribution of the charge density .

Most of the experimental investigation of many-body phenomena in quantum dots was performed at very low temperatures on the order of 100 mK and focused on the evolution of the ground state of the quantum dot as a function of increasing magnetic field and electron number . The thermal characteristics of these systems are largely unexplored, especially for large values of the magnetic field. We explore the thermal response of the broken-symmetry phase at constant magnetic field and electron number and demonstrate that transitions among a number of distinct, broken-symmetry phases occur as the temperature is increased.

In describing the ground-state and low-lying (intra-band) excitations of the $`N`$-electron semiconductor nanostructures, it is often sufficient to restrict consideration to the conduction band using the effective-mass approximation . We consider the problem of N electrons of effective mass $`m^{}`$ in a plane, $`(x,y)`$, confined by an external parabolic potential, $`V(r)=\frac{1}{2}m^{}\omega _0^2r^2`$, and subject to a strong magnetic field $`\stackrel{}{B}=B_0\stackrel{}{e}_z`$. We consider the Zeeman splitting but neglect the spin-orbit interaction. The Hamiltonian for such a system is

$`\widehat{H}`$ $`=`$ $`{\displaystyle \underset{i}{}}\left[{\displaystyle \frac{\left(\stackrel{}{p}\frac{e}{c}\stackrel{}{A}\right)^2}{2m^{}}}+V(r_i)+{\displaystyle \frac{g^{}\mu _B}{\mathrm{}}}\stackrel{}{B}\stackrel{}{S}_i\right]`$ (1)
$`+`$ $`{\displaystyle \underset{i<j}{}}{\displaystyle \frac{e^2}{\epsilon |\stackrel{}{r}_i\stackrel{}{r}_j|}},`$ (2)
where the vector potential is $`\stackrel{}{A}(\stackrel{}{r}_i)=(B_0/2)(y_i,x_i,0)`$, $`g^{}=0.54`$, $`\epsilon =12.9`$, $`m^{}=0.067m_e`$, and $`\mathrm{}\omega _0=3`$ meV.

We solve this equation at the finite-temperature Hartree–Fock level. This is the initial step to account for electron-electron correlations within the quantum dot and suffices for a qualitative discussion of phase transitions and thermal behavior. Hartree-Fock is implemented across nuclear and atomic physics as a first step towards the solution of the quantum many-body problem . Since the Hartree-Fock equations are nonlinear, the self-consistent potential obtained from their solution will not necessarily show the same symmetries as the original Hamiltonian. Of course, the exact wave will maintain the symmetries of the original Hamiltonian. As a point of reference, an estimate of the relative uncertainty in total angular momentum at high magnetic field strengths was shown to be roughly $`10\%`$ by the authors of Ref. .

The equations describing the static mean field (i.e., Hartree-Fock) at finite temperature for an $`N`$-electron system are

$$[\widehat{K}+\widehat{W}(\beta )E_\alpha ]|\varphi _\alpha (\beta )=0$$
(3)
where $`E_\alpha `$ are the single particle energies associated with the single-particle wave functions $`\varphi _\alpha `$, and $`\widehat{K}`$ is the one-body kinetic-energy operator. The one-body Hartree–Fock field, $`\widehat{W}(\beta )`$, is obtained in terms of the Hamiltonian (2) and the one-body density matrix, $`\widehat{\rho }(\beta )`$ as,

$$\widehat{W}(\beta )=\mathrm{Tr}[\widehat{\rho }(\beta )\widehat{H}].$$
(4)
The single-particle wave functions are denoted as $`|\varphi _\alpha (\beta )`$, and $`\beta =1/kT`$ is the inverse temperature. The representation of the one-body density is given in terms of the occupation numbers $`n_\alpha (\beta )`$, $`\rho (\beta )=_\alpha n_\alpha (\beta )\varphi _\alpha (\beta )\varphi _\alpha (\beta )`$, and the Fermi–Dirac occupation numbers at chemical potential $`\mu `$ are $`n_\alpha (\beta )=1/\left[e^{\beta (E_\alpha \mu )}+1\right]`$. Note that Tamura and Ueda  studied the number fluctuations $`\left(\delta N\right)^2`$ of a quantum dot as a function of $`N`$ for different field strengths and found that the fluctuations were rather small (on the order of 0.1) in finite-temperature Hartree–Fock applications. We expect qualitative agreement with experiment as Hartree–Fock is known to overestimate the magnetic field strength at which transitions occur . However, for a fixed field strength, experimental evidence of the thermal transitions should still be qualitatively visible.

We use a Fock–Darwin basis expansion to solve the finite-temperature Hartree-Fock equations. Since we use a high ($`12`$ T) magnetic field, we consider only angular momentum states with the $`n=0`$ principal quantum number. The electrons carry spin, and so our states are labeled by $`k=\{l_k,s_k\}`$, where $`l_k`$ is the angular momentum projection of the $`k`$-th state and $`s_k`$ is the spin of that state. We found convergence using fifty states for the $`N8`$ systems. We also checked our zero-temperature results with other publications  for various numbers of electrons in the dot and found satisfactory agreement.

We begin the discussion of our results by investigating the electron charge, angular momentum, and spin densities as a function of increasing temperature for the $`N=6`$ system at $`B_0=12.15`$ T. We show densities at representative temperatures of 3.87 K (the low temperature limit, $`\nu =6`$), 11.97 K (before the first phase transition, $`\nu =6`$), 13.65 K (in the second phase, $`\nu =5`$), and 14.32 K (in the third phase, $`\nu =4`$), where $`\nu `$ is the number of definable high-density regions (or vortices) in the charge density plots. The density, shown in Fig. 1a-d, begins as a fairly well-defined Wigner crystal at 3.87 K, which exhibits some degree of thermal broadening at 11.97 K. The $`\nu =5`$ and $`\nu =4`$ phases continue to show a similar amount of density in the remaining vortices, while the density of the thermally dissipated vortices have effectively spread through the entire dot. The angular momentum along the $`B`$-field direction, shown in Fig. 1e-h, exhibits well-defined structures at low temperatures which tend to decrease rapidly as one moves through the various phases. Although the high charge density regions in the $`\nu =4`$ phase are still well defined, the angular momentum in this high-temperature phase has nearly washed out. Finally, we show in Fig. 1i-l the spin density defined as $`\rho _s(x,y)=\left[\rho _{}(x,y)\rho _{}(x,y)\right]/\left[\rho _{}(x,y)+\rho _{}(x,y)\right]`$, where $`\rho _{}`$ ($`\rho _{}`$) refer to the spin up (down) density. At these temperatures, little appreciable spin depolarization occurs and the spin density remains above 0.8 for the entire region where there is appreciable charge density.

The suddenness of the phase transitions seen in Fig. 1 become quite evident when the internal energy of the quantum dot is plotted as a function of the temperature. We show the three phases of the dot in Fig. 2a. Note that the $`\nu =6`$ phase exists as an excited configuration when the most probable Hartree–Fock solution is the $`\nu =5`$ phase. Similarly, the $`\nu =4`$ configuration exists as a possible excited configuration of the system even at fairly low temperatures. The specific heat, $`C_v=dH/dT`$ (with T in units of eV), is shown in Fig. 2b. Clearly, when the energy undergoes a phase transition, the specific heat shows a sharp structure. This occurs since the energy is piece-wise continuous along the three phases.

These calculations suggest that the quantum dot exhibits a band structure of many-body levels. The low-temperature states all have the same intrinsic shape characteristic (the same vortex structure). As we increase the temperature, other $`\nu `$ phases become accessible. At the point when two bands of different intrinsic character cross in energy, we find a phase transition. Similar phenomena are found in nuclear physics, where at higher nuclear excitation energies the eigenstates of the system may be of a different intrinsic deformation when compared to states of the ground-state band .

We also find generally that the exchange energy decreases as a function of increasing temperature but with differing slopes in the different phases. In the region of phase transitions, the ratio of the exchange energy to the direct energy decreases from 0.38 (at 11.5 K)to 0.32 (at 13 K) in the $`\nu =6`$ phase. In the $`\nu =5`$ phase, this ratio decreases from 0.32 at $`T=13`$ K to 0.28 at $`T=13.8`$ K, and from 0.28 at $`T=13.8`$ K to 0.26 at $`T=14.5`$ K in the $`\nu =4`$ phase. The slopes of the decreasing ratio are different across the three phases, with $`\nu =6,5`$ being the largest and $`\nu =4`$ more gradual.

The occupation probabilities of the Fock–Darwin states, $`n_{FD}`$, for the temperature conditions of Fig. 1 are shown in Fig. 2c. In the low-temperature phase (3.9 K, solid line), the $`l=0`$ state is occupied, while an edge reconstruction has occurred for the remaining five electrons. As we increase the temperature to 12 K (dotted line), we see a spreading of the occupations in both the low and high angular momentum channels with about 0.8 particles in the $`l=0`$ state. The $`\nu =5`$ phase brings a dramatic decrease of occupation in the $`l=0`$ state and a shift to lower angular momentum for the reconstructed edge. This trend continues after the $`\nu =4`$ transition.

Signatures of the density transitions that we have seen also appear in the Hartree–Fock occupations $`n_{HF}`$ as shown in Fig. 2d. At low temperatures (3.9 K) the familiar step-function behavior is evident from the figure. At 12 K in the $`\nu =6`$ phase, we see a decrease of occupation to roughly 0.8 in the lowest six Hartree–Fock levels and a spreading to higher energy states. As the system undergoes the transition to $`\nu =5`$, we see only five Hartree–Fock levels significantly filled (with $`n_{HF}>0.7`$), and finally in the $`\nu =4`$ phase, only four Hartree-Fock levels are significantly filled. The occupation number-spreading, which is due to thermal excitation of the system, is enhanced significantly when the system undergoes a phase transition.

The phase transitions that we have shown in the preceding discussion have definite observable consequences. In Fig. 3a we plot the chemical potential, that is, the separation energy $`\mathrm{\Delta }(N,T)=E_N(T)E_{N1}(T)`$ to remove a particle from the quantum dot at a given temperature. Since this is an energy difference, $`\mathrm{\Delta }(N,T)`$ will be influenced by transitions within both the $`N=6`$ and $`N=5`$ systems, and we expect changes in slope at the transition points. The $`\nu =5`$ to $`\nu =4`$ transition in the $`N=5`$ dot occurs at roughly 12 K, causing a sharp rise in $`\mathrm{\Delta }(N=6)`$. A slope change in $`\mathrm{\Delta }(N=6)`$ is seen at $`13`$ K, where the $`\nu =6`$ to $`\nu =5`$ transition occurs in the $`N=6`$ system. The decrease from $`1414.3`$ K occurs when the $`N=6`$ system makes the transition from $`\nu =5`$ to $`\nu =4`$. A final change in slope occurs when the $`N=5`$ dot makes the transition from $`\nu =4`$ to $`\nu =3`$. For the brief temperature interval when the $`N=6`$ and $`N=5`$ dots are in the same $`\nu `$ phase, we see a decrease in the chemical potential.

We observe similar changes in the inverse compressibility, $`\mathrm{\Delta }_2(N,T)=E_{N+1}(T)2E_N(T)+E_{N1}(T)`$. This quantity has been measured for quantum dots in low magnetic fields  and studied in Hartree–Fock theory for ground-state properties . In our case, the $`N=7,6`$, and $`5`$ dots participate in the observable. We again notice strong effects as one passes through transition points in either of the three systems contributing to the observable. Figure 3b shows $`\mathrm{\Delta }_2(N=6)`$ as a function of temperature. Before transitions occur, $`\mathrm{\Delta }_2`$ remains fairly constant. A large decrease begins at 12 K, where the $`N=5`$ system undergoes its first transition. Interestingly, $`\mathrm{\Delta }_2`$ increases significantly when the $`N=6`$ and $`N=5`$ dots are in the $`\nu =4`$ phase.

In order to investigate the sensitivity of our results to the approximations made (namely the Hartree-Fock approximation) we extended our studies to include the second-order perturbative correction to the Hartree-Fock energy. This correction to the energy is given by

$$\mathrm{\Delta }E_2=\frac{1}{4}\underset{abϵ_F}{}\underset{rs>ϵ_F}{}\frac{ab\overline{v}rs^2}{\epsilon _a+\epsilon _b\epsilon _r\epsilon _s}$$
(5)
where we restrict the sums to be below and above the Fermi energy surface, $`ϵ_F`$. In this expression, the Hartree-Fock states are given by $`a,b,r,s`$, with associated single-particle energies $`\epsilon _a,\epsilon _b,\epsilon _r,\epsilon _s`$, and $`ab\overline{v}rs`$ are the antisymmetrized two-body Hartree-Fock matrix elements of the original two-body interaction of the Hamiltonian. We note that the absolute value of $`\mathrm{\Delta }E_2`$ in our calculations is always less than $`1\%`$ of the total energy, and in fact only slightly changes the observable quantities, as is evident from Fig. 3, where the second-order perturbative results are also shown for $`\mathrm{\Delta }`$ and $`\mathrm{\Delta }_2`$.

In addition to measurements of energy differences, one should be able to experimentally probe the thermal phase transitions using far-infrared spectroscopy and X-ray scattering. Far-infrared spectroscopy was used to investigate the excitations of InAs quantum dots as a function of the electron number per dot . Grazing incidence X-ray scattering was recently used to generate a full structural characterization of quantum dots, including information on the elastic form factor. While these experiments were carried out at very low temperatures, it is conceivable that one could study the thermal response of quantum dots using X-ray or far-infrared scattering. A Fourier transform of the charge density produced in our calculations gives the elastic form factor that can be used to characterize the dot. We show in Fig. 4 the form factor, $`\rho (q)`$, as a function of the momentum transfer vector $`q=\sqrt{q_x^2+q_y^2}`$. A well-defined minimum is apparent at approximately $`q=0.7`$ev<sup>-1</sup> in all three phases. This first minimum is related to the size of the dot and clearly does not change in the three phases. This is apparent also from a close inspection of Fig. 1: while the internal structure changes significantly as a function of increasing temperature, the dot size does not change. The height of the second maximum slightly increases as in phase $`\nu =5`$ before decreasing in the $`\nu =4`$ phase. The position of the second minimum increases as a function of $`q`$ significantly and would be a distinguishing feature in X-ray or far-infrared scattering experiments used to probe the thermal phases of a quantum dot.

We see in the above calculations that quantum dots exposed to high magnetic fields undergo various phase transitions as the temperature increases. Indeed, we have followed the $`N=6`$ system through the $`\nu =3`$ and $`2`$ phases as well. We have also studied these transitions at various magnetic-field strengths, and with differing numbers of electrons in the dot up to $`N=20`$, with the results presented here being characteristic of the general behavior. Our purpose in this paper has been to demonstrate that thermal excitations of a quantum dot induce phase transitions within the dot that should have observable consequences. Thus, we have shown that temperature may also be used as a control to tune the characteristics of quantum dots. Inclusion of the full many-body correlations may slightly dampen some of these effects , but signatures of the phase transitions should still appear in scattering experiments and measurements of the inverse compressibility and chemical potential of the quantum dot.

This research was sponsored by the Laboratory Directed Research and Development Program of Oak Ridge National Laboratory, under Contract No. DE-AC05-00OR22725 managed by UT-Battelle, LLC. by the U.S. Department of Energy.