text:
# 1 Symmetrization and the result on a local-in-time smoothness

## 1 Symmetrization and the result on a local-in-time smoothness

The result on a local in time existence of smooth solution $`𝐔(t,𝐱)`$ for symmetric hyperbolic systems, i.e. systems of the form

$$A^0(t,𝐱,𝐔)\frac{𝐔}{t}+\underset{k=1}{\overset{n}{}}A^k(t,𝐱,𝐔)\frac{𝐔}{𝐱^k}=g(t,𝐱,𝐔),$$
$`(1.1)`$
where the matrices $`A^j(t,𝐱,𝐔)`$ are symmetric, and additionally the matrix $`A^0(t,𝐱,𝐔)`$ is positive definite, is well known (,, the particular cases in , , )).

Namely, if the matrices $`A^j(t,𝐱,𝐩),g(t,𝐱,𝐩)`$ depending smoothly on their arguments, have continuous and bounded derivatives with respect to the variables $`(𝐱,𝐩)`$ up to order $`m+1`$ under bounded $`𝐩`$, the initial data $`𝐔_0=𝐔(0,𝐱)`$ and the function $`g(t,𝐱,\mathrm{𝟎})`$ belong to the Sobolev class $`H^m(^n),m>1+n/2,`$ for any fixed $`t0`$, then locally in time the corresponding Cauchy problem has a unique solution from the class $`_{j=0}^1C^j([0,T);H^{mj}(^n)),T>0.`$ Moreover,

$$\underset{tT0}{lim}sup(𝐔_L^{\mathrm{}}+_𝐱𝐔_L^{\mathrm{}})=+\mathrm{},$$
$`(1.2)`$
$`T=T(𝐔_0)`$ and $`\underset{𝐔_0_{H^m}0}{lim}T(𝐔_0)=+\mathrm{}.`$

(Note that less rigid requirements may be imposed on the initial data for the existence of smooth solution locally in time, f.e. local (in space), but uniform belonging to $`H^m(^n).`$ ,).

In many problems having the origin in physics the coefficients of system (1.1) depend only on solution, this simplifies significantly the formulation of the result. In the case for the local in time existence of the Cauchy problem in the Sobolev class described above, it is sufficient, beside of the coefficients smoothness, to require the implementation of the condition $`g(\mathrm{𝟎})=\mathrm{𝟎}.`$

As we deals with the solutions to system (1–3) such that $`inf\rho =infp=0,`$ we have to use the symmetrizartion firstly proposed in .

Involving the variable $`\mathrm{\Pi }=\kappa \left(p/2\right)^{\frac{\gamma 1}{2\gamma }},\kappa ={\displaystyle \frac{2\sqrt{\gamma }}{\gamma 1}},`$ we get the symmetric form as follows:

$$\mathrm{exp}(\frac{S}{\gamma })(_t+𝐕,)\mathrm{\Pi }+\frac{\gamma 1}{2}\mathrm{exp}(\frac{S}{\gamma })\mathrm{\Pi }div𝐕=0,$$
$`(1.3)`$
$$(_t+(𝐕,))𝐕+\frac{\gamma 1}{2}\mathrm{exp}(\frac{S}{\gamma })\mathrm{\Pi }\mathrm{\Pi }=𝐟_1(t,𝐱,\mathrm{\Pi },𝐕,S)=𝐟(t,𝐱,e^{\frac{S}{\gamma }}\left(\frac{\mathrm{\Pi }}{\kappa }\right)^{\frac{2}{\gamma 1}},𝐕,S),$$
$`(1.4)`$
$$(_t+(𝐕,))S=0.$$
$`(1.5)`$
Denote $`𝐔=(\mathrm{\Pi },𝐕,S).`$ As an immediate corollary of the general result on the symmetric hyperbolic systems we obtain the following

###### Theorem 1. 1

Let the initial data $`𝐔_0=(\mathrm{\Pi }_0,𝐕_0,S_0),`$ belong to the class $`H^m(^n),m>1+n/2.`$ Suppose the function $`𝐟_1(t,𝐱,𝐔))`$ have continuous and bounded derivatives with respect to space variables and the solution components up to order $`m+1`$ for bounded $`𝐔,`$ and $`𝐟_1(t,𝐱,\mathrm{𝟎})H^m(^n)`$ for any fixed $`t0`$.

Then the Cauchy problem for system (1 – 3) has locally in time a unique solution such that

$$(\mathrm{\Pi },𝐕,S,)_{j=0}^1C^j([0,T);H^{mj}(^n),$$

moreover

$$\underset{𝐔_0_{H^m}0}{lim}T=+\mathrm{}.$$

## 2 The interior solution

We call an interior solution globally in time smooth solution $`(\overline{\rho }(t,𝐱),\overline{𝐕}(t,𝐱),\overline{S}(t,𝐱))`$ to system (1–3) having the following property: another solution $`(\rho ,𝐕,S)`$ to the Cauchy problem (4) having sufficiently small norm

$$(\rho _0^{(\gamma 1)/2}(𝐱)\overline{\rho }^{(\gamma 1)/2}(0,𝐱),𝐕_0(𝐱)\overline{𝐕}(0,𝐱),S_0(𝐱)\overline{S}(0,𝐱)_{H^m(^n)}$$

is smooth globally in time as well, moreover

$$(\rho ^{(\gamma 1)/2}\overline{\rho }^{(\gamma 1)/2},𝐕\overline{𝐕},S\overline{S})_{j=0}^1C^j([0,\mathrm{});H^{mj}(^n)),m>1+n/2.$$

Note that the trivial solution is not interior if $`𝐟=\mathrm{𝟎}`$.

The set of interior solutions is not empty. In the paper  (the generalization of ) for $`𝐟=\mathrm{𝟎}`$ it is shown that the solution to system (1–3) $`(0,\overline{𝐕}(t,𝐱),const)`$ is interior, if $`\overline{𝐕}(t,𝐱)`$ is the solution to equation $`_t𝐕+(𝐕,)𝐕=0`$ such that $`SpecD\overline{𝐕}(0,𝐱)`$ is separated from the semi-axis $`_{}`$ and $`D\overline{𝐕}(0,𝐱)L^{\mathrm{}}(^n),D^2\overline{𝐕}(0,𝐱)H^{m1}(^n).`$ (We denote $`D^k`$ the vector of all spatial derivatives of order $`k`$.)

In particular, as it is noted in  in the role of such solution may be played by the solution with linear profile of velocity $`𝐕=A(t)𝐫,`$ where $`𝐫`$ is the radius-vector of point, $`A(t)`$ is the matrix with the coefficients depending on time such that $`\mathrm{Spec}A(0)_{}.`$ The solution may be found explicitly.

Note that the given interior solution (and the close solution as well) does not belong to the Sobolev space, since the velocity components grow at infinity. The result is clear from the physical point of view: the velocity field having the positive divergency ”spreads out” the initially concentrated sufficiently small mass, that prevent the singularity formation.

But the given solution with the zero density is not physical. It follows, in particular, in the case of the velocity Jacobian having initially negative spectrum at some value of $`x`$ , the appearance during a finite time of points with infinite negative divergency. If the solution with non-zero velocity, initially close to the described above, would behave like this, it would follow (according to (1)) the infinite value of density in some points.

The question arises: can one construct the interior solution with the density not close to zero in the Sobolev norm? Or to be interior the solution must describe the scattering of sufficiently small mass?

If from some reasons we succeed to find the velocity field, then from the linear with respect to $`\rho `$ and $`S`$ equations (1), (3) we can find the solution to the Cauchy problem (4), moreover, if the velocity field occurs smooth, then the rest of the components are smooth.

It needs to note, that the values $`\rho `$ and $`p`$ can become unbounded in some point of trajectory in a finite time $`T`$ if and only if $`_0^T𝑑iv𝐕𝑑t=\mathrm{}`$ on the trajectory (see (1)).

Before formulating the theorem we do some transformation under supposition that $`(\overline{\mathrm{\Pi }},\overline{𝐕}=A(t)𝐫,\overline{S})`$ is a globally in time smooth solution to system (1.3–1.5).

According to (1.3–1.5) the vector-function $`𝐮=𝐔\overline{𝐔}`$ (here $`𝐮=(\pi ,𝐯,s),`$$`𝐔=(\mathrm{\Pi },𝐕,S),`$$`\overline{𝐔}=(\overline{\mathrm{\Pi }},\overline{𝐕},\overline{S})`$) satisfies the system of equations

$$(_t+\overline{𝐕},)\pi +(𝐯,)\pi +\frac{\gamma 1}{2}(\pi +\overline{\mathrm{\Pi }})\mathrm{div}𝐕=$$

$$(𝐯,)\overline{\mathrm{\Pi }}\frac{\gamma 1}{2}\pi \mathrm{div}\overline{𝐕},$$
$`(2.1)`$
$$(_t+(\overline{𝐕},))𝐯+(𝐯,)𝐯+\frac{\gamma 1}{2}\mathrm{exp}(\frac{S}{\gamma })(\pi +\overline{\mathrm{\Pi }})\pi =$$

$$(𝐯,)\overline{𝐕}\frac{\gamma 1}{2}\mathrm{exp}(\frac{S}{\gamma })(\pi +\overline{\mathrm{\Pi }})\overline{\mathrm{\Pi }}+\frac{\gamma 1}{2}\mathrm{exp}(\frac{\overline{S}}{\gamma })\overline{\mathrm{\Pi }}\overline{\mathrm{\Pi }}+$$
$`(2.2)`$
$$+𝐟_1(t,𝐱,(\overline{𝐕}+𝐯),(\overline{S}+s),(\overline{\mathrm{\Pi }}+\pi ))𝐟_1(t,𝐱,\overline{𝐕},\overline{S},\overline{\mathrm{\Pi }}),$$

$$(_t+\overline{𝐕},)s+(𝐯,)s=(𝐯,)\overline{S}.$$
$`(2.3)`$
Further, following to , we carry out the nondegenerate change of variables such that the infinite semi-axis of time turns out to semi-interval $`[0,\sigma _{\mathrm{}})`$, to obtain the symmetric hyperbolic system with the coefficients allowing to apply the theorem on a local existence of smooth solution on the interval $`[0,\sigma _{}),`$$`\sigma _{}\sigma _{\mathrm{}}.`$

So, let $`(t_1,𝐲):=(t,A(t)𝐱)`$ be the new variables ($`A(t)=a_{ij}(t)`$ is a quadratic $`(n\times n)`$ matrix). Then $`_𝐱=A^{}_𝐲,\mathrm{div}_𝐱𝐕=\mathrm{div}_𝐲A𝐕,_{t_1}=_t+A𝐫_𝐲`$. Choose the nonnegative function $`\lambda (t)`$ such that the integral $`_0^{\mathrm{}}\lambda (\tau )𝑑\tau =\sigma _{\mathrm{}}<\mathrm{}`$ converges and set $`\sigma (t_1)=_0^{t_1}\lambda (\tau )𝑑\tau .`$ In that way, the semi-infinite axis of time goes to the semi-interval $`[0,\sigma _{\mathrm{}}).`$

Further, involve the variables

$$𝐖=A(t)𝐯\lambda ^1(t),P=\lambda ^q(t)\pi ,$$

the constant $`q`$ will be defined below.

After the transformation we get the system

$$(_\sigma +(𝐖,_𝐲))P+\frac{\gamma 1}{2}(P+\overline{P})\mathrm{div}_y𝐖=$$

$$(𝐖,_y)\overline{P}PQ_1(\sigma ),$$
$`(2.4)`$
$$(_\sigma +(𝐖,_y))𝐖+\frac{\gamma 1}{2}\mathrm{\Psi }(S,\sigma )(P+\overline{P})_yP=$$

$$\frac{\gamma 1}{2}\mathrm{\Psi }(S,\sigma )(P+\overline{P})_y\overline{P}+\frac{\gamma 1}{2}\mathrm{\Psi }(\overline{S},\sigma )\overline{P}_y\overline{P}𝐖Q_2(\sigma )+G,$$
$`(2.5)`$
$$(_\sigma +𝐖,_y)s=(𝐖,_y)\overline{S}.$$
$`(2.6)`$
Here we denote

$$Q_1(t(\sigma ))=\lambda (t)^1(\frac{\gamma 1}{2}\mathrm{tr}A(t)+q(\mathrm{ln}\lambda (t))^{}),$$

$$Q_2(t(\sigma ))=\lambda (t)^1((\mathrm{ln}\lambda (t))^{}E+A(t)((A^1)^{}(t)+E)),$$

$$\mathrm{\Psi }(S,\sigma )=\mathrm{exp}(\frac{S}{\gamma })R(\sigma )B(\sigma ),$$

$$B(t(\sigma ))=A(t)A^{}(t)(\mathrm{det}A(t))^{2/n}$$

$$R(t(\sigma ))=\lambda ^{2q2}(t)(\mathrm{det}A(t))^{2/n},$$

$$G=\lambda ^2(t)A(t)(𝐟_1(t,A(t)𝐱,\lambda ^q(t)(\overline{P}+P),\lambda (t)A^1(t)(\overline{𝐖}+𝐖),(\overline{S}+s))$$

$$𝐟_1(t,A(t)𝐱,\lambda ^q(t)\overline{P},\lambda (t)A^1(t)\overline{𝐖},\overline{S})),$$

taking into account that $`t=t(\sigma ).`$

Note that $`B`$ is a bounded invertible matrix.

After multiplying (2.5) by $`\mathrm{\Psi }^1(S,\sigma )`$ the system become symmetric hyperbolic, however, generally speaking, $`\mathrm{\Psi }^1(S,\sigma )0`$ as $`t\mathrm{}`$ or will be unbounded owing to the behaviour of $`R(\sigma (t))`$.

In another words, we can increase the time of existence of smooth solution to system (1–3) as long as we wish due to choosing the initial data small in the Sobolev $`H^m`$ norm, but we cannot guarantee the existence of the solution during the infinite time. Let

$$D^\alpha =\left(\frac{}{y_1}\right)^{\alpha _1}\mathrm{}\left(\frac{}{y_n}\right)^{\alpha _n},\underset{i}{}\alpha _{i=1}^n=\alpha .$$

Denote $`C_b^r(^n)`$ the space of continuous functions having continuous bounded in $`^n`$ derivatives up to order $`r=0,1,\mathrm{}`$. Further, denote $`𝐟_𝐕`$ the matrix $`\frac{f_i}{V_j}`$.

Now formulate the theorem taking into account all denotation involved.

###### Theorem 2. 1

Let the function $`𝐟_1(t,𝐱,\mathrm{\Pi },𝐕,S)=𝐟(t,𝐱,e^{\frac{S}{\gamma }}\left(\frac{\mathrm{\Pi }}{\kappa }\right)^{\frac{2}{\gamma 1}},𝐕,S)`$ have the derivatives with respect to all components of solution up to order $`m+1`$, continuous and bounded under bounded $`(\mathrm{\Pi },𝐕,S),`$ moreover $`f_1(t,𝐱,0,\mathrm{𝟎},0)H^m(^n)`$ at any fixed $`t0.`$ Suppose system (1–3) has the globally in time smooth solution $`\overline{U}=(\overline{\rho },\overline{𝐕},\overline{S})`$ with the linear profile of velocity $`\overline{𝐕}=A(t)𝐫`$ such that

a) $`\overline{\rho }^{\frac{\gamma 1}{2}}(0,𝐱)_{j=0}^{m+1}C_b^j(^n);`$ $`D\overline{S}(0,𝐱)_{j=0}^mC_b^j(^n);`$

b) $`\xi (t)=detA(t)>0`$ for $`tt_00;`$

and there exist a smooth real-valued function $`\lambda (t)`$, a constant $`q`$ and a skew symmetric matrix with the real-valued coefficients $`U_\varphi (t)`$ with the following properties:

c)$`_{t_0}^+\mathrm{}\lambda (\tau )𝑑\tau <\mathrm{},`$

d) $`_{t_0}^+\mathrm{}\lambda ^q(\tau )\xi ^{1/n}(\tau )𝑑\tau <\mathrm{},`$

e) functions

$$Q_1(t)R(t),$$
$`(2.7)`$
$$Q_2(t)\lambda ^1(A(t)𝐟_𝐕A^1(t)U_\varphi (t)),$$
$`(2.8)`$
$$(\mathrm{ln}R(t))^{}$$
$`(2.9)`$
are bounded on $`_+\times ^n`$.

Then the solution $`\overline{U}`$ is interior.

Proof. At first we consider the case of $`R(\sigma )`$ such that $`0<r_1=constR(\sigma )r_2=const<\mathrm{}`$. If the other conditions of Theorem 2.1 hold, then we can prove Theorem 1.1 immediately. The case is realized, for example, at $`\gamma =\frac{2}{n}+1,𝐟=0,`$ in the situation, considered in  ($`trA(t)\frac{n}{t},\xi (t)t^n(t\mathrm{}),\lambda (t)=(1+t)^2,q=\frac{n(\gamma 1)}{4}).`$

In the general case the proof is more complicated.

Suppose $`\sigma [0,\sigma _{})`$ and $`(X,𝐘,Z):R^nR\times R^n\times R`$ is a vector-function from $`L_2.`$ Following , involve the norm

$$[X,𝐘,Z]^2(\sigma ):=\underset{R^n}{}(X^2+Z^2+𝐘^{}\mathrm{\Psi }^1(S,\sigma )𝐘)𝑑𝐲.$$

For $`pN`$ define

$$E_p(\sigma )=\frac{1}{2}\underset{|\alpha |=p}{}[D^\alpha P,D^\alpha 𝐖,D^\alpha s]^2(\sigma ),$$

$$F_m(\sigma )=\underset{p=0}{\overset{m}{}}E_p(\sigma ).$$

If we should succeed to show that $`F_m(\sigma )`$ is bounded on $`[0,\sigma _{\mathrm{}}]`$ function, the theorem would be proved, as the impossibility to prolong the solution is connected with the going to the infinity of $`L_{\mathrm{}}`$ norms of the solution itself or its gradient (see (1.2)).

Showing that $`H^m`$ norm of the solution admits a majorization by means of some function, bounded on $`[0,\sigma _{\mathrm{}}]`$, we shall show that $`\sigma _{}=\sigma _{\mathrm{}}`$ and thus the theorem will be proved.

After involving of the new variable $`\mathrm{\Pi },`$ which is necessary for the symmetriization of the system (see above), globally smooth solution $`\overline{U}`$ corresponds to the globally smooth solution $`(\overline{\mathrm{\Pi }},\overline{𝐕},\overline{S})`$ with the initial data $`(\overline{\mathrm{\Pi }}_0,\overline{𝐕}_0,\overline{S}_0).`$ After the changing of variables it corresponds to the solution $`(\overline{P},\overline{𝐖},\overline{S})`$ with initial data $`(\overline{P}_0,\overline{𝐖}_0,\overline{S}_0).`$ The last solution is a smooth solution to system (2.4–2.6) on $`[0,\sigma _{\mathrm{}}]`$, and therefore its $`H^m`$ norm is bounded on the segment.

Compute

$$\frac{dE_q}{d\sigma }=\underset{|\alpha |=q}{}[\underset{R_n}{}(D^\alpha PD^\alpha _\sigma Pd𝐲+$$

$$\underset{R_n}{}D^\alpha sD^\alpha _\sigma sd𝐲+$$

$$\underset{R_n}{}D^\alpha 𝐖\mathrm{\Psi }^1(S,\sigma )D^\alpha _\sigma 𝐖d𝐲+$$

$$\underset{R_n}{}(D^\alpha 𝐖)^{}(\frac{\mathrm{\Psi }^1(S,\sigma )}{\sigma }\frac{1}{\gamma }\mathrm{\Psi }^1(S,\sigma )_\sigma S)D^\alpha 𝐖d𝐲]=$$

$$I_1+I_2+I_3+I_4$$

(further it is supposed the summation over repeated indices, the integrals $`I_k,k=1,2,3,4,`$ are numbered in consecutive order).

Estimate all integrals. Below we denote by $`c_i`$ the positive constants, depending only on initial data. Begin from $`I_4.`$ Before all we stress that according to (2.9) $`R^{}(\sigma )/R(\sigma )`$ is bounded as $`\sigma \sigma _{\mathrm{}},`$ therefore

$$|\frac{\mathrm{\Psi }^1(S,\sigma )}{\sigma }|=|\frac{\mathrm{\Psi }^2(S,\sigma )\mathrm{\Psi }(S,\sigma )R^{}}{R}+O(\mathrm{\Psi }^1(S,\sigma ))|c_1\mathrm{\Psi }^1(S,\sigma ).$$

Thus,

$$I_4c_2(1+𝐖_𝐲s_{\mathrm{}}+𝐖_𝐲\overline{S}_{\mathrm{}})E_p(c_2+c_3R^{1/2}F_m+c_4R^{1/2}F_m^{1/2})E_p.$$

Further,

$$I_1=\underset{R_n}{}D^\alpha PD^\alpha ((𝐖,_𝐲)P+\frac{\gamma 1}{2}(P+\overline{P})\mathrm{div}_𝐲𝐖+(𝐖,_𝐲)\overline{P}+PQ_1(\sigma ))𝑑𝐲.$$

The integral from the member of higher order

$$D^\alpha PD^\alpha ((𝐖,_𝐲)P+\frac{\gamma 1}{2}P\mathrm{div}_𝐲𝐖)$$

been integrated by parts can be reduced to

$$\underset{R_n}{}\left\{\frac{1}{2}\right|D^\alpha P|\mathrm{div}_𝐲𝐖+D^\alpha P\left((𝐖,_𝐲(D^\alpha P))D^\alpha (𝐖,_𝐲P)\right)\}𝑑𝐲$$

$$\frac{\gamma 1}{2}\underset{R_n}{}D^\alpha P\left\{(_𝐲P,D^\alpha 𝐖)+(P\mathrm{div}_𝐲D^\alpha 𝐖D^\alpha (P\mathrm{div}_𝐲𝐖))\right\}𝑑𝐲:=I_{11}.$$

Under sign of the integral there is the sum of members of the form $`^\alpha U_1^lU_2^{\alpha +1l}U_2,\mathrm{\hspace{0.17em}0}l\alpha `$ (by $`U_j,j=1,2,`$ we mean the corresponding components of the solution. At first let $`\alpha 0,1.`$ According to the Galiardo-Nirinberg inequality we have

$$|^iU_j|_{p_i}C|DU_j|_{\mathrm{}}^{12/p_i}|D^\alpha U_j|_2^{2/p_i}$$

for $`p_i=2\frac{\alpha 1}{i1}.`$ For $`l1,\alpha `$ it is true that $`\frac{1}{p_l}+\frac{1}{p_{\alpha l+1}}=\frac{1}{2}.`$ From the H$`\ddot{\mathrm{o}}`$lder inequality it follows that

$$\underset{R_n}{}|^\alpha U_1^lU_2^{\alpha +1l}U_2|𝑑𝐲$$

$$D^\alpha U_1_2^lU_2^{\alpha +1l}U_2_2$$

$$D^\alpha U_1_2^lU_2_{p_l}^{\alpha +1l}U_2_{p_{\alpha l+1}}CD^\alpha U_1_2D^\alpha U_2_2D^\alpha U_2_{\mathrm{}}.$$

For the other values of $`\alpha `$ and $`l`$ the last inequality is evident.

Thus,

$$I_{11}cD^pP_2(D^pP_2_𝐲𝐖_{\mathrm{}}+D^p𝐖_2_𝐲P_{\mathrm{}})cR^{1/2}F_m^{1/2}E_p.$$

(one can find the details of the inequalities applying in the book .)

Further, according to condition (2.7) of the Theorem

$$\left|\underset{R_n}{}D^\alpha PD^\alpha (𝐖_𝐲\overline{P})𝑑𝐲\right|c_5R^{1/2}F_m,$$

$$\left|\underset{R_n}{}D^\alpha PD^\alpha (PQ_1)𝑑𝐲\right|c_6Q_1RF_mc_7F_m.$$

The integral $`I_2`$ can be estimated analogously:

$$|I_2|c_8R^{1/2}F_m^{1/2}E_q+c_9R^{1/2}F_m.$$

At last,

$$I_3=\underset{R_n}{}(D^\alpha )^{}𝐖\mathrm{\Psi }^1(S,\sigma )\{(𝐖,_𝐲)D^\alpha 𝐖D^\alpha ((𝐖,_𝐲𝐖)\}d𝐲+$$

$$\frac{1}{2}\underset{R_n}{}\{(D^\alpha 𝐖)^{}\mathrm{\Psi }^1(S,\sigma )D^\alpha 𝐖)\mathrm{div}_𝐲𝐖\frac{𝐖_𝐲S}{\gamma }(D^\alpha 𝐖)^{}\mathrm{\Psi }^1(S,\sigma )D^\alpha 𝐖\}d𝐲+$$

$$\frac{\gamma 1}{2}\underset{R_n}{}\{(D^\alpha 𝐖)^{}\mathrm{\Psi }^1(S,\sigma )D^\alpha (P\mathrm{\Psi }(S,\sigma )_𝐲P+\overline{P}\mathrm{\Psi }(S,\sigma )_𝐲P$$

$$(P+\overline{P})\mathrm{\Psi }(S,\sigma )_𝐲\overline{P}+\overline{P}\mathrm{\Psi }(\overline{S},\sigma )_𝐲\overline{P})\}d𝐲+$$

$$\underset{R_n}{}(D^\alpha 𝐖)^{}\mathrm{\Psi }^1(S,\sigma )D^\alpha \left(𝐖Q_2+G\right)𝑑𝐲.$$

The first and second integrals can be estimated from above by values

$$c_{10}R^1D^p𝐖_2^2𝐖_{\mathrm{}}+c_{11}R^1D^p𝐖_2^2𝐖_2)c_{12}R^{1/2}F_m^{1/2}E_p,$$

and the third one by the value

$$c_{13}D^p𝐖_2(D^pP_2+D^ps_2+D^p\overline{P}_2+D^{p1}_𝐲\overline{S}_2)(1+P_{\mathrm{}}+_𝐲P_{\mathrm{}}+$$

$$_𝐲s+_{\mathrm{}}\overline{P}_{\mathrm{}}+_𝐲\overline{P}_{\mathrm{}}+_𝐲\overline{S}_{\mathrm{}})^{p+1}$$

$$c_{13}R^{1/2}(c_{14}+F_m^{(p+1)/2})E_p+c_{15}R^{1/2}(c_{16}+F_m^{(p+1)/2})E_p^{1/2},$$

and, at last, the fourth one taking into account conditions (2.8) by the value

$$c_{17}E_p.$$

In the last case note that

$$G=\lambda ^1(t)A(t)𝐟_𝐕(t,𝐱,\mathrm{\Pi },\overline{𝐕}+\theta 𝐯,S)A^1(t),\theta (0,1),$$

$$\underset{R_n}{}(D^\alpha 𝐖)^{}\mathrm{\Psi }^1(S,\sigma )D^\alpha (U_\varphi (t)𝐖)𝑑𝐲=0,$$

as $`(𝐰,U_\varphi (t)𝐰)=0`$ for any vector $`𝐰^n.`$

Uniting all the estimates we get

$$E_p^{}c_{18}E_p+c_{19}R^{1/2}E_p+c_{20}R^{1/2}F_m^{1/2}E_p+c_{21}R^{1/2}F_mE_p+$$

$$c_{22}R^{1/2}F_m^{(p+1)/2}E_p+c_{23}R^{1/2}F_m^{(p+1)/2}E_p^{1/2},$$

whence it follows after summation over $`p`$ that

$`F_m^{}c_{24}F_m+c_{25}R^{1/2}(F_m^{1/2}+F_m+F_m^{3/2}+F_m^{(p+2)/2}+F_m^{(p+3)/2}).`$

Set $`\mathrm{\Lambda }_m(\sigma )=\mathrm{exp}(c_{24}\sigma )F_m.`$ Then

$$\mathrm{\Lambda }_m^{}c_{26}(\mathrm{\Lambda }_m^{1/2}+\mathrm{\Lambda }_m+\mathrm{\Lambda }_m^{3/2}+\mathrm{\Lambda }_m^{(p+2)/2}+\mathrm{\Lambda }_m^{(p+3)/2})R^{1/2},$$
$`(2.10)`$
where the constant $`c_{26}`$ depends on $`c_{24},c_{25},\sigma _{\mathrm{}}`$.

We set

$$\mathrm{\Theta }(\stackrel{~}{g}):=\underset{\delta >0}{\overset{\stackrel{~}{g}}{}}\frac{dg}{g^{1/2}+g+g^{3/2}+g^{(p+2)/2}+g^{(p+3)/2}}.$$

The integral diverges in the zero, therefore $`\mathrm{\Theta }(0)=\mathrm{}.`$

Integrating inequalities (2.10) over $`\sigma `$ we obtain that

$$\mathrm{\Theta }(\mathrm{\Lambda }_m(\sigma ))\mathrm{\Theta }(\mathrm{\Lambda }_m(t_0))+c_{26}\underset{\sigma (t_0)}{\overset{\sigma }{}}R^{1/2}(\stackrel{~}{\sigma })𝑑\stackrel{~}{\sigma },$$

moreover, as follows from (2.8), the integral in the right-hand side of the last inequality converges as $`\sigma \sigma _{\mathrm{}}`$ to the constant $`C`$, depending only on the initial data.

Choosing $`\mathrm{\Lambda }_m(t_0),`$ (and together with the value the $`H_m`$ norm of the initial data) sufficiently small, one can obtain that $`\mathrm{\Theta }(\mathrm{\Lambda }_m(t_0))+C`$ is later then $`\mathrm{\Theta }(+\mathrm{}),`$ it signifies that $`\mathrm{\Lambda }_m(\sigma )`$ and $`F_m(\sigma )`$ are bounded from above for all $`\sigma [0,\sigma _{\mathrm{}})`$ and $`\sigma _{}=\sigma _{\mathrm{}}`$. So, Theorem 2.1 is proved.

Remark.  The solutions with the linear profile of velocity, satisfying the theorem conditions, exist. In the case described in  $`f=0,\gamma 1+\frac{2}{n},`$ $`A(t)=(E+A(0)t)^1A(0),\mathrm{Spec}A(0)_{},`$ $`\mathrm{tr}A(t)\frac{n}{t},`$$`\mathrm{det}A(t)\frac{1}{t^n},`$ $`Q_1=0,`$$`Q_2`$$`2(\frac{1}{n}((\mathrm{tr}A^1)EA^1),`$ as $`t\mathrm{}.`$ Here $`\lambda =(1+t)^2,`$$`q=\frac{n(\gamma 1)}{4},`$ $`U_\varphi (t)=0.`$

We shall get two corollaries from Theorem 2.1, basing on which we can assert that the solutions which will be constructed in the next section are interiors.

###### Corollary 2. 1

Let $`𝐟=\mathrm{𝟎}.`$ If system (1–3) has the globally smooth in time solution with the linear profile of velocity, described in the statement of Theorem 2.1 and $`A(t)\frac{\delta }{t}E,t\mathrm{},`$ where $`\delta `$ is a positive constant, then the solution is interior.

Proof. Choose $`\lambda (t)=t^{\delta +1},`$ $`q\frac{3}{2}\frac{n+1}{2(\delta +1)},`$ $`U_\varphi (t)=0.`$ In that way all conditions of Theorem 2.1 are satisfied.

###### Corollary 2. 2

Let $`𝐟=L𝐕,`$ where $`L`$ is a matrix with the smooth coefficients such that $`A_0LA_0^1U_1(t)=\mu E,`$ where $`A_0=\delta E+U_2`$ is a matrix with the positive determinant, the matrices $`U_1`$ and $`U_2`$ are skew symmetric, $`\delta `$ is a constant, $`\mu `$ is a positive constant. If system (1–3) has the globally smooth in time solution with the linear profile of velocity described in the statement of Theorem 2.1 and $`A(t)\frac{1}{t}A_0,t\mathrm{},`$ then the solution is interior.

Proof. To verify condition (2.8) one has to show the boundedness of the value

$$\lambda ^1(t)((\mathrm{ln}\lambda (t))^{}E+A(t)EA(t)LA^1(t)A^1(t)A^{}(t)U_\varphi )$$

$$\lambda ^1(t)((\mathrm{ln}\lambda (t))^{}E+\frac{1}{t}((\delta +1)E+U_2)U_1(t)+\mu EU_\varphi ),t\mathrm{}.$$

Choose $`\lambda (t)=t^{\delta +1}\mathrm{exp}\{\mu t\}.`$ Condition (2.7) is satisfied at $`q\frac{3}{2},`$ condition (2.9) can be verified elementary. It needs to choose the matrix $`U_\varphi (t)=\frac{1}{t}U_2U_1(t).`$

## 3 Constructing of interior solutions

Limit ourselves to the important case of $`n=2`$ and set $`𝐟=L𝐯,`$ where the matrix $`L=\left(\begin{array}{cc}\mu & l\\ l& \mu \end{array}\right),`$ $`\mu `$ is a nonnegative constant, $`l`$ is an arbitrary constant. Thus, we in the simplest way describe the Coriolis force and the Rayleigh friction, for example, in the meteorological model under neglecting the vertical processes and the Earth curvature.

It is easy to verify that on the smooth solutions to system (1–3), so quickly vanishing at infinity to assure the convergency of all integrals involved, there are the following conservation quantities: the mass $`m=\underset{^2}{}\rho 𝑑𝐱=const,`$ for $`\mu =0`$ the total energy $`E=\underset{^2}{}\left(\frac{\rho |𝐕|^2)}{2}+\frac{p}{\gamma 1}\right)𝑑𝐱=E_k(t)+E_p(t)=const`$ and the momentum of inertia $`J=\underset{^2}{}\rho ((𝐕_{},𝐫)+`$ $`\frac{l}{2}|𝐱|^2)d𝐱=const.`$ The conservation laws are true also for the solutions to system (1–3) with the density (and the pressure) quickly vanishing as $`|𝐱|\mathrm{}`$, and the velocity components may even increase.

We consider another integral functionals, characterizing the average properties of solutions, namely

$$G(t)=\frac{1}{2}\underset{^2}{}\rho |𝐫|^2𝑑𝐱,F_i(t)=\underset{^2}{}(𝐕,𝐗_i)\rho 𝑑𝐱,i=1,2,$$

where $`𝐗_1=𝐫=(x,y),𝐗_2=𝐫_{}=(y,x).`$

Note that $`G(t)>0.`$

On smooth solutions to system (1 – 3) the following relations take place

$$G^{}(t)=F_1(t),$$
$`(3.1)`$
$$F_2^{}(t)=lF_1(t)\mu F_2(t),$$
$`(3.2)`$
$$F_1^{}(t)=2(\gamma 1)E_p(t)+2E_k(t)lF_2(t)\mu F_1(t),$$
$`(3.3)`$
$$E^{}(t)=2\mu E_k(t).$$

Remark once more that at $`\mu =0`$ the quantity $`E(t)`$ of total energy is constant. The functions $`E_p(t)`$ and $`E_k(t)`$ generally cannot be expressed through $`G(t),F_1(t),F_2(t)`$.

We try to choose the form of velocity field in such way that from system (3.1 – 3.3) one could find it directly. Suppose, for example,

$$𝐕=A(t)𝐫,$$
$`(3.4)`$
where $`A(t)`$ is a (2 2) matrix , whose time-depending coefficients we have to find.

1) We get the simplest result if we choose the matrix $`A(t)`$ such that

$$𝐕=\alpha (t)𝐫+\beta (t)𝐫_{},$$

It is easy to see that in the case we have

$$F_1(t)=2\alpha (t)G(t),F_2(t)=2\beta (t)G(t),$$

$$E_k(t)=(\alpha ^2(t)+\beta ^2(t))G(t).$$

It is not difficult to obtain from the equation

$$_tp+(𝐕,p)+\gamma pdiv𝐕=0,$$
$`(3.5)`$
which is a corollary of system (1–3) and the state equation, that

$$E_p(t)=E_p(0)G^{\gamma 1}(0)\frac{1}{G^{\gamma 1}(t)}.$$

In that way, all functions involved in system (3.1 – 3.3) are expressed through $`G(t),\alpha (t),\beta (t).`$ For the sake of convenience we denote $`G_1(t)=1/G(t)`$ and get the system of equations

$$G_1^{}(t)=2\alpha (t)G_1(t),$$
$`(3.5)`$
$$\beta ^{}(t)=\alpha (t)(l2\beta (t))\mu \beta (t),$$
$`(3.6)`$
$$\alpha ^{}(t)=\alpha ^2(t)+\beta ^2(t)l\beta (t)\mu \alpha (t)+(\gamma 1)E_p(0)G_1^{1\gamma }(0)G_1^\gamma (t).$$
$`(3.7)`$
The solutions of the system of ordinary differential equations are smooth if they are can be prolonged up to all axis of time. The non-prolongation is connected with the escape of the components of solution to infinity during a finite time.

Note that for $`\gamma >1`$ the quantities $`\alpha (t)`$ and $`\beta (t)`$ are bounded, this follows immediately from the expression for the total energy of the system:

$$E(t)=(\alpha (t)^2+\beta (t)^2)G(t)+E_p(0)G^{\gamma 1}(0)\frac{1}{G^{\gamma 1}(t)}E(0),$$

therefore

$$\alpha ^2(t)+\beta ^2(t)E(0)G_1(t)E_p(0)G_1^{1\gamma }(0)G_1^\gamma (t)<+\mathrm{},$$

thus the quantities of density and pressure remain bounded for all solution we want to construct.

At $`\mu =0`$ system (3.5–3.7) can be integrated explicitly:

$$\alpha (G_1)=\pm \sqrt{KG_1^\gamma C^2G_1^2+(ElC)G_1l^2/4},$$

$$\beta (t)=CG_1(t)+l/2,$$
$`(3.8)`$
$$_{G_1(0)}^{G_1(t)}\frac{dG_1}{2G_1\alpha (G_1)}=t,$$

with the constants $`C=\frac{2\beta (0)l}{2G_1(0)},K=(\alpha ^2(0)+C^2G_1^2(0)(ElC)G_1(0)+l^2/4)/G_1^\gamma (0).`$

Due to (3.8) one can reduce the number of equations in the system and consider it on the phase plane $`(G_1,\alpha ),G_1>0`$. There exists the unique singular point: at $`l0`$ it is a center situated on the axis $`\alpha =0`$, at $`l=0`$ it is a complex equilibrium in the origin, the trajectories form the elliptic saddle point. In the last case the time of movement from any point of the phase plane to the origin is infinite. The equilibria, at $`l0`$ lying in the plane $`G_1=0`$, namely, the points $`\alpha =\beta =0`$ and $`\alpha =0,\beta =l,`$ are the centers under consideration in the plane $`(\alpha ,\beta ).`$

At $`\mu 0`$ one needs to integrate the system (3.5–3.7) numerically. However, the conclusion on the behaviour of the solution components at infinity we can do analytically as well. At $`\mu >0`$ the system has a unique stable equilibrium in the origin (at $`l=0`$ it is a knot, at $`l0`$ it is a focus.) Therefore we can find the solution as a formal asymptotic series by the negative powers of $`t`$. Limit ourselves to writing out the leading terms:

$$\alpha (t)a_1t^p,\beta (t)b_1t^r,G_1(t)c_1t^q,a_1,b_1,c_1=const0.$$

From (3.5) we have immediately that $`qc_1t^{q1}=2a_1c_1t^{p+q},`$ therefore $`p=1,q=2a_1.`$

From (3.6) we obtain that $`rb_1t^{r1}la_1t^p\mu b_1t^r2a_1b_1t^{p+r}.`$ If $`l0,\mu 0,`$ then $`r=p=1,`$ moreover, $`la_1=\mu b_1.`$

From (3.7) we get analogously that $`pa_1t^{p1}a_1^2t^{2p}+b_1^2t^{2r}lb_1t^r\mu a_1t^p+Kc_1^\gamma t^{\gamma q},`$ where it is denoted for short $`K=(\gamma 1)E_p(0)G_1^{1\gamma }(0),`$ that is $`a_1t^2a_1^2t^2+b_1^2t^2lb_1t^1\mu a_1t^1+Kc_1^\gamma t^{\gamma q}.`$

Let $`\gamma q<1,`$ that is the last term does not contain the senior degree. Then the rest of members of senior order must be eliminated, and $`lb_1=\mu a_1=\mu ^2b_1/l,`$ therefore $`l=\mu =0`$ in spite of the supposition. Consequently, $`q=1/\gamma .`$

So, at $`\mu 0,l0`$

$$\alpha (t)\frac{1}{2\gamma }t^1,\beta (t)\frac{\mu }{2l\gamma }t^1,G_1(t)\left(\frac{\mu ^2+l^2}{2K\gamma \mu }\right)^{1/\gamma }t^{1/\gamma },t\mathrm{}.$$

At $`l=0,`$ due to the existence of the integral $`\beta (t)=CG_1(t)\mathrm{exp}\{\mu t\},`$ system (3.5 – 3.7) can be reduced to the nonautonomous system of two equations:

$$G_1^{}(t)=2\alpha (t)G_1(t),$$

$$\alpha ^{}(t)=\alpha ^2(t)\mu \alpha (t)+C^2\mathrm{exp}\{2\mu t\}G_1^2(t)+KG_1^\gamma (t).$$
$`(3.9)`$
From (3.9) we have that $`pa_1t^{p1}a_1^2t^{2p}\mu a_1t^p+C^2c_1^2t^{2q}\mathrm{exp}\{2\mu t\}+Kc_1^\gamma t^{\gamma q},`$ where, remember, $`p=1,q=2a_1.`$ To compensate the senior term $`\mu a_1t^1`$ it is necessary that the equality $`q=1/\gamma `$ holds. Hence, we find the coefficient $`c_1.`$

Thus, at $`\mu 0,l=0`$

$$\alpha (t)\frac{1}{2\gamma }t^1,\beta (t)C\left(\frac{\mu }{2K\gamma }\right)^{1/\gamma }t^{1/\gamma }\mathrm{exp}\{\mu t\},G_1\left(\frac{\mu }{2K\gamma }\right)^{1/\gamma }t^{1/\gamma },t\mathrm{}.$$

At $`\mu =l=0,`$ as follows from the explicit formulas,

$$\alpha (t)\left(\sqrt{G_1(0)/E}+t\right)^1,G_1(t)\alpha ^2(t),\beta (t)=CG_1(t),t\mathrm{}.$$

Thus, as follows from Corollaries 2.1 and 2.2, in the cases $`\mu =l=0`$ and $`\mu >0`$ we have constructed the velocity field for the interior solution. More precisely, at $`\mu =l=0`$ $`A(t)\frac{1}{t}E,t\mathrm{},`$ Corollary 2.1 can be applied with $`\delta =1.`$ At $`\mu >0`$ one can apply Corollary 2.2. In the case $`\mu >0,l=0`$ we have $`U_1=U_2=0,`$ $`\delta =\frac{1}{2\gamma }`$. If $`\mu >0,l0,`$ then $`\delta =\frac{1}{2\gamma },`$ $`U_1=\left(\begin{array}{cc}0& l\\ l& 0\end{array}\right),`$ $`U_2=\frac{\mu }{2\gamma l}\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right).`$

As soon as $`\alpha (t)`$ and $`\beta (t)`$ are found, the rest of the components can be found elementary

$$\rho (t,|r|,\varphi )=\mathrm{exp}(2\underset{0}{\overset{t}{}}\alpha (\tau )𝑑\tau )\rho _0(|r|\mathrm{exp}(\underset{0}{\overset{t}{}}\alpha (\tau )𝑑\tau ),\varphi +\underset{0}{\overset{t}{}}\beta (\tau )𝑑\tau ),$$

$$S(t,|r|,\varphi )=S_0(|r|\mathrm{exp}(\underset{0}{\overset{t}{}}\alpha (\tau )𝑑\tau ),\varphi +\underset{0}{\overset{t}{}}\beta (\tau )𝑑\tau ).$$

$$p(t,|r|,\varphi )=\mathrm{exp}(2\gamma \underset{0}{\overset{t}{}}\alpha (\tau )𝑑\tau )p_0(|r|\mathrm{exp}(\underset{0}{\overset{t}{}}\alpha (\tau )𝑑\tau ),\varphi +\underset{0}{\overset{t}{}}\beta (\tau )𝑑\tau ).$$

From (2) and (3.7) we get that on the classical solution of the initial system must be satisfied the relation

$$p=(\gamma 1)G_1^{1\gamma }(0)E_p(0)G_1^\gamma (t)\rho 𝐫.$$

Hence it follows that the components of the initial data $`\rho _0`$ $`p_0`$ must be axisymmetric and compatible, i.e. connected as follows:

$$p_0=(\gamma 1)G_1(0)E_p(0)\rho _0𝐫.$$

However, in spite of the components of the pressure and the density must vanish as $`|𝐱|\mathrm{},`$ the entropy must even increase (remember that the conditions to Theorem 2.1 require only the boundedness of the entropy gradient). However, one can choose

$$p_0=\frac{1}{(1+|𝐫|^2)^a},a=const>3,\rho _0=\frac{2a}{(\gamma 1)G_1(0)E_p(0)}\frac{1}{(1+|𝐫|^2)^{a+1}},$$

then $`S_0=const+(a(\gamma 1)+\gamma )\mathrm{ln}(1+|𝐫|^2).`$

Remark. In the case $`G(0)0`$ the density cannot be compactly supported (in contrast $`G(0)=0`$). Really, since $`p=\pi ^{\frac{2\gamma }{\gamma 1}},\rho =\pi ^{\frac{2}{\gamma 1}}e^{\frac{S}{\gamma }},`$ then due to the compatibility conditions $`\pi \pi conste^{\frac{S}{\gamma }},|𝐱|c0,`$ where $`c`$ is a point of the support of $`\pi .`$ Therefore, for $`C^1`$ smooth $`\pi `$ it occurs that $`S+\mathrm{},|𝐱|c0,`$ and we cannot choose any smooth initial data. It is interesting that if one requires only $`C^0`$ smoothness of $`\pi (\pi const(c|𝐱|)^{1/2},|𝐱|c0,\pi =0,|𝐱|c>0),`$ the condition may be fulfilled. Moreover, $`\rho `$ and $`p`$ will be of the $`C^1`$ class of smoothness, however, neither the theorem on the local in time existence of the smooth solution, no Theorem 2.1 can be applied.

2) To consider the velocity field (3.4) with the matrix

$$A(t)=\left(\begin{array}{cc}a(t)& b(t)\\ c(t)& d(t)\end{array}\right),$$

involve the functions

$$G_x(t)=\frac{1}{2}\underset{^2}{}\rho x^2𝑑𝐱,G_y(t)=\frac{1}{2}\underset{^2}{}\rho y^2𝑑𝐱,G_{xy}(t)=\frac{1}{2}\underset{^2}{}\rho xy𝑑𝐱.$$

Then involve the auxiliary variables

$$G_1(t)=G_x(t)\mathrm{\Delta }^{(\gamma +1)/2}(t),G_2(t)=G_y(t)\mathrm{\Delta }^{(\gamma +1)/2}(t),G_3(t)=G_{xy}\mathrm{\Delta }^{(\gamma +1)/2}(t),$$
$`(3.10)`$
where $`\mathrm{\Delta }(t)=G_x(t)G_y(t)G_{xy}^2(t)`$ is a positive function on solutions to system (1–3). Note that the behaviour of $`\mathrm{\Delta }(t)`$ is governed by the equation

$$\mathrm{\Delta }^{}(t)=2(a(t)+d(t))\mathrm{\Delta }(t),$$
$`(3.11)`$
and a potential energy $`E_p(t)`$ is connected with $`\mathrm{\Delta }(t)`$ as follows:

$$E_p(t)=E_p(0)\mathrm{\Delta }^{(\gamma 1)/2}(0)\mathrm{\Delta }^{(\gamma +1)/2}(t).$$

To find $`G_1,G_2,G_3`$, and also the elements of the matrix $`A(t)`$ we get the system of equations

$$G_1^{}(t)=((1\gamma )a(t)(1+\gamma )d(t))G_1(t)+2b(t)G_3(t),$$

$$G_2^{}(t)=((1\gamma )d(t)(1+\gamma )a(t))G_2(t)+2c(t)G_3(t),$$

$$G_3^{}(t)=c(t)G_1(t)+b(t)G_2(t)\gamma (a(t)+d(t))G_3(t),$$

$$a^{}(t)=a^2(t)b(t)c(t)+lc(t)\mu a(t)+K_1G_2(t),$$
$`(3.12)`$
$$b^{}(t)=b(t)(a(t)+d(t))+ld(t)\mu b(t)K_1G_3(t),$$

$$c^{}(t)=c(t)(a(t)+d(t))la(t)\mu c(t)K_1G_3(t),$$

$$d^{}(t)=d^2(t)b(t)c(t)+lb(t)\mu d(t)+K_1G_1(t),$$

with $`K_1=\frac{\gamma 1}{2}(E_p(0)\mathrm{\Delta }^{(\gamma 1)/2}(0))`$.

###### Proposition 3. 1

In the case $`\mu =l=0`$ any solution to system (1–3) with a linear profile of velocity, satisfying the condition ) of Theorem 2.1, is interior.

Proof. Show that at $`\mu =l=0`$ there is the asymptotics $`A(t)\frac{\delta }{t}E,t\mathrm{},`$ that is according to Corollary 2.1 the corresponding solution with a linear profile of velocity is interior.

Go to the new variables $`a_1(t)=a(t)d(t),b_1(t)=b(t)+c(t),c_1(t)=b(t)c(t),d_1(t)=a(t)+d(t),G_4(t)=G_1(t)+G_2(t),G_5(t)=G_1(t)G_2(t).`$ In variables $`a_1,b_1,c_1,d_1,G_3,G_4,G_5`$ system (3.12) has a form:

$$a_1^{}(t)=a_1(t)d_1(t)K_1G_5(t),$$
$`(3.13)`$
$$b_1^{}(t)=b_1(t)d_1(t)2K_1G_3(t),$$
$`(3.14)`$
$$c_1^{}(t)=c_1(t)d_1(t),$$
$`(3.15)`$
$$d_1^{}(t)=\frac{1}{2}(a_1^2(t)+d_1^2(t))\frac{1}{2}(b_1^2(t)c_1^2(t))+K_1G_4(t),$$
$`(3.16)`$
$$G_3^{}(t)=\gamma d_1(t)G_3(t)+\frac{1}{2}b_1(t)G_4(t)\frac{1}{2}c_1(t)G_5(t),$$
$`(3.17)`$
$$G_4^{}(t)=\gamma d_1(t)G_4(t)+a_1(t)G_5(t)+2b_1(t)G_3(t),$$
$`(3.18)`$
$$G_5^{}(t)=\gamma d_1(t)G_5(t)+a_1(t))G_4(t)2c_1(t)G_3(t).$$
$`(3.19)`$
Let as $`t\mathrm{}`$ the asymptotics of functions involved in the system be the following: $`a_1(t)L_1t^{l_1},`$ $`b_1(t)L_2t^{l_2},c_1(t)L_3t^{l_3},`$$`d_1(t)L_4t^{l_4},`$$`G_3(t)Nt^q,G_4(t)M_1t^{p_1},`$ $`G_5(t)M_2t^{p_2},`$ where $`L_i,i=1,2,3,M_j,j=1,2,N`$ are some constants not equal to zero. Note that $`p_2p_1`$ and in virtue of $`\mathrm{\Delta }>0`$ the estimate $`qp_1`$ holds.

From (3.15) we get immediately that

$$l_3L_3t^{l_31}=L_3L_4t^{l_3+l_4},$$

hence $`l_4=1,L_4=l_3.`$

From (3.13) we have taking it into account

$$l_1L_1t^{l_11}=L_1L_4t^{l_11}2K_1M_2t^{p_2}.$$

The following variants are possible:

$`p_2<l_11,`$ hence $`l_1=L_4,`$ i.e. $`l_1=l_3,`$

$`p_2=l_11,`$ hence

$$l_1=L_4\frac{K_1M_2}{L_1}=l_3\frac{K_1M_2}{L_1}.$$
$`(3.20)`$
From (3.14) we have analogously

$$l_2L_2t^{l_21}=L_2L_4t^{l_21}2K_1Nt^q.$$

There are the variants:

$`q<l_21,`$ hence $`l_2=L_4,`$ i.e. $`l_2=l_3,`$

$`q=l_21,`$ hence

$$l_2=L_4\frac{2K_1N}{L_2}=l_3\frac{2K_1N}{L_2}.$$
$`(3.21)`$
From (3.16) we get

$$L_4t^2=\frac{1}{2}(L_1^2t^{2l_1}+L_4^2t^2)\frac{1}{2}(L_2^2t^{2l_2}L_3^2t^{2l_3})+K_1M_1t^{p_1}.$$
$`(3.22)`$
If $`l_i<1,i=1,2,3,`$ then at $`t\mathrm{}`$

$$a(t)\frac{1}{2}(L_1t^{l_1}+L_4t^1),b(t)\frac{1}{2}(L_2t^{l_2}+L_3t^{l_3}),$$

$$c(t)\frac{1}{2}(L_2t^{l_2}L_3t^{l_3}),d(t)\frac{1}{2}(L_4t^1L_1t^{l_1}),$$

that is $`A(t)\frac{L_4}{2t}E,L_4=const>0.`$

We shall show below that there are not others situations.

From (3.22) it follows the impossibility of the situation $`l_3=1,l_1<1,l_2<1,p_12.`$ Really, in the case $`l_3=L_4=1\sqrt{1+L_3^2}<1`$ or $`l_3=L_4=1\sqrt{1+L_3^2+2K_1M_1}<1,`$ it contradicts to the supposition.

Now let $`l_3>1.`$ Then from (3.22) $`l_1=l_3`$ or (and) $`l_2=l_3.`$ Suppose, for example, that $`l_1=l_3.`$ Then $`p_2<l_11,`$ from (3.19) it follows that

$$p_2M_2t^{p_21}=\gamma L_4M_2t^{p_21}+L_1M_1t^{p_1+l_1}+2L_2Nt^{q+l_2},$$
$`(3.23)`$
and therefore $`p_1+l_1=q+l_2,l_1ł_2,`$ i.e. $`l_2>1.`$ Besides,

$$L_1M_1=2L_2N.$$
$`(3.24)`$
From (3.17) we have

$$qNt^{q1}=\gamma L_4Nt^{q1}+\frac{1}{2}L_2M_1t^{p_1+l_2}\frac{1}{2}L_3M_2t^{p_2+l_3}.$$
$`(3.25)`$
As soon as $`q1p_11,`$ and $`p_11<p_1+l_2,`$ then $`p_1+l_2=p_2+l_3p_1+l_3,`$ and therefore $`l_2l_3,`$ and as $`l_1=l_3,`$ then $`l_1=l_2=l_3`$ and $`p_1=q.`$ Besides,

$$L_2M_1=L_3M_2.$$
$`(3.26)`$
Further, from (3.18) we obtain

$$p_1M_1t^{p_11}=\gamma L_4M_1t^{p_11}+L_1M_2t^{p_2+l_1}2L_3Nt^{q+l_3}.$$
$`(3.27)`$
As $`q+l_3>p_11,`$ then $`p_2+l_1=q+l_3,p_2=q,`$

$$L_2M_1=2L_3N.$$
$`(3.28)`$
From (3.24) and (3.28) we have $`\frac{M_1}{M_2}=\frac{L_2}{L_3},`$ and from (3.26) we obtain $`\frac{M_1}{M_2}=\frac{L_2}{L_3},`$ i.e. $`L_2^2=L_3^2,L_2=L_3=0`$ in spite of the supposition.

Suppose now that $`l_2=l_3.`$ Then $`q<l_21,p_2l_11.`$ From (3.17) we have that $`q1<l_2+p_1,`$ since $`l_2>1,`$ and therefore $`l_2+p_1=l_3+p_2,p_1=p_2,L_2M_1=L_3M_2.`$

Taking this into account from (3.23) we have that $`p_1+l_1=q+l_2p_1+l_2,l_1l_2,L_2M_1=2L_2N.`$

From (3.27) we get $`p_2+l_1=q+l_3p_1+l_2,p_2=q,l_1=l_2,L_1M_2=2L_3N.`$ In that way, we obtain the contradiction analogous to the previous one.

Now let $`l_31,`$ $`l_1>1`$ and (or) $`l_2>1.`$ Then from (3.22) it follows that $`p_1=2l_1`$ and (or) $`p_1=2l_2(p_1>2).`$ For example, if $`p_1=2l_1,`$ then from (3.24) we get $`p_21<p_1+l_1,p_1+l_1=q+l_2p_1+l_2,l_1l_2.`$ From (3.25) $`p_2+l_3p_11,l_2+p_1l_1+p_1>p_11,`$ therefore $`q1=l_2+p_1.`$ But $`q1p_11,`$ therefore $`l_21,`$ and $`l_21`$ in spite of the supposition.

If $`p_1=2l_2,`$ then from (3.25) we get $`l_3+p_2p_11,`$ therefore $`p_1+l_2>=p_2+l_3,q1=p_1+l_2,l_1l_2.`$ But $`q1p_11,`$ therefore $`l_21`$ in spite of the supposition.

So, it remains the unique possibility: $`l_31,l_11`$ and (or) $`l_21.`$ In the case, as follows from (3.22)

$$l_3=1\pm \sqrt{1(\delta _1L_1^2+\delta _2L_2^2\delta _3L_3^22\delta _4K_1M_1)},$$
$`(3.29)`$
where $`\delta _i=1,`$ if $`l_i=1,`$ and $`\delta _i=0`$ otherwise, $`i=1,2,3,`$ $`\delta _4=1,`$ if $`p_2=2`$ and $`\delta _4=0`$ otherwise. Hence

$$\delta _1L_1^2+\delta _2L_2^2\delta _3L_3^22\delta _4K_1M_11,$$
$`(3.30)`$
if inequality (3.30) is strict, then $`l_3<1.`$

We consider this case, that is $`l_3<1,l_1=1`$ and (or) $`l_2=1.`$

If $`l_1=1ł_3,`$ then $`p_2=l_11=2.`$ If $`l_2=1ł_3,`$ then $`q=l_21=2.`$ In that way, in any of these cases we have $`p_1=2.`$

At $`l_1=1,l_2<1,l_3<1,`$ from (3.23) we obtain $`p_2=\gamma L_4+\frac{L_1M_1}{M_2},`$

$$\frac{L_1M_1}{M_2}=(2+\gamma l_3),$$
$`(3.31)`$
from (3.27) we get $`p_1=\gamma L_4+\frac{L_1M_2}{M_1},`$

$$\frac{L_1M_2}{M_1}=(2+\gamma l_3),$$
$`(3.32),`$
and from (3.21)

$$\frac{K_1M_2}{L_1}=l_3+1.$$
$`(3.33)`$
From (3.31), (3.33) after excluding $`L_1`$ and $`M_2`$ we obtain

$$K_1M_1=(2+\gamma l_3)(l_3+1).$$
$`(3.34)`$
As $`K_1`$ and $`M_1`$ are positive, and $`l_31,`$ then

$$l_3>\frac{2}{\gamma }.$$
$`(3.35)`$
At $`\gamma >2`$ it goes already to the contradiction.

Further, from (3.32) and (3.33) we have $`\frac{L_1^2}{K_1M_1}=\frac{2+\gamma l_3}{l_3+1},`$ this fact together with (3.34) give

$$L_1^2=(2+\gamma l_3)^2.$$
$`(3.36)`$
From (3.31), (3.34), (3.35) and (3.36) we obtain the inequality

$$\frac{2}{\gamma }<1\sqrt{1(2+\gamma l_3)^22(2+\gamma l_3)(l_3+1)},$$
$`(3.37)`$
which cannot be true at $`\gamma >1.`$

If $`l_1<1,l_2=1,l_3<1,`$ then from (3.23) we get $`p_2=\gamma L_4+\frac{2L_2N}{M_2},`$

$$\frac{2L_2N}{M_2}=2\gamma l_3,$$
$`(3.38)`$
from (3.25) $`q=\gamma L_4+\frac{L_2M_1}{2N},`$

$$\frac{L_2M_1}{2N}=2\gamma l_3,$$
$`(3.39)`$
from (3.27)

$$p_2=\gamma L_4=\gamma l_3.$$

In that way, $`l_3=\frac{2}{\gamma },`$ this fact together with (3.38) or (3.39) contradicts to the fact that $`L_2,N`$ and $`M_1`$ are not equal to zero.

If $`l_1<1,l_2=1,l_3=1,`$ then from (3.23),(3.25) and (3.27) we get correspondingly

$$\frac{2L_2N}{M_2}+\frac{L_1M_1}{M_2}=2\gamma l_3,$$
$`(3.40)`$
$$\frac{L_2M_1}{2N}=2\gamma l_3,$$
$`(3.41)`$
$$\frac{L_1M_2}{M_1}=2\gamma l_3,$$
$`(3.42)`$
from (3.20), (3.21)

$$\frac{K_1M_2}{L_1}=\frac{2K_1N}{L_2}=l_3+1.$$
$`(3.43)`$
In that way, multiplying (3.40) by (3.42), taking into account (3.43) we get

$$(2+\gamma l_3)^2=L_1^2+\frac{M_2L_2^2}{M_1},$$

from (3.41), (3.42), (3.43)

$$M_1^2=M_2^2,$$
$`(3.42)`$
that is

$$(2+\gamma l_3)^2=L_1^2\pm L_2^2.$$
$`(3.45)`$
As above, from (3.31), (3.34), (3.35) and (3.45) we get the inequality (3.37) which cannot hold at $`\gamma >1.`$

It remains the unique possibility $`l_1=l_2=l_3=1.`$ In the case $`p_2<2,q<2`$ and, as follows from (3.23), $`p_1=p_2`$ or $`p_1=q,`$ that is $`p_1<2.`$ Therefore, as follows from (3.29)

$$L_1^2+L_2^2L_3^2=1.$$
$`(3.46)`$
If $`p_1=p_2`$, then $`|M_2|<M_1.`$ If $`p_1=p_2,q<p_1,`$ then from (3.23), (3.27) we obtain that $`p_2+\gamma =\frac{L_1M_1}{M_2},p_1+\gamma =\frac{L_1M_2}{M_1},`$ hence $`M_1^2=M_2^2,`$ it goes to the contradiction.

If $`p_1=q,p_2<p_1,`$ then from (3.23), (3.25), (3.27) we get

$$L_1M_1=2L_2N,\lambda =\frac{L_2M_1}{2N}=\frac{2L_3N}{M_1},$$
$`(3.47)`$
where $`\lambda =q+\gamma =p_1+\gamma .`$ It follows, in particular, that $`L_2^2=\lambda L_1,L_1L_3=\lambda L_2,\lambda ^2=L_2L_3.`$ Besides, if $`G_1(t)N_1t^{q_1},G_2(t)N_2t^{q_2},t\mathrm{},`$ where $`N_1,N_2`$ are some positive constants, $`q_1=q_2=p_1=q,N_1=N_2,M_1=2N_1,`$ and $`N_1^2N^2`$ in virtue of $`G_1(t)G_2(t)G_3^2(t)>0.`$ In such way, from (3.47) we get that $`L_2^2\lambda ^2,\lambda ^2L_3^2,L_1^2L_2^2,`$ and taking into account (3.46), $`L_1^21,L_2^21,\lambda ^21.`$ That is if $`\lambda >1(\gamma >3),`$ the conditions mentioned in the paragraph cannot hold together.

At last, consider the case $`p_1=p_2=q.`$ Then from (3.23), (3.25), (3.27) we have $`\lambda =\frac{L_1M_1}{M_2}+\frac{2L_2N}{M_2}=\frac{L_2M_1}{2N}\frac{L_3M_2}{2N}=\frac{L_1M_2}{M_1}\frac{2L_3N}{M_1},`$ that is the system of linear homogeneous equations with respect to the variables $`M_1,M_2,N`$

$$L_1M_1\lambda M_2+2L_2N=0,$$

$$L_2M_1L_3M_22\lambda N=0,$$

$$\lambda M_1L_1M_2+2L_3N=0.$$

For the existence of its nontrivial solution the determinant of the system must be equal to zero, i.e.

$$\lambda L_1^22\lambda L_2L_3+L_1L_2^2+L_1L_3^2\lambda ^3=0.$$

Taking into account (3.46), involve the function with respect to the variables $`L_2`$ and $`L_3`$, where $`\lambda `$ plays the role of parameter, namely

$$\mathrm{\Psi }_\lambda (L_2,L_3)=\lambda (1L_2^2+L_3^2)2\lambda L_2L_3+(L_2^2+L_3^2)\sqrt{1L_2^2+L_3^2}\lambda ^3.$$

By the standard methods one can show that the function is not equal to zero at $`\lambda >1,`$ i.e. at $`\gamma 3`$ (remember, that $`p_1<2`$).

So, it remains to investigate the case $`1<\gamma 3.`$

Consider equation (3.11), which can be written as

$$\mathrm{\Delta }^{}(t)=2d_1(t)\mathrm{\Delta }(t).$$

If we suppose that $`\mathrm{\Delta }(t)constt^m,t\mathrm{},m`$ is a constant, then $`m=2L_4=2.`$ From (3.11) we can get $`E_p(t)constt^{1\gamma }(0),t\mathrm{}.`$ Therefore, if we denote $`E`$ the quantity of the total energy of the system, then $`E_k(t)=EE_p(t)E(1C_1t^{1\gamma }),`$ here and further $`C_i`$ are some positive constants. But $`E_k(t)C_2G_i(t)t^2\mathrm{\Delta }^{\frac{\gamma +1}{2}}C_3G_kt^{\gamma 1},`$ where $`G_k`$ is at least one of functions $`G_1,G_2`$ or $`G_3.`$ That is $`G_kC_4(1C_1t^{1\gamma })t^{1\gamma }C_5t^{1\gamma }.`$ But then $`1\gamma p_1<2,\gamma >3,`$ and $`\gamma `$ do not belong to the interval under consideration.

So the proof of the proposition is over.

Remark 3.1. One can show more shortly, that in the physical case $`1<\gamma 2`$ the situation $`l_1=l_2=l_3=1`$ if impossible. Taking into account (3.10) we have $`G_1G_2G_3^2=\mathrm{\Delta }^\gamma constt^{2\gamma },t\mathrm{}.`$ But the degree of the leading term of the expression $`G_1G_2G_3^2`$ is not greater then $`2p_1,`$ therefore $`p_1\gamma ,`$ as $`p_1<2,`$ then $`\gamma >2.`$

Remark 3.2. Actually $`l_3=L_4=2,G_1G_2G_3^2constt^{4\gamma },t\mathrm{}.`$

Remark 3.3. In the case $`\mu >0,l=0`$ the asymptotic $`A(t)\frac{\delta }{t}E,t\mathrm{}`$ holds as well, that is according to Corollary 2.1 the corresponding solution with a linear profile of velocity is interior. If $`\mu >0,l0`$ it occurs that $`A(t)A_0\frac{1}{t},`$ but the matrix $`A_0`$ under arbitrary initial conditions $`A(0)`$ has not to be of form $`\delta E+U_2`$ (see the denotations in the statement of Corollary 2.2), i.e. the Corollary cannot be applied in the last case.

One can write out the conditions of compatibility of initial data of the density and the pressure for the case of solution with arbitrary linear profile of velocity.


format: text

text:
# Phase separation in the combined Falicov-Kimball and static Holstein model

## Abstract

The Falicov-Kimball model has been known to phase separate when the correlation strength is large and the temperature is low. We examine the stability of phase separation under the influence of an additional electron-phonon interaction that tries to force the electron density to be homogeneous. We find that the electron-phonon interaction eventually destroys the phase-separated state once it becomes larger than some critical value. The results of this work may have an influence on theories of collossal magentoresistance materials.

The phenomenon of strong electron correlations driving phase separation has been long conjectured in the Falicov-Kimball modelFreericks and Falicov (2001) and in the Hubbard modelEmery et al. (1990). Little progress has been made on rigorously establishing phase separation in the Hubbard model, but it has been proved to exist in the Falicov-Kimball model in one dimensionLemberger (1992), in infinite dimensionsFreericks et al. (1999), and in all other dimensionsFreericks et al. (2002a, b). In spite of these rigorous results, there has still been only limited calculations of phase diagrams and the transition temperature as a function of the parameters in the model. Essentially all finite-temperature calculations have been performed in the infinite-dimensional limitBrandt and Mielsch (1991); Letfulov (1999); Freericks et al. (1999); Freericks and Lemański (2000). Here we continue with this tradition, and illustrate how competing (finite) interactions can destroy the phase-separated state, even if the Falicov-Kimball interaction strength is infinitely large.

The system we choose to study is the harmonic version of the spinless Falicov-Kimball–static-Holstein modelFalicov and Kimball (1969); Holstein (1959). The Hamiltonian is

$`H`$ $`=`$ $`{\displaystyle \frac{t^{}}{2\sqrt{d}}}{\displaystyle \underset{i,j}{}}c_i^{}c_j+(E_f\mu ){\displaystyle \underset{i}{}}f_i^{}f_i+U{\displaystyle \underset{i}{}}c_i^{}c_if_i^{}f_i`$ (1)
$`+`$ $`{\displaystyle \frac{1}{2}}\kappa {\displaystyle \underset{i}{}}x_i^2+{\displaystyle \underset{i}{}}(gx_i\mu )(c_i^{}c_i\rho _e),`$ (2)
where $`c_i^{}`$ ($`c_i`$) is the creation (annihilation) operator for an itinerant electron at site i, $`f_i^{}`$ ($`f_i`$) is the creation (annihilation) operator for a localized electron at site i, $`t^{}/2\sqrt{d}`$ is the hopping integral between nearest-neighbor sites (each pair is summed over twice for hermiticity), $`E_f`$ is the local f-electron site energy, $`U`$ is the Falicov-Kimball interaction, $`x_i`$ is the static phonon coordinate, $`\kappa `$ is the spring constant (chosen to be 1 here), and $`g`$ is the electron-phonon interaction. The symbol $`\rho _e`$ denotes the average itinerant-electron concentration $`c^{}c`$. We examine the system in the infinite-dimensional limitMetzner and Vollhardt (1989) and take $`t^{}`$ as our energy unit. A common chemical potential $`\mu `$ is used to fix the total electron concentration (itinerant plus localized). Note that we have chosen the phonons to interact with the local charge fluctuation, which is a common choice, but not the only one that could be made. Since the term that is subtracted depends on the localized electron concentration it is not just a constant that can be included or removed; it has a nontrivial effect on the value of the Helmholz free energy, and thereby on the phase separation. For concreteness we solve our system of equations on an infinite-coordination Bethe lattice, which has a semicircular electronic density of states

$$\rho (ϵ)=\frac{1}{2\pi }\sqrt{4ϵ^2}.$$
(3)
The system is solved by employing dynamical mean-field theory. The procedure simply combines the strategy used for the Falicov-Kimball modelBrandt and Mielsch (1989) with that used for the static Holstein modelMillis et al. (1995, 1996). We summarize the procedure here to standardize our notation. Our interest is in the $`U\mathrm{}`$ limit. In that case, the itinerant electrons and the localized electrons are not allowed to sit at the same lattice site. The dynamical mean-field theory is straightforward to evaluate. The partition function becomes

$`Z`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑x2e^{\frac{\beta }{2}(gx\mu )}{\displaystyle \underset{n}{}}{\displaystyle \frac{G_0^1(i\omega _n)gx}{i\omega _n}}`$ (4)
$`\times `$ $`\mathrm{exp}\{\beta [(gx\mu )\rho _e+{\displaystyle \frac{1}{2}}\kappa x^2]\}`$ (5)
$`+`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑x2e^{\frac{\beta }{2}(gx\mu +U)}{\displaystyle \underset{n}{}}{\displaystyle \frac{G_0^1(i\omega _n)gxU}{i\omega _n}}`$ (6)
$`\times `$ $`\mathrm{exp}\{\beta [(gx\mu )\rho _e+E_f\mu +{\displaystyle \frac{1}{2}}\kappa x^2]\},`$ (7)
with $`\beta =1/T`$ the inverse temperature, $`G_0(i\omega _n)=1/(i\omega _n+\mu \lambda _n)`$ the effective medium, $`i\omega _n=i\pi T(2n+1)`$ the Fermionic Matsubara frequency, and $`\lambda (i\omega _n)=\lambda _n`$ the dynamical mean field. Once the partition function is determined (in the limit $`U\mathrm{}`$), then one can calculate the single-particle electronic Green’s function via

$$G(i\omega _n)=G_n=_{\mathrm{}}^{\mathrm{}}𝑑x\frac{w_0(x)}{G_0^1(i\omega _n)gx},$$
(8)
with $`w_0(x)`$ the phonon distribution function (when no localized-electron occupies the lattice site), which satisfies

$`w_0(x)`$ $`=`$ $`2e^{\frac{\beta }{2}(gx\mu )}{\displaystyle \underset{n}{}}{\displaystyle \frac{G_0^1(i\omega _n)gx}{i\omega _n}}`$ (9)
$`\times `$ $`\mathrm{exp}\{\beta [(gx\mu )\rho _e+{\displaystyle \frac{1}{2}}\kappa x^2]\}/Z;`$ (10)
the integral of $`w_0(x)`$ satisfies a sum rule

$$1\rho _f=_{\mathrm{}}^{\mathrm{}}𝑑xw_0(x).$$
(11)
\[The sum rule is not one, because we have taken the limit $`U\mathrm{}`$. If $`U`$ is finite, we would have a second phonon distribution function $`w_1(x)`$ to consider when calculating the Green’s function, which enters when the site is occupied by a localized-electron. In the $`U=\mathrm{}`$ limit, we can always rewrite $`w_0(x)=P(x)w_1(x)`$ where the integral of $`P(x)`$ is one and can be interpreted as the $`U=\mathrm{}`$ phonon distribution function, and $`w_1(x)`$ integrates to the localized electron concentration. But in this work, we never distinguish these objects, so there is no need to elaborate further on them.\] To complete the loop for the self-consistent solution of the problem, we need to note that Dyson’s equation is

$$G_0^1(i\omega _n)=G_n^1+\mathrm{\Sigma }_n,$$
(12)
and

$$G_n=_2^2𝑑ϵ\rho (ϵ)\frac{1}{i\omega _n+\mu \mathrm{\Sigma }_nϵ},$$
(13)
with $`\mathrm{\Sigma }(i\omega _n)=\mathrm{\Sigma }_n`$ the itinerant-electron self energy. This provides a closed set of equations that is needed to solve the many-body problem.

The computational algorithm that we follow is the iterative procedure introduced by JarrellJarrell (1992): (i) we begin with the self energy set equal to zero; (ii) Eq. (13) is used to determine $`G_n`$; (iii) Eq. (12) is used to extract the effective medium $`G_0`$; (iv) Eqs. (7), (10), and (8), then determine $`G_n`$; and (v) Eq. (12) is employed to extract the new self energy. Then steps (ii-v) are repeated until convergence is reached. Once the equations have converged on the imaginary axis, one can take the function $`w_0(x)`$ and the analytic continuation of Eqs. (13), (12), and (8) to determine the Green’s function on the real axis as well.

In order to study phase separation, we also need to calculate the free energy. This is accomplished by following the procedure of Brandt and MielschBrandt and Mielsch (1991). The Helmholz free energy satisfies

$`F`$ $`=`$ $`T[\mathrm{ln}Z+{\displaystyle }_n\mathrm{ln}G_n`$ (14)
$``$ $`{\displaystyle \underset{n}{}}{\displaystyle _2^2}dϵ\rho (ϵ)\mathrm{ln}{\displaystyle \frac{1}{i\omega _n+\mu \mathrm{\Sigma }_nϵ}}]`$ (15)
$``$ $`(E_f\mu )\rho _f.`$ (16)
It is important to include a sufficient number of Matsubara frequencies in the summation to achieve convergence of the summation; indeed this ultimately limits the accuracy of our results.

It is well known that if both $`\rho _e`$ and $`\rho _f`$ lie between zero and one, then the ground state is phase separatedFreericks et al. (1999). The key point in this work is that we fix $`E_f`$. If the chemical potential lies below $`E_f`$ then the ground state contains only itinerant electrons and is not phase separated. As the chemical potential is increased, we reach a critical value of $`\mu `$ where the average localized electron filling becomes nonzero, and the system goes into a phase-separated state. It is straightforward to determine the value of this critical chemical potential when $`g=0`$. The interacting density of states is just the semicircle in Eq. (3) when there are no localized electrons and it becomes

$$\rho _{int}(ϵ)=\frac{1}{2\pi }\sqrt{4(1\rho _f)ϵ^2}$$
(17)
when there are $`\rho _f`$ interacting electrons. So we must compare the energy with no localized electrons $`E_{gs}=(4\mu ^2)^{3/2}/6\pi `$ to the mixture of states with all localized particles and no itinerant particles (with weight $`\rho _f`$) and the state with no localized particles but with an electron filling of $`\rho _e=\rho _{tot}(1\rho _{tot})\rho _f`$ (with weight $`1\rho _f`$) in the limit where $`\rho _f`$ approaches zero. Note that one must adjust the chemical potential $`\mu ^{}`$ of the mixed state, to have the correct electron density. After performing some straightforward algebra, we find a transcendental equation for the critical value of the chemical potential

$`E_f`$ $`=`$ $`{\displaystyle \frac{(4\mu ^2)^{\frac{3}{2}}}{6\pi }}\mu ({\displaystyle \frac{1}{2}}+{\displaystyle \frac{1}{\pi }}\mathrm{sin}^1{\displaystyle \frac{\mu }{2}}`$ (18)
$`+`$ $`{\displaystyle \frac{\mu }{4\pi }}\sqrt{4\mu ^2}).`$ (19)
As $`E_f`$ ranges from $`2`$ to 0, the critical total electron density for phase separation ranges from 0 to 1; for $`E_f>0`$ the system does not phase separate as $`T0`$.

We will examine the case with $`E_f=0.5`$ here. Eq. (19) is solved by $`\mu =0.145`$ and $`\rho _e+\rho _f=0.454`$ for the lower boundary where phase separation begins. Hence phase separation enters when $`\rho _{tot}>0.454`$. To determine the phase separation, we must perform a Maxwell construction. If we plot the electron filling as a function of the chemical potential, we will find a multivalued function in the region of phase separation. One way to perform the Maxwell construction is by drawing a line that has equal areas above and below, with the points of intersection corresponding to the densities of the phases that compose the mixture of the phase-separated state (the other way is to construct the convex hull of the free-energy curve; both methods yield the same results). The double-valued function and the Maxwell construction are plotted for $`g=1.2`$ and $`T=0.04`$ in Fig. 1. It is nontrivial to generate such a curve. The technical difficulty arises from the fact that for a fixed value of $`E_f`$, there are multiple solutions of $`\rho _f`$ which have the same $`\rho _{tot}`$ when one is in the phase separated region. Hence, we sometimes need to either adjust $`\rho _f`$ and calculate $`E_f`$ once the equations have converged, or employ a two-dimensional rootfinder for $`\mu `$ and $`\rho _f`$ to get the equations to converge. Even so, we find regions of parameter space where we are unable to stabilize the algorithm.

We examine the transition temperature for different electron-phonon coupling strengths in Fig. 2. The phase-separated region hardly changes at zero temperature as the electron-phonon coupling increases, but the transition temperatures are sharply reduced. We run into numerical problems with the stability of our iterative solutions when $`g`$ is increased to 1.2 and above. We are unable to find any strategy to solve for the phase separation at low temperature in this region of parameter space. We expect, however, that the phase separation continues to be suppressed, and conjecture that there is a critical value of the electron-phonon coupling where the phase separation disappears, but we cannot determine it’s value with our numerical algorithms. In addition to the binodal (first-order) phase transition temperature, we also include the spinodal decomposition temperature (below which the system cannot be supercooled) for the $`g=0`$ case. This is found by determining where the $`k=0`$ charge susceptibility diverges. Of course the binodal and spinodal temperatures must agree at the peak of the curve.

The disappearance of the phase separation occurs when the electron-phonon interaction energy becomes similar in magnitude to the kinetic-energy gain of the phase separated state. This occurs at a finite value of $`g`$ because the kinetic energy gain is finite, and does not diverge as $`U\mathrm{}`$.

To understand this phenomenon further, we examine in detail the system in the low-temperature limit. We already saw that when $`g=0`$ there is a minimum total electron density that must be reached before phase separation occurs. Similarly, when $`g0`$, we find that such a critical density continues to hold for small $`g`$. This occurs because, as $`T`$ is lowered, the function $`w_0(x)`$ becomes sharply peaked at $`x=0`$, eventually becoming a delta function at $`T=0`$. Hence, the system remains a Fermi gas in this regime of parameter space. As $`g`$ is increased to a critical value, $`w_0(x)`$ develops a local minimum at $`x=0`$ and two maxima, one at $`x<0`$ and one at $`x>0`$. At this point, the ground state becomes a non-Fermi-liquid, and we believe that the non-Fermi-liquid character suppresses the phase separation, ultimately eliminating it. To be more quantitative, we compute

$$T\frac{d\mathrm{ln}w_0(x)}{dx}=g(\rho _e\frac{1}{2})\kappa xgT\underset{n}{}\frac{1}{G_0^1(i\omega _n)gx}.$$
(20)
When this derivative vanishes, we have an extremum of the phonon distribution function. When $`\rho _f=0`$, this equation always has a solution at $`x=0`$ because the sum over Matsubara frequencies equals $`\rho _e1/2`$. If $`g`$ is increased, we reach a point where the second derivative becomes positive, and the $`x=0`$ solution becomes a local minimum instead of a maximum. At this point two nonzero $`x`$ solutions exist, and the system becomes a non-Fermi-liquid. The second derivative (at $`x=0`$) is

$$T\frac{d^2\mathrm{ln}w_0(x=0)}{dx^2}=\kappa g^2T\underset{n}{}G_0^2(i\omega _n),$$
(21)
which becomes

$`T{\displaystyle \frac{d^2\mathrm{ln}w_0(x=0)}{dx^2}}`$ $`=`$ $`\kappa g^2{\displaystyle _2^\mu }𝑑ϵ\rho (ϵ)\mathrm{Re}G(ϵ)`$ (22)
$`=`$ $`\kappa +g^2{\displaystyle \frac{(4\mu ^2)^{\frac{3}{2}}}{12\pi }}`$ (23)
in the low-temperature limit. For $`\mu =0.145`$ the critical value of $`g`$ to enter the non-Fermi-liquid state is 2.18.

It is useful to discuss the impact of this work on models for the collossal magnetoresistance manganite materials. There are a number of different elements that appear to be important in describing the CMR materials. All models start with double-exchangeFurukawa (1994), but then need something else to enhance the metal-insulator transition. Three main ideas have emerged to do this: a strong electron-phonon interactionMillis et al. (1995, 1996); Edwards (2002); strong disorder scatteringAllub and Alascio (1996); Varma (1996); Allub and Alascio (1997); and proximity to phase separationDagotto et al. (2001). The results of this investigation indicate that the three elements are not fully independent. Indeed, strong disorder scattering can lead to phase separation, and strong electron-phonon scattering can suppress the phase separation. One of the consequences of this work is that one needs to carefully investigate the implications of all three of these physical effects, and how they interrelate, to construct a meaningful theory for the manganites.

In conclusion, we have analyzed the suppression of the phase separation transition in a combined Falicov-Kimball–static-Holstein model. We see that increasing the electron-phonon coupling forces the system to be more homogeneous and suppresses the phase separation. We conjecture that a finite value of $`g`$ will be sufficient to remove phase separation entirely, but are unable to access that region of the phase diagram due to numerical difficulties. Our results could have an impact on CMR materials that have both disorder and electron-phonon interactions. Any phase separation induced by the disorder may be suppressed by the electron-phonon coupling.

Acknowledgments: B.M.L. acknowledges support from the Project for Supporting of Scientific Schools, N 00-15-96544. J.K.F. acknowledges support from the National Science Foundation under Grant No. DMR-9973225.


text:
# The most distant structure of galaxies known: a protocluster at 𝑧=4.1Based on observations carried out at the European Southern Observatory, Paranal, Chile; program LP167.A-0409

## 1 Introduction

Studies of high redshift (proto)clusters of galaxies ($`z>2`$) can directly constrain theories of galaxy evolution and cosmological models (e.g. Bahcall & Fan 1998), but the detection of (proto)clusters at these redshifts using conventional optical and X-ray techniques is difficult. High redshift radio galaxies (HzRGs, $`z>2`$) can help: they are among the most evolved and most massive galaxies in the early Universe (De Breuck et al. 2002) and are most likely located in dense (proto)cluster environments (e.g., Le Fèvre et al. 1996; Pascarelle et al. 1996; Röttgering et al. 1996; Carilli et al. 1997; van Ojik et al. 1997; Pentericci et al. 2000b). In addition, HzRGs have properties that would be expected of forming central cluster galaxies. Their extremely clumpy morphologies revealed by Hubble Space Telescope images (Pentericci et al. 1999) are strikingly similar to simulations of forming brightest cluster galaxies, based on hierarchical models (e.g., Aragón-Salamanca, Baugh, & Kauffmann 1998).

A pilot project on the Very Large Telescope (VLT) to search for Ly$`\alpha `$-emitting galaxies around the clumpy radio galaxy PKS 1138–262, resulted in the discovery of 14 galaxies and a QSO at approximately the same redshift as the radio galaxy (Kurk et al. 2000; Pentericci et al. 2000a). If the structure found is virialized, the total mass of the protocluster would be $`10^{14}`$ $`M_{\mathrm{}}`$. Motivated by this result, we started a large program at the VLT to search for forming clusters (protoclusters) around HzRGs at redshifts 2 and higher. The first radio galaxy field we observed was TN J1338–1942 at a redshift of 4.1 (De Breuck et al. 1999). This HzRG is one of the brightest and most luminous in Ly$`\alpha `$ known. Both its Ly$`\alpha `$ profile and radio structure are very asymmetric (De Breuck et al. 1999), which indicates strong interaction with dense gas, and the rest-frame radio luminosity is comparable to that of the most luminous 3CR sources (P$`{}_{178\mathrm{M}\mathrm{H}\mathrm{z}}{}^{}4\times 10^{35}`$ erg s-1 Hz-1 sr-1). Here we report on the discovery of a substantial overdensity of Ly$`\alpha `$ emitters around this radio galaxy at $`z4.1`$.111Throughout this Letter, magnitudes are in the AB system and a $`\mathrm{\Lambda }`$-dominated cosmology with H$`{}_{0}{}^{}=65`$ km s-1 Mpc-1, $`\mathrm{\Omega }_M=0.3`$, and $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ is assumed.

## 2 Observations and candidate selection

### 2.1 VLT imaging and selection of candidate Ly$`\alpha `$ emitters

We carried out narrow- and broadband imaging on 2001 March 25 and 26 with the 8.2 m VLT Kueyen (UT2), using the imaging mode of the FOcal Reducer/low dispersion Spectrograph 2 (FORS2). At $`z=4.10`$ the Ly$`\alpha `$ emission line is redshifted to 6202 Å, which falls in our custom narrowband filter with a central wavelength of 6195 Å and FWHM of 60 Å. The broadband R filter had a central wavelength of 6550 Å and FWHM of 1650 Å. The detector was a SiTE CCD with 2048$`\times `$2048 pixels, with a scale of 0$`\stackrel{}{\mathrm{.}}`$2 per pixel and a field of view of 6$`\stackrel{}{\mathrm{.}}`$8$`\times `$6$`\stackrel{}{\mathrm{.}}`$8. We took 18 separate 1800 s exposures and one 900 s exposure in the narrow band and 21 exposures of 300 s in $`R`$, shifted by $`15`$″ with respect to each other to minimize flat-fielding problems and to handle cosmic rays. The nights were photometric and the average seeing in both narrowband and broadband images was 0$`\stackrel{}{\mathrm{.}}`$65. The 1 $`\sigma `$ limiting magnitude per square arcsecond was 28.6 for the narrowband and 29.2 for the broadband image. For the flux calibration the spectrophotometric standard star GD 108 (Oke 1990) was used. The final images have sizes of 6$`\stackrel{}{\mathrm{.}}`$4 $`\times `$ 6$`\stackrel{}{\mathrm{.}}`$2 (39.7 arcmin2). The total volume probed at $`z=4.1`$ by our narrowband filter is 7315 Mpc3.

NTT imaging data, taken on 1998 April 27 and 29 under nonphotometric conditions with SuSI2, were also used to provide $`B`$\- and $`I`$-band magnitudes for candidate emitters where possible. The 1 $`\sigma `$ limiting magnitude per square arcsecond was about 27.8 in both bands.

For the detection and photometry of objects in the images, we used the object detection and classification program SExtractor (Bertin & Arnouts 1996). Detected objects had at least 9 connected pixels with values larger than the rms sky noise on the narrowband image. A photometric analysis was then carried out on both the narrowband and the broadband image. In total 2407 objects were extracted. Based on the statistics of our detection of Ly$`\alpha `$ emitters around PKS 1138–262, we selected objects with a rest-frame equivalent width (EW$`{}_{0}{}^{})`$ greater than 15 Å (or $`m_{\mathrm{BB}}`$ $`m_{\mathrm{NB}}>0.84`$) and significance $`\mathrm{\Sigma }>3`$ as candidate Ly$`\alpha `$-emitting galaxies, with $`\mathrm{\Sigma }`$ the ratio of continuum-subtracted counts in the narrowband to the combined noise in the broadband and narrowband (Bunker et al. 1995). Of the 2407 objects, 34 objects satisfied these criteria, including the radio galaxy. Of these 34 objects, 31 were detected in the $`R`$-band image, 9 were detected in the NTT $`I`$ image and 5 in the NTT $`B`$ image. The 5 objects with a detection in B have a color of $`BR`$ $``$ 0–1, inconsistent with the colors expected for galaxies at $`z=4`$ (Steidel et al. 1999). Three of these objects are located in the halo of the radio galaxy, which affects the narrowband photometry. The other two are likely to be foreground objects with another line falling in the narrowband filter, e.g. \[O III\] $`\lambda 5007`$ or \[O II\] $`\lambda 3727`$. Excluding the radio galaxy, the resulting 28 objects were our Ly$`\alpha `$-emitting candidates for follow-up spectroscopy.

### 2.2 VLT spectroscopy

On 2001 May 20, 21, and 22 we carried out spectroscopy using FORS2 in the mask multiobject spectroscopy mode with standard resolution. The nights were photometric with an average seeing of 1″. The spectra were obtained with the 600RI grism with a dispersion of 1.32 Å pixel-1 and a wavelength range from 5300 Å to 8000 Å. This grism was preferred because of its high throughput (peak efficiency is 87%). Two different masks were used to observe 23 of the candidate Ly$`\alpha `$ emitters and the radio galaxy, with slit sizes of 1″, resulting in a resolution of 6 Å, which corresponds to 290 km s-1 at $`z=4.1`$. The total exposure time was 31,500 s for the first mask and 35,100 s for the second mask. Wavelength calibration was obtained from exposures of He, Ar, Ne, and HgCd arc-lamps. The accuracy of the wavelength calibration was better than 0.05 Å. For the flux calibration, long slit exposures with a 5″ slit were used of the spectrophotometric standard stars EG 274 and LTT 7987 (Stone & Baldwin 1983; Baldwin & Stone 1984).

## 3 Results

### 3.1 Line emitters in the field

#### 3.1.1 Line identifications

Of the 23 candidates observed, 20 show an emission line with a peak between $`\lambda =6187`$ Å and $`\lambda =6216`$ Å. The signal-to-noise of these lines is at least 10. Ten randomly selected spectra from our sample are shown in Figure 1. Two nondetections were very faint ($`m_{\mathrm{NB}}=25.9`$ and $`m_{\mathrm{NB}}=26.1`$ respectively). The third nondetection had a very low surface brightness. The success rate of our selection criteria was therefore 87%.

The first question is whether the detected lines are indeed due to Ly$`\alpha `$ at the redshift of the radio galaxy. Intervening neutral hydrogen will absorb emission blueward of the Lya line (Steidel et al. 1999). This discontinuity of the continuum over the Lya emission line is observed in one of the spectra, but the continuum emission of the other emitters was too faint to be detected (R $`27`$). Identification of the lines with \[O III\] $`\lambda 5007`$ at $`z0.24`$ can be excluded, because of the lack of confirming lines \[O III\] $`\lambda 4959`$ and H$`\beta `$ in the spectra. The position of the emitter within the slit mask determines the wavelength coverage of the resulting spectrum. In 9 cases this coverage was suitable to exclude the identification with \[O II\] $`\lambda 3727`$ at $`z0.66`$ on similar grounds. The detected emitters are anyway unlikely to be foreground \[O II\] galaxies. First of all, if one of the lines would be \[O II\], then the rest-frame equivalent width would be at least 70 Å. A survey of nearby field galaxies, conducted by Jansen et al. (2000), gives a mean EW0 of the \[O II\] line of $`30`$ Å for galaxies with $`M_B=16`$, which roughly compares to R $`27`$ at $`z=0.66`$. Only 2 galaxies out of 159 galaxies with \[O II\] in emission have EW0(\[O II\]) $`>70`$ Å. The number of \[O II\] emitters expected in our field, using another study (Hogg et al. 1998), is $`7`$. Again, only a few percent of the \[O II\] emitters observed by Hogg et al. have an EW$`{}_{0}{}^{}>70`$ Å. Therefore, from our sample of 20 emitters, $`<1`$ is expected to be an \[O II\] emitter.

Additional evidence that the observed lines are predominantly Ly$`\alpha `$ lines at $`z4.1`$ associated with the radio galaxy is provided by their velocity distribution, which has a dispersion of $`326\pm 73`$ km s-1 and a FWHM of $`768\pm 172`$ km s-1 (Fig. 2), which is 4 times smaller than the FWHM of the narrowband filter. Further, the velocity distribution peaks within 200 km s-1 of the radio galaxy Ly$`\alpha `$ peak, corrected for absorption (De Breuck et al. 1999).

For all these reasons we interpret the observed emission lines as Ly$`\alpha `$. To estimate the redshifts, flux densities, and widths (FWHM) of the lines, a Gaussian function was fitted to each of the one-dimensional spectra. Details will be provided in a future paper (B.P. Venemans et al., in preparation).

#### 3.1.2 Significance and properties of the overdensity

The next question to be addressed is to what extent the statistics of our detections represent a significant overdensity of galaxies in the field. A “blank-field” study of Ly$`\alpha `$ emitters at approximately the same redshift is the Large-Area Lyman Alpha survey (LALA survey, Rhoads et al. 2000). Preliminary results of this survey indicate a number density of $`4000\pm 460`$ deg-2 $`\mathrm{\Delta }z^1`$ for objects with EW$`{}_{0}{}^{}15`$Å and line + continuum $`>2.6\times 10^{17}`$ erg s-1 cm-2. For our field the second criterion corresponds to $`m_{\mathrm{NB}}<24.55`$ and we expect $`2.3\pm 0.3`$ such Ly$`\alpha `$ emitters within our probed volume. However, 9 of the confirmed emitters satisfy the above criteria (including the radio galaxy). In our cosmology the comoving volume density of the LALA survey is $`n_{LALA}=(3.1\pm 0.4)\times 10^4`$ Mpc-3. The volume density in our field is $`n_{1338}=9/7315`$ Mpc$`{}_{}{}^{3}=(12\pm 4)\times 10^4`$ Mpc-3. The difference in number density is $`n_{1338}`$/$`n_{LALA}=4.0\pm 1.4`$. However, the FWHM of the velocity distribution is approximately 4 times smaller than the FWHM of the filter (see Fig. 2). This implies that our radio galaxy field is overdense in Ly$`\alpha `$ emitters by a factor of 15 compared with a blank field.

The spatial distribution of the emitters is not homogeneous (Fig. 3). The structure appears to be bound in the northwest but our FOV is not large enough to show a boundary in the south. The size of the structure is therefore at least 6′$`\times `$ 4′, which corresponds to greater than 2.7 Mpc $`\times `$ 1.8 Mpc. Remarkably, the radio galaxy, probably the most massive system in the structure, does not appear to be at the center of the protocluster, in disagreement with models for the formation of dominant cluster galaxies (e.g., West 1994).

### 3.2 Radio galaxy halo

The radio galaxy is located close to the apparent north-west boundary of the galaxy overdensity structure. The radio emission is dominated by two components, separated by 5$`\stackrel{}{\mathrm{.}}`$5 (De Breuck et al. 1999). The brightest component coincides with the optical emission, while the other is in the south-east, towards the center of the galaxy overdensity. A spectacular feature of the radio galaxy, visible in the narrowband image, is the large Ly$`\alpha `$ halo (Fig. 4). Although giant Ly$`\alpha `$ halos are a common feature of HzRGs, they are usually fairly symmetrically extended around the radio galaxy. In the case of TN J1338–1942, the halo is highly asymmetric and extends for $`15`$″ (110 kpc) to the north-west, i.e., away from the center of the overdensity structure. In its narrowness and asymmetry the TN J1338–1942 halo is unique among known distant radio galaxies. Possible mechanisms to produce this structure include cooling flows in colliding sub-structures and buoyancy effects (Gisler 1976) and will be considered in a subsequent paper.

## 4 Discussion

### 4.1 Nature of the overdensity

Could the structure that we have detected be a protocluster at $`z4.1`$ that will evolve into a rich cluster of galaxies in the local Universe? At $`z=4.1`$ the Universe is $`1.6`$ Gyr old, too short for the structure to have virialized since the mean crossing time for galaxies at the observed velocity dispersion is at least 4 Gyr. Thus the observed velocities are probably infall velocities of the galaxies accreting onto a large overdensity. The total mass of this structure can be estimated by using $`M=\overline{\rho }V(1+\delta _m)`$ with $`\overline{\rho }`$ the mean density of the Universe and $`\delta _m`$ the mass overdensity within our volume $`V`$ (Steidel et al. 1998). The mass overdensity is related to the galaxy overdensity $`\delta _{\mathrm{gal}}`$ through $`1+b\delta _m=C(1+\delta _{\mathrm{gal}})`$, where $`C`$ takes into account the redshift space distortions caused by peculiar velocities and $`b`$ is the bias parameter. From the statistics of redshift “spikes”, Steidel et al. (1998) argue that $`b4`$. Taking $`b`$ in the range 3–5 and assuming that the structure is just breaking away from the Hubble expansion, $`\delta _m`$ is estimated to be 0.5–2.3, giving a mass of our structure of (1–2) $`\times 10^{15}`$ $`M_{\mathrm{}}`$. This is comparable to the mass of the Lyman break galaxy overdensity (spike) found by Steidel et al. and to that of the Coma cluster (e.g., Fusco-Femiano & Hughes 1994).

### 4.2 Relation to overdensity spikes

It is instructive to compare the number of luminous radio sources at $`z3`$ to the number of redshift spikes. Steidel et al. (1998) estimate that about one velocity spike is detected in each 9′$`\times `$ 18′ field in their spectroscopic survey, corresponding to $`9\times 10^5`$ redshift spikes in the whole sky. Their survey is sensitive to redshifts between 2.7 and 3.4, corresponding to a cosmic evolution time of 0.6 Gyr.

How many luminous radio sources are there at $`z3`$? Using the pure luminosity evolution model of Dunlop & Peacock (1990) to describe the steep spectrum radio luminosity function, we estimate that in the redshift range $`2.7<z<3.4`$ there should be $`1.2\times 10^4`$ radio sources with luminosities exceeding $`10^{33}`$ erg s-1 Hz-1 sr-1 at 2.7 GHz (“Cygnus A type” radio sources: P2.7GHz(Cygnus A) $`2\times 10^{33}`$ erg s-1 Hz-1 sr-1, Becker, White, & Edwards 1991). Assuming that HzRGs are only once active for $`10^7`$ yr (Blundell & Rawlings 1999), we expect the number of (previously) active radio sources in this redshift range to be $`7\times 10^5`$.

Hence, the luminosity functions and lifetimes of luminous radio sources are consistent with every velocity “spike” in the space densities of Lyman break galaxies being associated with a massive galaxy that has been or will become a luminous radio source once. Note that West (1994) presented similar statistical evidence to argue that distant powerful radio galaxies are the precursors of cD galaxies at the centers of galaxy clusters.

## 5 Conclusion

We have found a structure of 20 Ly$`\alpha `$-emitting galaxies around the high redshift radio galaxy TN J1338–1942. The overdensity of this protocluster is on the order of 15 compared to field samples. Our results demonstrate that by $`z=4.1`$, megaparsec-scale structure had already formed.

Together with our previous data, this implies that the most luminous radio sources are tracers of regions of galaxy overdensity in the early Universe. The estimated masses of $`10^{14}`$$`10^{15}`$ $`M_{\mathrm{}}`$ are consistent with the overdensities being ancestors of rich clusters of galaxies in the local Universe.

We thank the staff on Paranal, Chile and Gero Rupprecht at ESO for their splendid support, and William Grenier of Andover Corporation for his efforts in ensuring the narrowband filter was manufactured in time for the observations. We acknowledge useful discussions with P. Martin (University of Toronto) and M. Davis (University of California at Berkeley) on radio galaxy halos and protoclusters, and their possible role in constraining cosmological models. We thank the referee, Adam Stanford, for improving this Letter. The work by W. v. B. was performed under the auspices of the US Department of Energy, National Nuclear Security Administration by the University of California, Lawrence Livermore National Laboratory under contract No. W-7405-Eng-48. The NRAO is operated by associated universities inc., under cooperative agreement with the NSF. This work was supported by the European Community Research and Training Network “The Physics of the Intergalactic Medium”.