text:
# Finite dimensional integrable Hamiltonian systems associated with DSI equation by Bargmann constraints

## 1 Introduction

The Davey-Stewartson I (DSI) equation is a famous 2+1 dimensional integrable equation which describes the motion of water wave . This equation has localized soliton solutions and has been studied in various ways, such as inverse scattering , binary Darboux transformation , nonlinearization to 1+1 dimensional problems  etc.

For 1+1 dimensional integrable systems, the nonlinearization procedure, both mono-nonlinearization  and binary nonlinearization , reduces them to finite dimensional (1+0 dimensional) integrable Hamiltonian systems . Therefore, it transforms a partial differential equation to a systems of ordinary differential equations. This greatly simplifies the procedure of getting solutions, at least numerical solutions. Some important explicit solutions, especially periodic or quasi-periodic solutions have been obtained in this way.

The nonlinear constraint method has also been applied to some 2+1 dimensional equations like the KP, MKP, N-wave equations etc. . For the DSI equation, we have already found its new Lax system (1) by nonlinearization in which all the derivatives are separated. Each pair of equations in this system is 1+1 dimensional. Hence the derived Lax system is looked as essentially 1+1 dimensional because we can use 1+1 dimensional method to solve it. It is possible to nonlinearize this essentially 1+1 dimensional system again to get finite dimensional Hamiltonian systems.

In the present paper, we show that there are Bargmann constraints which reduce the DSI equation to finite dimensional Hamiltonian systems. These Hamiltonian systems have a full set of integrals of motions and these integrals of motion are functionally independent in a dense open subset of the phase space. Therefore, these Hamiltonian systems are completely integrable in Liouville sense.

## 2 Nonlinearization

We consider the following Lax system

$$\begin{array}{c}\mathrm{\Phi }_x=W^x\mathrm{\Phi }=\left(\begin{array}{ccc}\text{i}\lambda & 0& \text{i}f\\ 0& \text{i}\lambda & \text{i}g\\ \text{i}\overline{f}& \text{i}\overline{g}& 0\end{array}\right)\mathrm{\Phi }\hfill \\ \mathrm{\Phi }_y=W^y\mathrm{\Phi }=\left(\begin{array}{ccc}\text{i}\lambda & u& \text{i}f\\ \overline{u}& \text{i}\lambda & \text{i}g\\ \text{i}\overline{f}& \text{i}\overline{g}& 0\end{array}\right)\mathrm{\Phi }\hfill \\ \mathrm{\Phi }_t=W^t\mathrm{\Phi }\hfill \\ =(\begin{array}{c}2\text{i}\lambda ^2+\text{i}|u|^2+\text{i}v_1\hfill \\ 2\overline{u}\lambda +\text{i}\overline{u}_y\hfill \\ 2\text{i}\overline{f}\lambda +2\overline{f}_y\hfill \end{array}\begin{array}{c}2u\lambda +\text{i}u_y\\ 2\text{i}\lambda ^2\text{i}|u|^2\text{i}v_2\\ 2\text{i}\overline{g}\lambda +2\overline{g}_y\end{array}\begin{array}{c}\hfill 2\text{i}f\lambda 2f_y\\ \hfill 2\text{i}g\lambda 2g_y\\ \hfill 2\text{i}(|f|^2|g|^2)\end{array})\mathrm{\Phi }.\hfill \end{array}$$
(1)
Here $`u`$, $`f`$ and $`g`$ are complex functions, $`v_1`$ and $`v_2`$ are real functions.

Its integrability conditions $`\mathrm{\Phi }_{xy}=\mathrm{\Phi }_{yx}`$, $`\mathrm{\Phi }_{xt}=\mathrm{\Phi }_{tx}`$ and $`\mathrm{\Phi }_{yt}=\mathrm{\Phi }_{ty}`$ consist of the following three parts.

(1) DSI equation

$$\begin{array}{c}\text{i}u_t=u_{xx}+u_{yy}+2|u|^2u+2(v_1+v_2)u\hfill \\ v_{1,x}v_{1,y}=v_{2,x}+v_{2,y}=(|u|^2)_x.\hfill \end{array}$$
(2)
(2) Standard Lax pair of the DSI equation

$$\begin{array}{c}F_y=\left(\begin{array}{ccc}1& 0& \\ 0& 1& \end{array}\right)F_x+\left(\begin{array}{ccc}0& u& \\ \overline{u}& 0& \end{array}\right)F\hfill \\ F_t=2\text{i}\left(\begin{array}{ccc}1& 0& \\ 0& 1& \end{array}\right)F_{xx}+2\text{i}\left(\begin{array}{ccc}0& u& \\ \overline{u}& 0& \end{array}\right)F_x\hfill \\ +\text{i}\left(\begin{array}{ccc}|u|^2+2v_1& u_x+u_y& \\ \overline{u}_x+\overline{u}_y& |u|^22v_2& \end{array}\right)F\hfill \end{array}$$
(3)
where $`F=(f,g)^T`$.

(3) Nonlinear constraint

$$FF^{}=\frac{1}{2}\left(\begin{array}{ccc}v_1& u_x& \\ \overline{u}_x& v_2& \end{array}\right).$$
(4)
Hence any solution of (2)–(4) gives a solution of the DSI equation.

Notice that if $`\mathrm{\Phi }`$ is a vector solution of (1) for real $`\lambda `$, then $`\mathrm{\Psi }=\text{i}\overline{\mathrm{\Phi }}`$ is a solution of the adjoint equations

$$\begin{array}{c}\mathrm{\Psi }_x=(W^x)^T\mathrm{\Psi }\mathrm{\Psi }_y=(W^y)^T\mathrm{\Psi }\hfill \\ \mathrm{\Psi }_t=(W^t)^T\mathrm{\Psi }\hfill \end{array}$$
(5)
where each entry of $`\overline{\mathrm{\Phi }}`$ is the complex conjugation of the corresponding entry of $`\mathrm{\Phi }`$.

In order to obtain the finite dimensional Hamiltonian systems, we first nonlinearize the $`y`$-equation of (1) in the following way.

Consider the pair $`\mathrm{\Phi }_y=W^y\mathrm{\Phi }`$ and $`\mathrm{\Phi }_t=W^t\mathrm{\Phi }`$. Let $`w=(u,\overline{u},\text{i}f,\text{i}\overline{f},\text{i}g,\text{i}\overline{g})`$ containing all the variables in $`W^y`$ . Then the recursion relations of this AKNS system can be expressed in Lenard form

$$JG_l=KG_{l1}(l=1,2,\mathrm{})$$
(6)
where $`(J,K)`$ is the Lenard pair ($`J`$ is a non-degenerate constant matrix) which was given by  (here $`\alpha _1,\alpha _2,\alpha _3,\beta _1,\beta _2,\beta _3`$ in  are $`\alpha _1=\text{i}`$, $`\alpha _2=\text{i}`$, $`\alpha _3=0`$, $`\beta _1=2\text{i}`$, $`\beta _2=2\text{i}`$, $`\beta _3=0`$) and $`\{G_0,G_1,G_2,\mathrm{}\}`$ is the Lenard sequence. The first element of this Lenard sequence is given by  as

$$G_0=2(\overline{u},u,\text{i}\overline{f},\text{i}f,\text{i}\overline{g},\text{i}g)$$
(7)
which is a kernel of $`K`$.

On the other hand, the variation of spectral parameter can be computed by the general formula

$$\begin{array}{cc}\hfill \frac{\delta \lambda }{\delta w}& =C_0\text{ tr}\left(\mathrm{\Phi }\mathrm{\Psi }^T\frac{W^y}{w}\right)\hfill \\ & =C_0(\text{i}\overline{\varphi }_1\varphi _2,\text{i}\overline{\varphi }_2\varphi _1,\text{i}\overline{\varphi }_1\varphi _3,\text{i}\overline{\varphi }_3\varphi _1,\text{i}\overline{\varphi }_2\varphi _3,\text{i}\overline{\varphi }_3\varphi _2)\hfill \end{array}$$
(8)
where $`C_0`$ is a constant.

Now let $`\lambda _1,\mathrm{},\lambda _N`$ be $`N`$ distinct non-zero real numbers, $`(\varphi _{1\alpha },\varphi _{2\alpha },\varphi _{3\alpha })^T`$ be the corresponding solution of the Lax system (1) for $`\lambda =\lambda _\alpha `$. Let $`\mathrm{\Lambda }=\text{ diag}(\lambda _1,\mathrm{},\lambda _N)`$, $`\mathrm{\Phi }_j=(\varphi _{j1},\mathrm{},\varphi _{jN})^T`$.

By the general formulation of nonlinearization, we impose the nonlinear constraint

$$G_0=2\underset{j=1}{\overset{N}{}}\frac{\delta \lambda _j}{\delta w}$$
(9)
which gives the relations

$$\mathrm{\Phi }_2,\mathrm{\Phi }_1=\text{i}u\mathrm{\Phi }_3,\mathrm{\Phi }_1=f\mathrm{\Phi }_3,\mathrm{\Phi }_2=g$$
(10)
where $`V_1,V_2=V_1^{}V_2`$ for any two vectors $`V_1`$ and $`V_2`$. These are the Bargmann constraints between $`(u,f,g)`$ and $`(\mathrm{\Phi }_1,\mathrm{\Phi }_2,\mathrm{\Phi }_3)`$.

###### Remark 1

If we consider the pair $`\mathrm{\Phi }_y=W^y\mathrm{\Phi }`$ and $`\mathrm{\Phi }_x=W^x\mathrm{\Phi }`$, $`G_0`$ is different. In that case, we can not obtain Bargmann constraints, but only Neumann constraints .

Let

$$\begin{array}{c}L(\lambda )=\left(\begin{array}{ccc}1& & \\ & 1& \\ & & 0\end{array}\right)\hfill \\ +\underset{\alpha =1}{\overset{N}{}}\frac{1}{\lambda \lambda _\alpha }\left(\begin{array}{ccc}\overline{\varphi }_{1\alpha }\varphi _{1\alpha }& \overline{\varphi }_{2\alpha }\varphi _{1\alpha }& \overline{\varphi }_{3\alpha }\varphi _{1\alpha }\\ \overline{\varphi }_{1\alpha }\varphi _{2\alpha }& \overline{\varphi }_{2\alpha }\varphi _{2\alpha }& \overline{\varphi }_{3\alpha }\varphi _{2\alpha }\\ \overline{\varphi }_{1\alpha }\varphi _{3\alpha }& \overline{\varphi }_{2\alpha }\varphi _{3\alpha }& \overline{\varphi }_{3\alpha }\varphi _{3\alpha }\end{array}\right).\hfill \end{array}$$
(11)
###### Lemma 1

The Lax equations

$$L_x=[W^x,L]L_y=[W^y,L]L_t=[W^t,L]$$
(12)
hold if and only if (10) holds.

Proof. Let $`F_\alpha =(\varphi _{1\alpha },\varphi _{2\alpha },\varphi _{3\alpha })^T`$, then

$$L(\lambda )=C+\underset{\alpha =1}{\overset{N}{}}\frac{1}{\lambda \lambda _\alpha }F_\alpha F_\alpha ^{}.$$
(13)
Since $`F_{\alpha ,y}=W^y(\lambda _\alpha )F_\alpha `$ and $`(W^y(\lambda _\alpha ))^{}=W^y(\lambda _\alpha )`$, we have

$$\begin{array}{cc}\hfill L_y(\lambda )& =\underset{\alpha }{}\frac{1}{\lambda \lambda _\alpha }[W^y(\lambda _\alpha ),F_\alpha F_\alpha ^{}]\hfill \\ & =[W^y(\lambda ),L(\lambda )C]\hfill \\ & \underset{\alpha }{}\frac{1}{\lambda \lambda _\alpha }[W^y(\lambda )W^y(\lambda _\alpha ),F_\alpha F_\alpha ^{}]\hfill \\ & =[W^y(\lambda ),L(\lambda )C]\text{i}[C,\underset{\alpha }{}F_\alpha F_\alpha ^{}].\hfill \end{array}$$
(14)
Hence $`L_y(\lambda )=[W^y(\lambda ),L(\lambda )]`$ if and only if

$$\text{i}[C,\underset{\alpha }{}F_\alpha F_\alpha ^{}]=[C,W^y(\lambda )].$$
(15)
Written in the components, it becomes (10). If (10) holds, the other two equations in (12) are obtained similarly. The lemma is proved.

From the above constraints, we have

$$\begin{array}{c}f_x=\text{i}\mathrm{\Phi }_3,\mathrm{\Lambda }\mathrm{\Phi }_1+\text{i}f(\mathrm{\Phi }_3,\mathrm{\Phi }_3\mathrm{\Phi }_1,\mathrm{\Phi }_1)ug\hfill \\ g_x=\text{i}\mathrm{\Phi }_3,\mathrm{\Lambda }\mathrm{\Phi }_2\text{i}g(\mathrm{\Phi }_3,\mathrm{\Phi }_3\mathrm{\Phi }_2,\mathrm{\Phi }_2)\overline{u}f\hfill \\ f_y=\text{i}\mathrm{\Phi }_3,\mathrm{\Lambda }\mathrm{\Phi }_1+\text{i}f(\mathrm{\Phi }_3,\mathrm{\Phi }_3\mathrm{\Phi }_1,\mathrm{\Phi }_1)\hfill \\ g_y=\text{i}\mathrm{\Phi }_3,\mathrm{\Lambda }\mathrm{\Phi }_2+\text{i}g(\mathrm{\Phi }_3,\mathrm{\Phi }_3\mathrm{\Phi }_2,\mathrm{\Phi }_2)\hfill \\ u_x=2f\overline{g}\hfill \\ u_y=2\mathrm{\Phi }_2,\mathrm{\Lambda }\mathrm{\Phi }_1+\mathrm{\Phi }_2,\mathrm{\Phi }_1(\mathrm{\Phi }_1,\mathrm{\Phi }_1\mathrm{\Phi }_2,\mathrm{\Phi }_2)\hfill \\ v_1=2\mathrm{\Phi }_1,\mathrm{\Phi }_3\mathrm{\Phi }_3,\mathrm{\Phi }_1v_2=2\mathrm{\Phi }_2,\mathrm{\Phi }_3\mathrm{\Phi }_3,\mathrm{\Phi }_2\hfill \end{array}$$
(16)
With the constraints (10), the system (1) is changed to a system of ordinary differential equations

$$\begin{array}{c}\mathrm{\Phi }_{1,x}=\text{i}\mathrm{\Lambda }\mathrm{\Phi }_1+\text{i}f\mathrm{\Phi }_3\mathrm{\Phi }_{2,x}=\text{i}\mathrm{\Lambda }\mathrm{\Phi }_2+\text{i}g\mathrm{\Phi }_3\hfill \\ \mathrm{\Phi }_{3,x}=\text{i}\overline{f}\mathrm{\Phi }_1+\text{i}\overline{g}\mathrm{\Phi }_2\hfill \\ \mathrm{\Phi }_{1,y}=\text{i}\mathrm{\Lambda }\mathrm{\Phi }_1+u\mathrm{\Phi }_2+\text{i}f\mathrm{\Phi }_3\hfill \\ \mathrm{\Phi }_{2,y}=\overline{u}\mathrm{\Phi }_1\text{i}\mathrm{\Lambda }\mathrm{\Phi }_2\text{i}g\mathrm{\Phi }_3\hfill \\ \mathrm{\Phi }_{3,y}=\text{i}\overline{f}\mathrm{\Phi }_1\text{i}\overline{g}\mathrm{\Phi }_2\hfill \\ \mathrm{\Phi }_{1,t}=(2\text{i}\mathrm{\Lambda }^2+\text{i}|u|^2+\text{i}v_1)\mathrm{\Phi }_1+(2u\mathrm{\Lambda }+\text{i}u_y)\mathrm{\Phi }_2\hfill \\ +(2\text{i}f\mathrm{\Lambda }2f_y)\mathrm{\Phi }_3\hfill \\ \mathrm{\Phi }_{2,t}=(2\overline{u}\mathrm{\Lambda }+\text{i}\overline{u}_y)\mathrm{\Phi }_1+(2\text{i}\mathrm{\Lambda }^2\text{i}|u|^2\text{i}v_2)\mathrm{\Phi }_2\hfill \\ +(2\text{i}g\mathrm{\Lambda }2g_y)\mathrm{\Phi }_3\hfill \\ \mathrm{\Phi }_{3,t}=(2\text{i}\overline{f}\mathrm{\Lambda }+2\overline{f}_y)\mathrm{\Phi }_1+(2\text{i}\overline{g}\mathrm{\Lambda }+2\overline{g}_y)\mathrm{\Phi }_2\hfill \\ 2\text{i}(|f|^2|g|^2)\mathrm{\Phi }_3\hfill \end{array}$$
(17)
where $`u`$, $`u_y`$, $`v_1`$, $`v_2`$, $`f`$, $`g`$, $`f_y`$, $`g_y`$ are given by (10) and (16).

$`\text{ Re}\varphi _{j\alpha }`$ and $`\text{ Im}\varphi _{j\alpha }`$ $`(j=1,2,3;\alpha =1,\mathrm{},N)`$ form a system of coordinates of the phase space $`\text{R}^{6N}`$. For simplicity, we use the complex coordinates $`\varphi _{j\alpha }`$ and $`\overline{\varphi }_{j\alpha }`$ instead of $`\text{ Re}\varphi _{j\alpha }`$ and $`\text{ Im}\varphi _{j\alpha }`$. $`\text{R}^{6N}`$ has the standard symplectic form

$`\omega `$ $`=`$ $`2{\displaystyle \underset{j=1}{\overset{3}{}}}{\displaystyle \underset{\alpha =1}{\overset{N}{}}}\text{d}\text{ Im}(\varphi _{j\alpha })\text{d}\text{ Re}(\overline{\varphi }_{j\alpha })`$ (18)
$`=`$ $`\text{i}{\displaystyle \underset{j=1}{\overset{3}{}}}{\displaystyle \underset{\alpha =1}{\overset{N}{}}}\text{d}\overline{\varphi }_{j\alpha }\text{d}\varphi _{j\alpha }.`$
The corresponding Poisson bracket of two functions $`\xi `$ and $`\eta `$ is

$$\{\xi ,\eta \}=\text{i}\underset{j=1}{\overset{3}{}}\underset{\alpha =1}{\overset{N}{}}\left(\frac{\xi }{\varphi _{j\alpha }}\frac{\eta }{\overline{\varphi }_{j\alpha }}\frac{\xi }{\overline{\varphi }_{j\alpha }}\frac{\eta }{\varphi _{j\alpha }}\right).$$
(19)
By direct computation, we have

###### Lemma 2

(17) is equivalent to three Hamiltonian equations

$$\begin{array}{cc}\text{i}\varphi _{j\alpha ,x}=\frac{H^x}{\overline{\varphi }_{j\alpha }}\hfill & \text{i}\overline{\varphi }_{j\alpha ,x}=\frac{H^x}{\varphi _{j\alpha }}\hfill \\ \text{i}\varphi _{j\alpha ,y}=\frac{H^y}{\overline{\varphi }_{j\alpha }}\hfill & \text{i}\overline{\varphi }_{j\alpha ,y}=\frac{H^y}{\varphi _{j\alpha }}\hfill \\ \text{i}\varphi _{j\alpha ,t}=\frac{H^t}{\overline{\varphi }_{j\alpha }}\hfill & \text{i}\overline{\varphi }_{j\alpha ,t}=\frac{H^t}{\varphi _{j\alpha }}\hfill \end{array}$$
(20)
where

$$\begin{array}{cc}\hfill H^x=& \mathrm{\Phi }_1,\mathrm{\Lambda }\mathrm{\Phi }_1\mathrm{\Phi }_2,\mathrm{\Lambda }\mathrm{\Phi }_2|\mathrm{\Phi }_1,\mathrm{\Phi }_3|^2\hfill \\ & +|\mathrm{\Phi }_2,\mathrm{\Phi }_3|^2\hfill \\ \hfill H^y=& \mathrm{\Phi }_1,\mathrm{\Lambda }\mathrm{\Phi }_1+\mathrm{\Phi }_2,\mathrm{\Lambda }\mathrm{\Phi }_2|\mathrm{\Phi }_1,\mathrm{\Phi }_3|^2\hfill \\ & |\mathrm{\Phi }_2,\mathrm{\Phi }_3|^2|\mathrm{\Phi }_1,\mathrm{\Phi }_2|^2\hfill \\ \hfill H^t=& 2\mathrm{\Phi }_1,\mathrm{\Lambda }^2\mathrm{\Phi }_12\mathrm{\Phi }_2,\mathrm{\Lambda }^2\mathrm{\Phi }_2\hfill \\ & +4\text{ }\text{Re}(\mathrm{\Phi }_1,\mathrm{\Phi }_3\mathrm{\Phi }_3,\mathrm{\Lambda }\mathrm{\Phi }_1)\hfill \\ & +4\text{ }\text{Re}(\mathrm{\Phi }_2,\mathrm{\Lambda }^{}\mathrm{\Phi }_3\mathrm{\Phi }_3,\mathrm{\Lambda }\mathrm{\Phi }_2)\hfill \\ & +4\text{ }\text{Re}(\mathrm{\Phi }_1,\mathrm{\Phi }_2\mathrm{\Phi }_2,\mathrm{\Lambda }\mathrm{\Phi }_1)\hfill \\ & +2(\mathrm{\Phi }_3,\mathrm{\Phi }_3\mathrm{\Phi }_1,\mathrm{\Phi }_1)|\mathrm{\Phi }_1,\mathrm{\Phi }_3|^2\hfill \\ & 2(\mathrm{\Phi }_3,\mathrm{\Phi }_3\mathrm{\Phi }_2,\mathrm{\Phi }_2)|\mathrm{\Phi }_2,\mathrm{\Phi }_3|^2\hfill \\ & (\mathrm{\Phi }_1,\mathrm{\Phi }_1\mathrm{\Phi }_2,\mathrm{\Phi }_2)|\mathrm{\Phi }_1,\mathrm{\Phi }_2|^2.\hfill \end{array}$$

## 3 Integrability

###### Lemma 3

For any two complex numbers $`\lambda `$, $`\mu `$ and any positive integers $`k`$, $`l`$,

$$\{\text{ }\text{tr}L^k(\lambda ),\text{ }\text{tr}L^l(\mu )\}=0.$$
(21)
Proof. Denote $`e_j=(0,\mathrm{},0,\underset{j}{1},0,\mathrm{},0)^T`$, then

$$\frac{L(\lambda )}{\overline{\varphi }_{j\alpha }}=\frac{1}{\lambda \lambda _\alpha }F_\alpha e_j^T\frac{L(\lambda )}{\varphi _{j\alpha }}=\frac{1}{\lambda \lambda _\alpha }e_jF_\alpha ^{}.$$
(22)
Since

$$\begin{array}{c}\text{ tr}\left(\frac{L^k(\lambda )}{\varphi _{j\alpha }}\right)=\text{ tr}\underset{r=1}{\overset{k}{}}L^{r1}(\lambda )\frac{L(\lambda )}{\varphi _{j\alpha }}L^{kr}(\lambda )\hfill \\ =k\text{ tr}\left(L^{k1}(\lambda )\frac{L(\lambda )}{\varphi _{j\alpha }}\right),\hfill \\ \text{ tr}\left(\frac{L^k(\lambda )}{\overline{\varphi }_{j\alpha }}\right)=k\text{ tr}\left(L^{k1}(\lambda )\frac{L(\lambda )}{\overline{\varphi }_{j\alpha }}\right),\hfill \end{array}$$
(23)
we have

$$\begin{array}{c}\frac{\text{i}}{kl}\{\text{ tr}L^k(\lambda ),\text{ tr}L^l(\mu )\}\hfill \\ =\underset{j,\alpha }{}(\text{ tr}(L^{k1}(\lambda )\frac{L(\lambda )}{\varphi _{j\alpha }})\text{ tr}(L^{l1}(\mu )\frac{L(\mu )}{\overline{\varphi }_{j\alpha }})\hfill \\ \text{ tr}(L^{l1}(\mu )\frac{L(\mu )}{\varphi _{j\alpha }})\text{ tr}(L^{k1}(\lambda )\frac{L(\lambda )}{\overline{\varphi }_{j\alpha }}))\hfill \\ =\underset{a,b,j,\alpha }{}\frac{1}{\lambda \lambda _\alpha }(L^{k1}(\lambda ))_{aj}\overline{\varphi }_{a\alpha }\frac{1}{\mu \lambda _\alpha }(L^{l1}(\mu ))_{jb}\varphi _{b\alpha }\hfill \\ \underset{a,b,j,\alpha }{}\frac{1}{\mu \lambda _\alpha }(L^{l1}(\mu ))_{aj}\overline{\varphi }_{a\alpha }\frac{1}{\lambda \lambda _\alpha }(L^{k1}(\lambda ))_{jb}\varphi _{b\alpha }\hfill \\ =\underset{a,b,\alpha }{}\frac{1}{\lambda \lambda _\alpha }\frac{1}{\mu \lambda _\alpha }\overline{\varphi }_{a\alpha }\varphi _{b\alpha }[L^{k1}(\lambda ),L^{l1}(\mu )]_{ab}\hfill \\ =\underset{a,b,\alpha }{}\frac{1}{\mu \lambda }(\frac{1}{\lambda \lambda _\alpha }\frac{1}{\mu \lambda _\alpha })\overline{\varphi }_{a\alpha }\varphi _{b\alpha }\hfill \\ [L^{k1}(\lambda ),L^{l1}(\mu )]_{ab}\hfill \\ =\frac{1}{\mu \lambda }\text{ tr}\left((L(\lambda )L(\mu ))[L^{k1}(\lambda ),L^{l1}(\mu )]\right)=0.\hfill \end{array}$$
(24)
The lemma is proved.

Using this lemma, we can construct involutive integrals of motion from $`\text{ tr}L^k`$ $`(k=1,2,\mathrm{})`$.

For any complex number $`\xi `$, let

$$det(\xi L(\lambda ))=\xi ^3p_1(\lambda )\xi ^2+p_2(\lambda )\xi p_3(\lambda )$$
(25)
then $`p_k(\lambda )`$ is the sum of all the determinants of the principal submatrices of $`L(\lambda )`$ of order $`k`$. Suppose the eigenvalues of $`L(\lambda )`$ are $`\nu _1(\lambda ),\nu _2(\lambda ),\nu _3(\lambda )`$, then

$$\begin{array}{c}\text{ tr}L^k(\lambda )=\underset{j=1}{\overset{3}{}}\nu _j^k(\lambda )\hfill \\ p_k(\lambda )=\underset{1j_1<\mathrm{}<j_k3}{}\nu _{j_1}(\lambda )\mathrm{}\nu _{j_k}(\lambda ).\hfill \end{array}$$
(26)
Since $`p_k(\lambda )`$ can be expressed as a polynomial of $`\text{ tr}L^l(\lambda )`$ $`(l=1,2,\mathrm{})`$, $`\{p_j(\lambda ),p_k(\mu )\}=0`$ for any two complex numbers $`\lambda `$ and $`\mu `$ and two integers $`j,k0`$.

Expand $`p_k(\lambda )`$ as a Laurant series of $`\lambda `$

$$p_k(\lambda )=\underset{m=1}{\overset{\mathrm{}}{}}E_m^{(k)}\lambda ^{m1}$$
(27)
which is convergent when $`|\lambda |>\underset{1\alpha N}{\mathrm{max}}|\lambda _\alpha |`$, then

$$\begin{array}{c}E_m^{(1)}=\mathrm{\Phi }_1,\mathrm{\Lambda }^m\mathrm{\Phi }_1+\mathrm{\Phi }_2,\mathrm{\Lambda }^m\mathrm{\Phi }_2+\mathrm{\Phi }_3,\mathrm{\Lambda }^m\mathrm{\Phi }_3\hfill \\ E_m^{(2)}=\mathrm{\Phi }_1,\mathrm{\Lambda }^m\mathrm{\Phi }_1+\mathrm{\Phi }_2,\mathrm{\Lambda }^m\mathrm{\Phi }_2\hfill \\ +\underset{1i<j3}{}\underset{l=1}{\overset{m}{}}\left|\begin{array}{cc}\mathrm{\Phi }_i,\mathrm{\Lambda }^{l1}\mathrm{\Phi }_i& \mathrm{\Phi }_j,\mathrm{\Lambda }^{ml}\mathrm{\Phi }_i\\ \mathrm{\Phi }_i,\mathrm{\Lambda }^{l1}\mathrm{\Phi }_j& \mathrm{\Phi }_j,\mathrm{\Lambda }^{ml}\mathrm{\Phi }_j\end{array}\right|\hfill \\ E_m^{(3)}=\mathrm{\Phi }_3,\mathrm{\Lambda }^m\mathrm{\Phi }_3\hfill \\ \underset{l=1}{\overset{m}{}}\left|\begin{array}{cc}\mathrm{\Phi }_1,\mathrm{\Lambda }^{l1}\mathrm{\Phi }_1& \mathrm{\Phi }_3,\mathrm{\Lambda }^{ml}\mathrm{\Phi }_1\\ \mathrm{\Phi }_1,\mathrm{\Lambda }^{l1}\mathrm{\Phi }_3& \mathrm{\Phi }_3,\mathrm{\Lambda }^{ml}\mathrm{\Phi }_3\end{array}\right|\hfill \\ +\underset{l=1}{\overset{m}{}}\left|\begin{array}{cc}\mathrm{\Phi }_2,\mathrm{\Lambda }^{l1}\mathrm{\Phi }_2& \mathrm{\Phi }_3,\mathrm{\Lambda }^{ml}\mathrm{\Phi }_2\\ \mathrm{\Phi }_2,\mathrm{\Lambda }^{l1}\mathrm{\Phi }_3& \mathrm{\Phi }_3,\mathrm{\Lambda }^{ml}\mathrm{\Phi }_3\end{array}\right|\hfill \\ +\underset{\genfrac{}{}{0pt}{}{i+j+k=m2}{i,j,k0}}{}|\begin{array}{c}\mathrm{\Phi }_1,\mathrm{\Lambda }^i\mathrm{\Phi }_1\\ \mathrm{\Phi }_1,\mathrm{\Lambda }^i\mathrm{\Phi }_2\\ \mathrm{\Phi }_1,\mathrm{\Lambda }^i\mathrm{\Phi }_3\end{array}\begin{array}{c}\mathrm{\Phi }_2,\mathrm{\Lambda }^j\mathrm{\Phi }_1\\ \mathrm{\Phi }_2,\mathrm{\Lambda }^j\mathrm{\Phi }_2\\ \mathrm{\Phi }_2,\mathrm{\Lambda }^j\mathrm{\Phi }_3\end{array}\begin{array}{c}\mathrm{\Phi }_3,\mathrm{\Lambda }^k\mathrm{\Phi }_1\\ \mathrm{\Phi }_3,\mathrm{\Lambda }^k\mathrm{\Phi }_2\\ \mathrm{\Phi }_3,\mathrm{\Lambda }^k\mathrm{\Phi }_3\end{array}|.\hfill \end{array}$$
(28)
The sums are zero if the lower bound is greater than the upper bound. These $`E_m^{(k)}`$’s are in involution.

Let

$$\begin{array}{c}\mathrm{\Omega }_1=\mathrm{\Phi }_1,\mathrm{\Phi }_1=(E_0^{(1)}E_0^{(2)}+E_0^{(3)})/2\hfill \\ \mathrm{\Omega }_2=\mathrm{\Phi }_2,\mathrm{\Phi }_2=(E_0^{(1)}+E_0^{(2)}+E_0^{(3)})/2\hfill \\ \mathrm{\Omega }_3=\mathrm{\Phi }_3,\mathrm{\Phi }_3=E_0^{(3)}\hfill \end{array}$$
(29)
then

$$\begin{array}{c}H^x=E_1^{(1)}E_1^{(3)}(\mathrm{\Omega }_1\mathrm{\Omega }_2)\mathrm{\Omega }_3\hfill \\ H^y=E_1^{(2)}\mathrm{\Omega }_1\mathrm{\Omega }_2\mathrm{\Omega }_1\mathrm{\Omega }_3\mathrm{\Omega }_2\mathrm{\Omega }_3\hfill \\ H^t=2E_2^{(2)}+2(\mathrm{\Omega }_1+\mathrm{\Omega }_2)E_1^{(1)}\hfill \\ +(\mathrm{\Omega }_1+\mathrm{\Omega }_22\mathrm{\Omega }_3)H^x+(\mathrm{\Omega }_1\mathrm{\Omega }_2)H^y\hfill \end{array}$$
(30)
which have all been expressed as polynomials of $`E_m^{(k)}`$’s.

Hence the following lemma holds.

###### Lemma 4

$$\{H^x,H^y\}=\{H^x,H^t\}=\{H^y,H^t\}=0,$$
(31)
$$\{H^x,E_m^{(k)}\}=\{H^y,E_m^{(k)}\}=\{H^t,E_m^{(k)}\}=0$$
(32)
and

$$\{E_m^{(k)},E_p^{(l)}\}=0$$
(33)
for all $`k,l=1,2,3`$; $`m,p=0,1,2,\mathrm{}`$.

Now we consider the independence of $`E_m^{(k)}`$’s.

###### Lemma 5

$`E_m^{(k)}`$ $`(1k3;\mathrm{\hspace{0.17em}0}mN1)`$ are functionally independent in a dense open subset of $`\text{R}^{6N}`$.

Proof. Let $`P_0`$ in $`\text{R}^{6N}`$ be given by $`\varphi _{j\alpha }=ϵ`$ $`(j=1,2,3;\alpha =1,2,\mathrm{},N)`$ where $`ϵ`$ is a small real constant. Then, at $`P_0`$,

$$\frac{E_m^{(k)}}{\overline{\varphi }_{j\alpha }}=\epsilon \lambda _\alpha ^mb_{kj}+O(\epsilon ^3)$$
(34)
where

$$(b_{kj})=\left(\begin{array}{ccc}1& 1& 1\\ 1& 1& 0\\ 0& 0& 1\end{array}\right).$$
(35)
The Jacobian determinant

$$\begin{array}{c}J\frac{(E_0^{(1)},\mathrm{},E_{N1}^{(1)},E_0^{(2)},\mathrm{},E_{N1}^{(2)},E_0^{(3)},\mathrm{},E_{N1}^{(3)})}{(\overline{\varphi }_{11},\mathrm{},\overline{\varphi }_{1N},\overline{\varphi }_{21},\mathrm{},\overline{\varphi }_{2N},\overline{\varphi }_{31},\mathrm{},\overline{\varphi }_{3N})}\hfill \\ =2ϵ^{3N}\left(\underset{1\alpha <\beta N}{}(\lambda _\beta \lambda _\alpha )\right)^3+O(ϵ^{3N+2}).\hfill \end{array}$$
(36)
$`J0`$ at $`P_0`$ if $`ϵ0`$ is small enough. Since $`J`$ is a real analytic function of

$$(\varphi _{11},\mathrm{},\varphi _{3N},\overline{\varphi }_{11},\mathrm{},\overline{\varphi }_{3N}),$$

$`J0`$ in a dense open subset of $`\text{R}^{6N}`$. Therefore, the Jacobian determinant of

$$(E_0^{(1)},\mathrm{},E_{N1}^{(1)},E_0^{(2)},\mathrm{},E_{N1}^{(2)},E_0^{(3)},\mathrm{},E_{N1}^{(3)})$$

to $`6N`$ real coordinates

$$(\text{ Re}(\varphi _{11}),\mathrm{},\text{ Re}(\varphi _{3N}),\text{ Im}(\varphi _{11}),\mathrm{},\text{ Im}(\varphi _{3N}))$$

is of full rank $`3N`$. The lemma is proved.

In summary, from Lemma 2–5, we have

###### Theorem 1

The Hamiltonian systems given by the Hamiltonians (2) are involutive and completely integrable in Liouville sense. The integrals of motion are $`E_m^{(k)}`$ $`(k=1,2,3;m=0,\mathrm{},N1)`$ given by (28), which are in involution and functionally independent in a dense open subset of the phase space $`\text{R}^{6N}`$. Moreover, each solution of these Hamiltonian systems gives a solution of the DSI equation.

## Acknowledgements

This work was supported by the Chinese Major State Basic Research Project “Nonlinear Science”, the City University of Hong Kong (Grant No. 7001041), the Research Grants Council of Hong Kong (Grant No. 9040395,9040466), the Doctoral Program Foundation and the Key Project for Young Teachers of the Ministry of Education of China. The first author (Z. X. Zhou) is grateful to the Department of Mathematics of the City University of Hong Kong for the hospitality.


format: text

text:
# Quantitative Analysis of the Void Size Distribution in the LCRS and in CDM Models

## 1 Voids in the LCRS

The appearence of large voids in the galaxy distribution of extended galaxy surveys as the Las Campanas Redshift Survey (LCRS, Shectman et al. ) was one of the most spectacular findings characterizing the large-scale matter distribution of the universe. We expect large voids due to the matter outflow from positive primordial potential perturbations . A second mechanism for the lack of galaxies in voids is the lower probability of galaxy formation in underdense regions known as biasing. Recently, Friedmann and Piran  presented a simple model for void formation based on the probability distribution of spherical underdense regions in the universe, and the less abundance of density peaks therein. There it was aimed to explain that in recent redshift surveys voids have radii between 13 and 30 $`h^1`$Mpc, occupy 50% of the volume, and that the probability of larger voids falls off exponentially. Since the theoretical conclusions are based on a simple analytical model, and since galaxy surveys before LCRS were not extended enough to get reliable void statistics, we have undertaken a new quantitative study of voids in the LCRS and in numerical simulations . Here this analysis is evaluated and compared with theoretical expectations.

In our analysis  we have used automatic void finders in 2D. Voids are identified as connected empty cells on a grid, where first square base voids are identified, and then boundary layers are added if they are empty and larger than 2/3 of the previous extension. The void finder starts with the search for the largest voids. The algorithm was first described and tested in 2D using slices through the CfA galaxy redshift survey by Kauffmann & Fairall .

From the LCRS we have selected 14 volume limited subsets in different radial ranges of the six narrow slices, projected them into the central planes and corrected for variations in the sampling fraction. In Fig. 1 the void distribution of such sets is shown which illustrates that we got a large number of voids in each data set, and that the voids cover smoothly the surface of the slice, with no obvious inhomogeneities or biasing by boundary effects or radial gradients. To emphasize the importance of large voids, we constructed the size distribution as cumulative coverage of the survey area with voids. Results are shown in Fig. 2. In the left pannel we overplot the void distributions of all 14 samples in normalizing them with the median void sizes from each sample. The thick lines show the results from the best sampled part of the survey. Altogether we get a nice convergence to a mean behavior shown by a smooth parametric fit. Also the derivative is shown in comparison with the distribution of voids in random point distributions which shows the main result of the gravitational clustering of galaxies in overdense regions: The relative abundence both of large voids and of small voids in the data is larger than in the random point distributions (the dashed lines).

## 2 Comparison with simulations

We have compared the LCRS voids with cold dark matter simulations in spatially flat $`\mathrm{\Lambda }`$-term cosmology ($`\mathrm{\Lambda }`$CDM, $`\mathrm{\Omega }_0=0.3`$) and in open models (OCDM, $`\mathrm{\Omega }_0=0.5`$). PM simulations were performed in 500 $`h^1`$Mpc boxes using $`300^3`$ particles and a threshold bias prescription to identify ‘galaxies’ with single particles. In addition we analyzed P3M simulations in 280 $`h^1`$Mpc boxes and with $`256^3`$ particles, a friend-of-friend halo finder, and a procedure to split off unvirialized halos that result from particle clumps merged by numerical effects. The mock samples in the different simulations fit well the 2-point correlation function of the LCRS as already shown in Tucker et al. . In the right pannel of Fig. 3, we compare the real space correlation function of the mock samples with the reconstructed 3D correlation function of APM galaxies . The mock samples provide good fits to the observations especially on large scales.

The cumulative void size distribution of ten OCDM mock samples is shown in the right pannel of Fig. 2 together with the variance range. Obviously this model has difficulty in matching the void distribution for a large range from 15 to 35 $`h^1`$Mpc. The $`\mathrm{\Lambda }`$CDM fits better, but it also seems to underestimate the size of the largest voids.

## 3 Scaling Relation

The 14 analyzed volume-limited data sets show variations in the galaxy densities which are due to different absolute magnitude cuts and different sampling fractions. We already showed  that differences in the magnitude ranges have negligible influence on the size distribution of voids. So we ascribed the differences in the parameters of the void distributions to variations in the galaxy sampling densities. We measured it by the mean galaxy separation $`\lambda `$. In analyzing both data and mock samples, we found a characteristic scaling relation for the mean void sizes that cover 25%, 50% and 75% of the survey area. This scaling is a simple linear relation $`D=D_0+\nu \lambda `$, with size $`D_0=(6,12,17)h^1\text{Mpc}`$ and slope $`\nu =(0.9,1.1,1.5)`$ for 25%, 50% and 75% of void coverage, cp. Table 1.

In Fig. 3 we show that the scaling relations in simulated OCDM mock samples are steeper, $`\nu _{model}>\nu _{data}`$. It is even steeper for the $`\mathrm{\Lambda }`$CDM model as seen in Table 1. According to our analysis, the median size of voids is larger than reported from earlier surveys . The largest voids have a size of $`D_{max}50h^1\text{Mpc}`$ not met in simulations. Tests with mock samples indicate that the 2D-void size distributions in LCRS slices have a fixed ratio to the 3D-voids. Further, we have shown that large voids increase in size by 10 % in redshift space, but small voids remain unchanged, cp. also . Our mock samples were always analyzed in redshift space.

## 4 Discussion

The void statistics provide important characteristics of the galaxy distribution in low to medium dense regions where the correlation function is less sensitive. The self-similarity of the size distribution of samples with different mean galaxy separations culled by random dilutions and the scaling relation of the percentile sizes represent a genuine feature of hierarchical clustering. The discrepancy of the simulated void distribution to observations requires a stronger biasing of simulated galaxies and a higher concentration of them in superclusters, as has also been seen by analyzing overdensity regions with the minimal spanning tree technique illustrated in Fig. 1, cp. Doroshkevich et al. .


text:
# Proportion of frozen local polarization in relaxor ferroelectrics

## Abstract

A Landau-type phenomenological cluster theory was presented to model the freezing process of local polarization in relaxors. Based on the theory, the proportion of frozen polarization in Pb(Mg1/3Nb2/3)O3-PbTiO3 (PMNT) was calculated from the experiment of dielectric nonlinearity. The local polarization was shown to freeze continuously in a cooling process. The amount of frozen polarization increases with increasing the measuring frequency.

Relaxor ferroelectrics (relaxors) have been showing great importance both for the fundamental solid state science and the applications to advanced technology. Relaxors experience no macroscopic phase transition at zero electric field until very low temperatures. However, local polarization exists at much higher temperatures, which is widely believed to be responsible for the special characteristics of relaxors. The concept of “polar micro region” in relaxors was firstly presented by Smolenskii in the chemical inhomogeneity theory. The presence of polar microregions was later confirmed by the experiments on the non-linear variation of the optical refractive index, the thermal strain, and the thermal expansion coefficient. The deviation from the paraelectric Curie-Weiss behavior of the dielectric permittivity also suggested the presence of polar regions. Cross presented in a superparaelectric model that the polar microregions (represented by independent dipoles) are able to switch between the equivalent orientation states and gave an appropriate picture of the nature of the micro polarization at high temperatures. When the temperature decreases, the coupling between the polar microregions controls the kinetics of the fluctuations and the system is frozen into a polar-glassy state. Gui et al. used the Monte Carlo method to simulate the freezing process in relaxors and showed that some dipoles are slowed or frozen due to the interactions. A special phase transition of ergodic space shrinking in succession was also proposed according to the freezing process. Kleemann et al., however, proposed that the polarization is frozen due to the quenched random electric fields.

Some works have been done to quantitatively analyze the freezing process of local polarization. Qian and Bursill developed a phenomenological theory to describe the interaction of polar domains and to simulate the dielectric relaxation and phase transition of relaxors. Nambu and Sugimoto proposed a Landau type mean field theory by considering a gradual condensation of local polarization and confirm a picture of diffuse phase transition in relaxors. In this paper, we propose a Landau-type phenomenological cluster theory of relaxors and calculate the proportion of frozen polarization from the experiments of dielectric nonlinearity.

A phenomenological free energy is defined as

$$F=2^1N\alpha (T)P^2+4^1N\alpha _{11}P^4+2^1\gamma P^2\underset{i}{}n_iP_i^2+2^1\underset{i}{}n_i\alpha _i(T)P_i^2+4^1\alpha _{11}^{}\underset{i}{}n_iP_i^4,$$
(1)
where $`P`$ is the uniform polarization and $`P_i`$ is the $`i`$-th frozen local polarization. $`N`$ is the total number of lattice sites, and $`n_i`$ is the lattice-site number of the $`i`$-th local polarization. $`\alpha `$ and $`\alpha _i`$ are written as

$`\alpha (T)=(TT_0)/\epsilon _0C,`$ (2)
$`\alpha _i(T)=(TT_i)/\epsilon _0C_i,`$ (3)
and $`\alpha _{11}`$ and $`\alpha _{11}^{}`$ are constants as the conventional Landau theory. The couplings between the global polarization and the local ones are written as $`\gamma P^2P_i^2`$ from the symmetry consideration.

The equilibrium values of the frozen polarization $`\{P_i\}`$ are determined by minimizing the free energy $`F`$ in Eq. (1). Minimizing $`F`$ with respect to $`\{P_i\}`$ gives

$$F/P_i=0,$$
(4)
i.e.,

$$\gamma P^2n_iP_i+n_i\alpha _i(T)P_i+\alpha _{11}^{}n_iP_i^3=0.$$
(5)
Thus $`P_i`$ is solved from Eq. (5) as

$$P_i^2=\alpha _i(T)/\alpha _{11}^{}\gamma P^2/\alpha _{11}^{}.$$
(6)
Define $`P_{i0}`$ as

$$P_{i0}^2=\alpha _i(T)/\alpha _{11}^{}=(T_iT)/\epsilon _0C_i\alpha _{11}^{},$$
(7)
and then $`P_i`$ can be rewritten as

$$P_i^2=P_{i0}^2\gamma P^2/\alpha _{11}^{}.$$
(8)
When relaxors stays in a unpolar state ($`P`$=0), the frozen polarization $`P_i`$ is equal to $`P_{i0}`$. Substituting Eqs (7) and (8) into Eq. (1) yields the free energy

$$\overline{F}=2^1NP^2\left[\alpha (T)+\gamma \frac{1}{N}\underset{i}{}n_iP_{i0}^2\right]+4^1NP^4\left[\alpha _{11}\frac{\gamma ^2}{\alpha _{11}^{}}\frac{1}{N}\underset{i}{}n_i\right]+\mathrm{const}.$$
(9)
Introduce

$$q(T)=\frac{1}{N}\underset{i}{}n_iP_{i0}^2$$
(10)
and

$$n(T)=\frac{1}{N}\underset{i}{}n_i,$$
(11)
and then Eq. (9) can be expressed as

$$\overline{F}=2^1Na(T)P^2+4^1Nb(T)P^4+\mathrm{const},$$
(12)
where

$`a(T)=\alpha (T)+\gamma q(T),`$ (13)
$`b(T)=\alpha _{11}{\displaystyle \frac{\gamma ^2}{\alpha _{11}^{}}}n(T).`$ (14)
For relaxors, the sum of local polarizations is equal to zero, while the sum of the square of local polarizations is not zero that can be measured in the nonlinear behaviors of optical refractive index and the thermal strain, etc. So $`q(T)`$ is the order parameter to represent the appearence of local polarizations. $`n(T)`$ is the proportion of the lattice site of frozen polarization. It should be noted that Eq. (7) is valid only for $`T<T_i`$ (for $`T>T_i`$, $`P_{i0}`$=0, i.e., the $`i`$-th local polarization is not frozen), so the summarization in Eqs. (9-11) can be conducted only on the nonzero local polarization ($`P_{i0}^2>0`$), and $`n(T)`$ and $`q(T)`$ vary with temperatures.

The dielectric susceptibility is easily calculated from the thermodynamic relation between the dielectric field,

$$E=\frac{1}{N}\frac{\overline{F}}{P}=a(T)P+b(T)P^3,$$
(15)
and the uniform polarization $`P`$, i.e.,

$$\frac{1}{\epsilon (T)}=\frac{\epsilon _0E}{P}|_{E=0}=\epsilon _0a(T)=\frac{TT_0}{C}+\epsilon _0\gamma q(T).$$
(16)
It can be seen that the dielectric susceptibility deviates from the Curie-Weiss behavior due to the existence of the local order parameter $`q(T)`$. Schmitt and Kirsch gave the similar relation of Eq. (16) by simple phenology in (Pb,La)(Zr,Ti)O3 (PLZT) system. Nambu and Sugimoto derived the same relation on the basis of a more general mean field theory and explained the difference of the phase transitions in PMN and Pb(Sc1/2Ta1/2)O3 (PST). The purpose of this works was to investigate $`n(T)`$, but not $`q(T)`$. So we consider the dielectric nonlinearity, i.e., the susceptibility corresponding to a nonzero electric field $`E`$,

$$\epsilon _E(T)=\frac{P}{\epsilon _0E}=\frac{1}{\epsilon _0}\frac{1}{a(T)+b(T)P^2}\frac{1}{\epsilon _0a(T)}\frac{b(T)}{\epsilon _0a^2(T)}P^2=\epsilon (T)\epsilon _0^3E^2\epsilon ^2(T)\epsilon _E^2(T)b(T).$$
(17)
According to Eqs (17) and (14), $`n(T)`$ can be obtained as

$$n(T)=\frac{\alpha _{11}^{}}{\gamma ^2}\left[\alpha _{11}+\frac{1}{\epsilon _0^3E^2}\frac{\epsilon _E(T)\epsilon (T)}{\epsilon ^2(T)\epsilon _E^2(T)}\right].$$
(18)
The value of $`\alpha _{11}`$ can be obtained from the experiment data at high temperatures by setting $`n(T)=0`$ in Eq. (18) since there is no frozen polarization at high temperatures. So, by using Eq. (18), the proportion of frozen polarization, $`n(T)`$, can be determined except a constant coefficient $`\alpha _{11}^{}/\gamma ^2`$ from the experiment of dielectric nonlinearity.

We take a typical and well-known relaxors, Pb(Mg1/3Nb2/3)O3-PbTiO3 (PMNT), as an example to analyze the proportion of frozen polarization. The samples were prepared by two-stage calcination method. The magnesium niobate is synthesized by calcining MgCO$`{}_{3}{}^{}`$Mg(OH)$`{}_{2}{}^{}`$6H2O and Nb2O5 at 1000$`{}_{}{}^{}C`$ for 6 h. The columbite phase was then mixed with lead oxide and titanium oxide to form the composition 0.96PMN-0.04PT. Excess PbO (0.3 wt%) was added to compensate for PbO loss during heat treatment. The mixture was ball-milled and calcined again at 800$`{}_{}{}^{}C`$ for 2 h. A 10 mm uniaxial steel die was employed to produce green pellets from the calcined powder, using a pressure of 100 MPa. The pressed pellets were subsequently sintered in a covered alumina crucible at 1200$`{}_{}{}^{}C`$ for 1 h. After sintering, the pellets were polished to a thickness of 0.6 mm and silver paste was applied and fired at 600$`{}_{}{}^{}C`$ to achieve a conductive and adherent coating. The dielectric susceptibility was measured at various frequencies with an LCR precision meter (Model HP 4284A, Hewlett Packard, Palo Alto, CA) remotely controlled through a desktop computer. A temperature chamber (Model 2300, Delta Design, San Diego, CA) was interfaced to the computer to allow measurement of dielectric properties at various temperatures. The amplitude of the ac measurement field was 0.05 and 0.25 kV/cm.

An illustration of the temperature dependence of the dielectric response at a low frequency 100 Hz is shown in Fig. 1(a), and the corresponding calculated proportion of frozen polarization is given in Fig. 1(b). The curves in Fig. 1(a) demonstrate typical dielectric nonlinear behavior with the magnitude of the susceptibility increasing with increasing the field amplitude and the maximum shifting to lower temperatures, which is consistent with the previous observations. The curve of the proportion of frozen polarization, $`n(T)`$, is determined with an arbitrary factor. It can be seen that $`n(T)`$ is equal to zero at high temperatures, which means that the local polarization is all dynamic. When the temperature decreases to approach the temperature of the susceptibility maximum ($`T_m`$), $`n(T)`$ starts to increase, i.e., some local polarization is frozen. After crossing $`T_m`$, $`n(T)`$ increases rapidly with decreasing temperature, and an extrapolation of the slope at the inflection point to zero yielded a critical temperature $`T_p=8.3`$$`{}_{}{}^{}C`$. The monotonous increasing of $`n(T)`$ with decreasing temperature in Fig. 1(b) shows that the calculation is physically realistic.

In the works of other researchers, the calculated quantity is the local order parameter $`q(T)`$ that is determined from the temperature dependence of zero-field susceptibility \[see Eq. (16)\]. The increase of $`q(T)`$ with decreasing temperature is ensured by the decrease of susceptibility at low temperatures. However, the quantity investigated here, $`n(T)`$, is calculated from the nonlinear effect \[see Eq. (18)\]. It should be noted that the dielectric nonlinearity is the strongest when the temperature is near $`T_m`$ while the proportion of frozen polarization is not large at that temperature range. It seems to be conflict with the thought that the frozen polarization is the origin of the dielectric nonlinearity. The key to answer this problem is that the nonlinearity is approximately proportional to the forth power of the susceptibility \[see Eq. (17)\]. The susceptibility reaches the maximum at $`T_m`$, so the nonlinear is strong near $`T_m`$. This can also explain the weakening of the nonlinearity when the proportion of frozen polarization increases in further decreasing the temperature.

Another feature of Fig. 1 is that the frozen polarization, i.e., the local polarization which experiences a phase transition and loses the ergodicity, does not appear suddenly, but increases continuously, which is accompanied with the decreasing of the susceptibility. It is a kind of phase transition where the ergodic space shrinks in succession.

A few sentences can be presented here to explore the origin of the effect of dielectric amplitude on the susceptibility. From Eq. (8) one knows that the magnitude of frozen polarization decreases with increasing external field. In other words, the driving forces on the polarization is enhanced when the ac field amplitude increases, so the frozen polarization is forced to flip faster, and some frozen polarization is unfrozen and gives contribution to the polarization process. Thus the susceptibility increases with increasing the field amplitude.

Figure 2 shows the curves of frozen polarization proportion at different measuring frequencies. When the frequency increases, the curve of the proportion shifts slightly towards higher temperatures. It implies that the time scale of polarization flipping is shortened when the frequency increases, so more polarization cannot reach the equilibrium states in the observation time, i.e., more polarization is frozen. The slight increase of frozen polarization would result in the decreasing of susceptibility at higher frequency, which is known as the frequency dispersion in relaxors. The curve of $`T_p`$ \[$`T_p`$ is defined in Fig. 1(b)\] is depicted in Fig. 3 together with the curve of $`T_m`$ for comparison. It clearly demonstrates the increasment of $`T_p`$ with increasing frequency. And we can see that $`T_m`$ increases more rapidly than $`T_p`$ does.

There are several points that should be discussed here. Firstly, the couplings between local polarization are not included in Eq. (1). However, if the couplings between polarization are considered, through the similar solving procedure as that in Ref. 15, the results presented in this paper is still valid by redefining the parameters. Secondly, we only consider the forth-power terms in the free energy and ignore the influence of the higher terms in this works. And the coefficients $`\alpha _{11}`$, $`\alpha _{11}^{}`$ and $`\gamma `$ are considered as temperature-independent. All these defects restrict the applications of the theory. At last, in this paper the local properties are included within the Landau expansion by introducing the local changing quantities $`n_i`$. Consequently, one expects the appearance of a distribution for $`n_i`$ which is missed in the model. As the result all quantities should be averaged on such a distribution.

In summary, a Landau-type phenomenological cluster theory of relaxors is proposed in this paper to evaluate the influence of frozen polarization. The proportion of frozen polarization in PMNT is calculated from the experimental data of dielectric nonlinearity. It is shown that the frozen polarization starts to appear when the temperature decreases near to $`T_m`$, and increases rapidly after crossing $`T_m`$. When the ac field frequency increases, the amount of frozen polarization arises too.

This work was supported by the Chinese National Science Foundation (Grant NO. 59995520) and State Key Program of Basic Research Development (Grant No. G2000067108).