text:
# Spectral Gaps for Periodic Elliptic Operators with High Contrast: an Overview

## Introduction

We consider periodic elliptic partial differential operators of the second order with high contrast in (some of) the coefficients. The classes of operators discussed here include Schrödinger operators, magnetic Schrödinger and Pauli operators, the acoustic operator, Laplace-Beltrami operators, and Maxwell operators. We are mainly interested in the question whether high contrast in the coefficients leads to weak coupling between the period cells and to spectral concentration and/or spectral gaps. For most of the paper, we assume a simple geometry like the one shown in Figure 1, for which the strongest results are obtained. We include some remarks on more general situations, like sponges in dimensions $`3`$ and higher. Our main theme is to show that simple topological assumptions lead to spectral consequences for a wide class of—rather different—elliptic operators. We do not discuss the question of absolute continuity here.

The paper is organized as follows. Section 1 contains some basic definitions, assumptions and background results concerning monotone convergence of quadratic forms. In Section 2, we briefly report on some “classical” results on the convergence of Schrödinger operators with high barriers. Similar results for strong magnetic barriers in the Schrödinger and Pauli case are discussed in Section 3. Section 4 is devoted to divergence type operators $`(1+\lambda \chi _\mathrm{\Omega })`$, where $`\mathrm{\Omega }`$ is a periodic, open subset of $`^m`$, with $`m2`$. There is some superficial similarity between divergence type operators and Laplace-Beltrami operators. The existence of spectral gaps for periodic Laplace-Beltrami operators with a conformal factor in the metric tensor is discussed in Section 5. We also refer to constructions of a purely “geometric” nature that yield Laplace-Beltrami-operators with gaps. We conclude in Section 6 with some references to the work of A. Figotin, P. Kuchment and others on the Maxwell operator and applications to photonic crystals.

## 1. Preliminaries and basic assumptions

Throughout the present paper, we will always assume that $`\mathrm{\Omega }^m`$ is open and periodic w.r.t. some discrete lattice $`\mathrm{\Gamma }`$; for most purposes it is enough to consider $`\mathrm{\Gamma }=^m`$. $`Q`$ denotes a fundamental period cell of $`\mathrm{\Gamma }`$; we may take $`Q=(0,1]^m`$ if $`\mathrm{\Gamma }=^m`$. The complement of $`\mathrm{\Omega }`$ is denoted as $`M`$, a closed and periodic set. We will study two-component media with material parameter $`\lambda =1`$ on $`M`$, and $`\lambda 1`$ on $`\mathrm{\Omega }`$. The strongest results are obtained in the “model situation” of Figure 1 where $`\mathrm{\Omega }`$ contains the boundary of the period cell $`Q`$.

Our basic Hilbert spaces are $`=L_2(^m)`$ and $`^1(^m)`$, the Sobolev space of $`L_2`$-functions that have weak derivatives in $`L_2`$.

All our operators are semi-bounded, self-adjoint operators $`T`$ acting in some $`L_2`$-space; they are obtained from some closed, semi-bounded quadratic form $`𝔱`$ via the usual representation theorem (\[K, RS-I\]) and satisfy $`D(T)D(𝔱)`$ and the relation

$$𝔱[u,v]=Tu,v,uD(T),vD(𝔱).$$
$`(1.1)`$
As usual, we say $`T_2T_1`$ iff the associated quadratic forms satisfy $`D(𝔱_2)D(𝔱_1)`$ and $`𝔱_2[u,u]𝔱_1[u,u]`$, for all $`uD(𝔱_2)`$. We also write $`𝔱[u]:=𝔱[u,u]`$ in the sequel.

For a sequence of non-negative, self-adjoint operators $`T_n`$ in $``$ we say that $`T_n`$ converges to a limit $`T_{\mathrm{}}`$ in the *strong resolvent sense* iff

$$(T_n+I)^1f(T_{\mathrm{}}+I)^1f,f;$$
$`(1.2)`$
for convergence in the *norm resolvent sense* we require uniform convergence for all $`f`$ s.th. $`f1`$.

Most of our results hinge upon the following convergence theorem for quadratic forms (\[RS-I\], \[S\]):

###### 1.1 Theorem.

Let $`(T_n)`$ be a monotonic sequence of non-negative, self-adjoint operators. Then there exists a self-adjoint operator $`T_{\mathrm{}}0`$ s.th. $`T_nT_{\mathrm{}}`$ in strong resolvent sense.

For an *increasing* sequence of quadratic forms $`0𝔱_n𝔱_{n+1}`$ with domains $`𝒬_n`$ the limiting form domain is given by $`𝒬_{\mathrm{}}=\{u_{n=1}^{\mathrm{}}𝒬_n;sup𝔱_n[u]<\mathrm{}\}`$. Here it may happen that $`T_{\mathrm{}}`$ acts in a strictly smaller Hilbert space than the $`T_n`$. In this case, $`(T_{\mathrm{}}+I)^1`$ has to be complemented by the zero operator on a suitable subspace of $``$. For a *decreasing* sequence it is shown in \[RS-I, S\] that there exists a largest closable form that is smaller than $`𝔱_n`$, for all $`n`$.

In a periodic situation, strong resolvent convergence often implies convergence of spectral densities. For the definition of the *integrated density of states* (i.d.s.) we refer to \[RS-IV\], \[PF\] or \[DIT\]. In our context, the i.d.s. is a monotonically increasing function $`F:`$ that tends to zero at $`\mathrm{}`$. Such functions have at most a countable number of discontinuities. The associated Borel measure $`\mu `$ is then called the *density of states measure*. We will say that a sequence $`(\mu _n)`$ of density of states measures converges weakly to $`\mu _{\mathrm{}}`$ if $`F_n(t)F_{\mathrm{}}(t)`$ at all points of continuity $`t`$ of $`F_{\mathrm{}}`$. The proof of the subsequent proposition is elementary if one uses the definition of \[RS-IV\] for the i.d.s.

###### 1.2 Proposition.

In addition to the assumptions of Theorem 1.1, suppose that the operators $`T_n`$ are periodic (with respect to the the same lattice $`\mathrm{\Gamma }`$), and that a (finitely-valued) i.d.s. $`F_n`$ exist for $`T_n`$ as well as for $`T_{\mathrm{}}`$. Then the corresponding measures $`\mu _n`$ converge weakly to $`\mu _{\mathrm{}}`$.

It is well-known that *norm resolvent convergence* implies spectral convergence (cf. \[RS-I; Ch. VIII\]) on any compact interval of the real line. In particular, we have the following simple consequence for the existence of spectral gaps: suppose $`T_nT_{\mathrm{}}`$ in norm-resolvent sense. If $`T_{\mathrm{}}`$ has a gap around $`E`$, then the $`T_n`$ will also have a gap around $`E`$, for large $`n`$. More precisely, suppose that $`(a,b)\sigma (T_{\mathrm{}})=\mathrm{}`$ and let $`\epsilon >0`$. Then $`(a+\epsilon ,b\epsilon )\sigma (T_n)=\mathrm{}`$, for all $`nn_\epsilon `$.

## 2. Schrödinger operators

As a warm-up, we illustrate the above convergence theorems by some simple results on Schrödinger operators of the form $`\mathrm{\Delta }+\lambda V(x)`$, where $`V`$ is non-negative. It has been a “folk theorem” for quite some time that $`\mathrm{\Delta }+\lambda \chi _B`$ converges in strong resolvent sense to the Dirichlet Laplacian on the complement $`^mB`$, under suitable assumptions on $`B`$.

Suppose $`V(x)0`$ is a periodic, real-valued (continuous or measurable and bounded) potential s.th. $`\{x^mV(x)0\}=\mathrm{\Omega }`$. We let

$$H_\lambda =\mathrm{\Delta }+\lambda V(x),\lambda 0,$$
$`(2.1)`$
with associated quadratic forms $`𝔥_\lambda [u]=u^2+\lambda V|u|^2`$, for $`u^1`$.

To describe the limit, we will need here (and in the following sections) the *Dirichlet Laplacian* on the closed set $`M=^m\mathrm{\Omega }`$, denoted as $`\mathrm{\Delta }_M`$ (\[Hed, HZ\]). The operator $`\mathrm{\Delta }_M`$ is constructed from the quadratic form $`u^2`$ with domain

$$\stackrel{~}{}_0^1(M):=\{u^1(^m)u(x)=0\text{ a.e. in }M^C\}.$$
$`(2.3)`$
For sufficiently regular boundary of $`M`$, $`\mathrm{\Delta }_M`$ agrees with the usual Dirichlet Laplacian. In the general case, it is shown in \[HZ\] that there exists a Borel set $`M^{}`$, $`M^{}M^{}M`$, such that $`\mathrm{\Delta }_M`$ acts as a self-adjoint operator in the Hilbert-space $`L_2(M^{})`$. It is now easy to see that

$$\{u^1(^m)\underset{\lambda 0}{sup}𝔥_\lambda [u]<\mathrm{}\}=\stackrel{~}{}_0^1(M),$$
$`(2.4)`$
and so Theorem 1.1 implies that $`H_\lambda `$ converges in strong resolvent sense to $`\mathrm{\Delta }_M`$. (To be more precise, the resolvent of $`\mathrm{\Delta }_M`$ acts in $`L_2(M^{})`$ and has to be complemented by the zero operator on $`L_2(^mM^{})`$ in Eqn. (1.1).) Up to this point, periodicity was not required. If we now assume that $`\mathrm{\Omega }`$ is periodic, then for all operators involved the i.d.s. is well-defined, and we may also conclude that the i.d.s. converges.

As was observed in \[HH-1\], one can actually do better and upgrade to norm resolvent convergence. The proof combines uniform exponential decay estimates for the resolvent kernels with Schur’s Lemma and local compactness of $`\mathrm{\Delta }`$.

The above convergence results do not require $`\mathrm{\Omega }`$ to be connected. In our “model case” of Figure 1 where $`\mathrm{\Omega }`$ contains a neighborhood of $`Q`$, $`\mathrm{\Delta }_M`$ is a countable direct sum of copies of $`\mathrm{\Delta }_{MQ}`$. Therefore, the spectrum of $`\mathrm{\Delta }_M`$ is a discrete set with each point in the spectrum an eigenvalue of infinite multiplicity. We therefore see that the spectrum of $`H_\lambda `$ concentrates at a discrete set of points, and an arbitrarily large number of spectral gaps opens up, as $`\lambda \mathrm{}`$.

## 3. Magnetic Schrödinger and Pauli operators

Here we briefly mention some results of \[HH-1,2\], \[HN\], \[B-1,2\] on strongly coupled periodic magnetic fields.

(1) Let $`\stackrel{}{a}C^1(^m,^m)`$ be a magnetic vector potential and $`B=d\stackrel{}{a}`$ the associated magnetic field. We assume that $`B`$ is periodic w.r.t. the lattice $`\mathrm{\Gamma }`$, and that $`B(x)=0`$ for all $`xM`$, while $`B(x)0`$ for all $`x\mathrm{\Omega }`$. Note that periodicity of $`\stackrel{}{a}`$ is *not* required. The magnetic Schrödinger operator is then given as \[CFKS\]

$$H(\lambda \stackrel{}{a})=(\mathrm{i}\lambda \stackrel{}{a})^2,$$
$`(3.1)`$
and obtained via the closure of the forms

$$H(\lambda \stackrel{}{a})u,u=(\mathrm{i}\lambda \stackrel{}{a}(x))u^2,uC_c^{\mathrm{}}(^m).$$

Evidently, there is *no* monotonicity w.r.t. $`\lambda `$. However, in any open subset $`U`$ of $`\mathrm{\Omega }`$ where a fixed entry $`B_{ij}(x)`$ of the magnetic tensor is bounded away from zero, the Avron-Herbst-Simon bound \[AHS\] gives a *local* lower bound for the quadratic form, expressing the barrier-like effect of a strong magnetic field. If the vector potential $`\stackrel{}{a}`$ vanishes on $`M`$, the above results for $`\mathrm{\Delta }+\mu \chi _U`$, $`\mu \mathrm{}`$, can be combined with the Feynman-Kac-Itô-formula to obtain strong resolvent convergence of $`H(\lambda \stackrel{}{a})`$ to $`\mathrm{\Delta }_M`$, as $`\lambda \mathrm{}`$ (\[HH-1\]). For given magnetic field $`B`$, it depends on the topology of the set $`\mathrm{\Omega }`$ and on flux conditions whether a gauge exists s.th. $`d\stackrel{}{a}=B`$ and $`\stackrel{}{a}(x)=0`$ a.e. on $`M`$.

Under the above conditions, we also obtain an upgrade to norm-resolvent convergence, as before, and get the same spectral consequences as for $`\mathrm{\Delta }+\lambda \chi _\mathrm{\Omega }`$.

(2) The situation where $`\stackrel{}{a}(x)`$ is allowed to have non-zero values on $`M`$ has been analyzed in detail in \[HN\], using some basic homology theory and de Rham’s Theorem. It can be shown that, under suitable assumptions, the norm difference of the resolvents of $`H(\lambda \stackrel{}{a})`$ and $`H_M(\lambda \stackrel{}{a})`$ tends to zero, as $`\lambda \mathrm{}`$, where $`H_M(\lambda \stackrel{}{a})=(\mathrm{i}\lambda \stackrel{}{a})^2`$, acting in $`L_2(M^{})`$ with Dirichlet boundary conditions. As a consequence, one can conclude that the spectrum of $`H(\lambda \stackrel{}{a})`$ approaches a periodic or a quasi-periodic function of $`\lambda `$, as $`\lambda \mathrm{}`$. The precise result involves the flux of $`B`$ through the connected components of $`(\overline{M^{}})^C`$, the complement of the closure of the interior of $`M`$.

(3) For the Pauli operator in $`^2`$, we have the pair of operators

$$H_\pm (\lambda \stackrel{}{a})=H(\lambda \stackrel{}{a})B,$$
$`(3.2)`$
both acting in $`L_2(^2)`$. Here supersymmetry (cf. \[CFKS\]) implies

$$\sigma (H_+(\lambda \stackrel{}{a}))\{0\}=\sigma (H_{}(\lambda \stackrel{}{a}))\{0\}.$$
$`(3.3)`$
Note that the results in (1) and (2) remain unchanged if we include an electric background potential $`V`$ (e.g., $`V`$ bounded and measurable), while supersymmetry is destroyed upon addition of an electric potential.

Under suitable assumptions, one can show (\[B-1,2\]) that the operators $`H_\pm (\lambda \stackrel{}{a})`$ have a common gap: in the simplest case, we assume $`B0`$. A more delicate result requires that $`B`$ is positive in a neighborhood of $`Q`$ and that there exists a periodic gauge s.th. $`d\stackrel{}{a}=B`$. (By the divergence theorem, this gauge condition is equivalent to the *flux condition* $`_QB(x)𝑑x=0`$).

## 4. Periodic divergence type operators

Here we report on some results of \[HL-1,2\] concerning periodic two-component media with high contrast. In the following we always assume $`m2`$. The operators dicussed below occur in acoustics, heat conduction, and propagation of electro-magnetic waves in photonic crystals.

The quadratic form is given by the *Dirichlet integral*

$$𝔱_\lambda [u]:=(1+\lambda \chi _\mathrm{\Omega }(x))|u(x)|^2𝑑x,u^1,\lambda 0;$$
$`(4.1)`$
this defines a monotonically increasing family of forms. The associated self-adjoint operators $`T_\lambda `$ can formally be written as $`T_\lambda =(1+\lambda \chi _\mathrm{\Omega })`$.

It is clear that only such functions $`u`$ from the form domain survive taking $`\lambda \mathrm{}`$ that are constant on each connected component of $`\mathrm{\Omega }`$. If $`\mathrm{\Omega }`$ is connected, this implies that only such $`u`$’s will remain in the limiting form domain $`𝒬_{\mathrm{}}`$ that are identically zero on $`\mathrm{\Omega }`$. We have:

###### 4.1 Theorem.

\[HL-1\] Suppose $`\mathrm{\Omega }`$ is open, periodic, and connected. Then $`T_\lambda `$ converges to the Dirichlet-Laplacian $`\mathrm{\Delta }_M`$ on $`M=^m\mathrm{\Omega }`$ in strong resolvent sense. Furthermore, the associated density of states measures converge weakly.

The standard Floquet-decomposition gives a more detailed picture. Let $`T_\lambda ^{(\theta )}`$ denote the operator $`(1+\lambda \chi _\mathrm{\Omega })`$, acting in $`L_2(Q)`$, with $`\theta `$-periodic boundary conditions on $`Q`$, for $`\theta (\pi ,\pi ]^m`$; cf. \[RS-IV, Ku\]. Then Theorem 1.1 can be applied in each fiber, i.e., for each fixed $`\theta `$. As each $`T_\lambda ^{(\theta )}`$ has compact resolvent, we even obtain norm resolvent convergence in the fibers. The limiting operators $`T_{\mathrm{}}^{(\theta )}`$ have form domains $`𝒬_{\mathrm{}}^{(\theta )}`$ consisting of those functions $`u^1(Q)`$ that are constant in each connected component of $`\mathrm{\Omega }Q`$ and satisfy $`\theta `$-periodic boundary conditions. If $`\mathrm{\Omega }`$ is connected and intersects each face of the period cell $`Q`$, it follows that such $`u`$ must vanish on $`\mathrm{\Omega }Q`$, provided $`\theta (0,\mathrm{},0)`$. Pursuing this observation in detail, one obtains

###### 4.2 Theorem.

(\[HL-1\]) Suppose $`\mathrm{\Omega }^m`$ is open, periodic, connected and contains $`Q`$. We also assume that $`M^{}`$ (defined in Section 2) is not of measure zero. Let $`0<\delta _1\delta _2\mathrm{}`$ denote the eigenvalues of $`\mathrm{\Delta }_{MQ}`$. We then have:

$`(i)`$ The first spectral band of $`T_\lambda `$ converges to the interval $`[0,\delta _1]`$.

$`(ii)`$ For $`\lambda `$ large, $`T_\lambda `$ has a spectral gap following the first band.

$`(iii)`$ The i.d.s. measure of $`T_\lambda `$ concentrates at the set $`\{\delta _k;k\}`$, as $`\lambda \mathrm{}`$.

It follows from $`(i)`$ that there is *no* norm resolvent convergence of the family $`(T_\lambda )`$ (The only possible candidate for the limit is $`\mathrm{\Delta }_M`$ which has no spectrum below $`\delta _1>0`$. But norm resolvent convergence implies spectral convergence.)

For smoothly bounded $`M`$ a much more detailed picture of the process of spectral concentration is given in \[F\].

###### 4.3 Remark.

More general results can be found in \[HL\]. First of all, without additional effort one can handle elliptic operators of the type $`(𝐚(x)+\lambda 𝐛(x))`$, under suitable assumptions on the (symmetric) matrices $`𝐚(x)`$ and $`𝐛(x)`$, like $`a_{ij},b_{ij}L_{\mathrm{}}(^m)`$, $`𝐚(x)`$ uniformly positive definite, $`𝐛(x)=\mathrm{𝟎}`$ on $`M`$, and positive definite for $`x\mathrm{\Omega }`$.

Second, \[HL\] also gives criteria for the existence of higher gaps: if the eigenspace of the $`k`$-th Dirichlet eigenvalue $`\delta _k`$ of $`\mathrm{\Delta }_{MQ}`$ contains an eigenfunction $`d_k`$ s.th. $`d_k0`$ then $`\sigma (T_\lambda )`$ will have a gap for $`\lambda `$ large that converges to $`(\delta _k,\nu _k)`$, for some $`\nu _k>\delta _k`$.

###### 4.4 Remark.

There are natural and simple situations where one is led to consider period cells that are more general than $`(0,1]^m`$. In fact, all of the above works if $`\mathrm{\Gamma }`$ is a discrete sublattice of $`^m`$ that spans $`^m`$ and if $`Q`$ is a contractible subset of $`^m`$ with piecewise smooth boundary tessalating $`^m`$ via the translations $`\gamma \mathrm{\Gamma }`$. Here we think, in particular, of the class of results that depend on the assumption that $`\mathrm{\Omega }`$ contains the boundary of the period cell.

###### 4.5 Remark.

In $`^m`$, $`m3`$, it is possible to have both $`\mathrm{\Omega }`$ and $`M`$ connected. A *sponge* is a typical example (albeit, in reality, not periodic but random!). Here our analysis would yield that, for high contrast, the i.d.s. is approximately given by the i.d.s. of $`T_{\mathrm{}}`$. The latter can be further analyzed by Floquet-decomposition if we assume that $`\mathrm{\Omega }`$ intersects each face of $`Q`$ (think of the case where $`M`$ is connected from one cell to the neighboring cells through thin filaments only).

###### 4.6 Remark.

It is natural to ask what happens for a *decreasing* sequence of operators like

$$S_\lambda =\frac{1}{1+\lambda \chi _M},0\lambda \mathrm{}.$$
$`(4.2)`$
The associated quadratic forms $`𝔰_\lambda `$ are defined by $`𝔰_\lambda [u]=(1+\lambda \chi _M)^1|u|^2𝑑x`$ on $`^1(^m)`$. Evidently, the largest form that is smaller than $`𝔰_\lambda `$, for all $`\lambda 0`$, has domain $`^1(^m)`$ and takes the values $`_\mathrm{\Omega }|u|^2𝑑x`$. Let us call this form $`𝔰_{\mathrm{}}`$. It is easy to see that $`𝔰_{\mathrm{}}`$ is closable. Let $`S_{\mathrm{}}`$ denote the self-adjoint operator associated with the closure of $`𝔰_{\mathrm{}}`$. It follows from Simon’s theorem \[RS-I; Thm. S. 16\] that the operators $`S_\lambda `$ converge to $`S_{\mathrm{}}`$ in strong resolvent sense.

Let us assume, for the moment, that any $`u^1(\mathrm{\Omega })`$ can be extended to a function in $`^1(^m)`$ and that $`C_c^{\mathrm{}}(M^{})`$ is dense in $`L_2(M)`$. Then it is not difficult to show that the closure of the form $`𝔰_{\mathrm{}}`$ has domain $`^1(\mathrm{\Omega })L_2(M)`$, and takes the values $`𝔰_{\mathrm{}}[u]=_\mathrm{\Omega }|u(x)|^2𝑑x`$. Under these assumptions $`S_{\mathrm{}}`$ is the direct sum of the Neumann Laplacian on $`\mathrm{\Omega }`$ and the zero operator on $`L_2(M)`$. In particular, the i.d.s. of $`S_{\mathrm{}}`$ would be infinite on the positive real line, and it seems rather hopeless to extract any useful spectral information on $`S_\lambda `$. (Of course, the rôles of $`M`$ and $`\mathrm{\Omega }`$ are interchangeable here if the boundary of $`M`$ is smooth enough.)

###### 4.7 Remark.

In the theory of photonic crystals, one is led to study the “high-contrast/thin-wall limit” of the acoustic operator in $`^2`$. In typical applications, the optically dense medium occupies a connected, periodic set $`M`$ consisting of thin walls, cf. \[FK-1,2\]. In our notation, we would then consider the limit of $`T_\lambda =(1+\lambda \chi _\mathrm{\Omega })`$, as $`\lambda \mathrm{}`$. Note that now $`\mathrm{\Omega }`$ is *not* connected. Our results yield a strong resolvent limit and a Floquet-decomposition analogous to what we have discussed above. The limiting operators $`T_{\mathrm{}}^{(\theta )}`$ in the fibers now have a form domain consisting of functions that are constant on $`\mathrm{\Omega }Q`$, where the constant value will be non-zero, in general. Therefore, it is much harder to determine the spectrum of the $`T_{\mathrm{}}^{(\theta )}`$ than in the case of Eqn. (4.1). If $`\mathrm{\Omega }Q`$ exhausts “most of” $`Q`$ (thin wall asymptotics), then the analysis can proceed again \[FK-1,2\], \[KK\].

###### 4.8 Remark.

Random Media and the convergence of the i.d.s. for the operators of Eqn. (4.1) with a random geometry have been studied by K. Lienau \[L\]. Since we need a connected $`\mathrm{\Omega }`$ for our approach to work well, the ergodic process that generates the random set $`M`$ has to respect this property. Here various methods from stochastic geometry (\[SKM\]) can be employed.

## 5. Periodic Laplace-Beltrami operators

In this section we report on some related results on periodic Laplace-Beltrami operators, cf. \[P-1,2\]. Let $``$ be a non-compact Riemannian manifold of dimension $`d2`$. We call $``$ a *periodic* or *covering manifold* iff $`\mathrm{\Gamma }=^m`$ acts properly discontinuously and isometrically on $``$. Under these conditions, the quotient $`/\mathrm{\Gamma }`$ is again a Riemannian manifold which we assume to be compact.

We call $`Q`$ a *period cell* or *fundamental domain* of $``$ iff $`Q`$ is open, $`Q`$ does not intersect any other translate $`\gamma +Q`$ except for $`\gamma =0`$, and the union of all translates of $`\overline{Q}`$ covers the whole manifold $``$. Note that we assume here that $`Q`$ is an open set. One can apply Floquet decomposition in the same way as before; here, $`\vartheta `$-periodicity means that $`u(\gamma +x)=e^{i\vartheta \gamma }u(x)`$ for all $`x\overline{Q}`$ and $`\gamma \mathrm{\Gamma }`$ such that $`\gamma +x\overline{Q}`$ (plus an analogous condition for the first derivatives).

Denote by $`(g_{ij})`$ the metric tensor of $``$ in some chart. We consider conformal deformations of the given periodic metric, i.e., we set

$$g_{ij}(x;\lambda )=\varrho _\lambda ^2(x)g_{ij}(x).$$

Here, $`\varrho _\lambda `$ is a strictly positive, periodic, smooth function converging pointwise to the indicator function of a closed periodic set $`M`$ as $`\lambda \mathrm{}`$. In particular, we assume that

(1) $`M_0:=MQ`$ is compactly contained in $`Q`$ with smooth $`M`$;

(2) $`\varrho _\lambda (x)=1`$ for all $`xM`$;

(3) $`\varrho _\lambda (x)=\frac{1}{\lambda }`$, for all $`x\mathrm{\Omega }`$ with dist$`(x,M)\lambda ^m`$.

To simplify the notation, we restrict ourselves to the conformally flat case, i.e., we assume that $`=^m`$ with the Euclidean metric tensor $`g_{ij}(x)=\delta _{ij}`$, defined on $`^m`$. We assume that the boundary $`M`$ has a tubular neighborhood where we can introduce local coordinates $`(r,y)`$ given by the distance $`r`$ from $`M`$ and $`yM`$. For technical reasons we need the additional assumption that these coordinates can be extended to $`\mathrm{\Omega }Q`$; e.g., think of a ball $`M_0`$ centered in the period cell $`Q=(0,1)^m`$. The Laplace-Beltrami-operator $`B_\lambda `$ corresponding to the deformed metric is given by

$$B_\lambda =\frac{1}{\sqrt{g_\lambda }}_j\left(g^{ij}(;\lambda )\sqrt{g_\lambda }\right)_i$$

(using the summation convention), where $`g_\lambda =detg_{ij}(,\lambda )=\rho _\lambda ^{2m}`$ denotes the square of the volume density of the deformed metric tensor. This operator is defined as a non-negative, self-adjoint operator via the quadratic form

$$𝔟_\lambda [u]:=g^{ij}(,\lambda )_iu_j\overline{u}\sqrt{g_\lambda }dx=|u|^2\varrho _\lambda ^{m2}𝑑x,u^1(^m)$$
$`(5.1)`$
in the Hilbert space $`_\lambda :=L_2(^m)`$ with inner product $`u\overline{v}\sqrt{g_\lambda }𝑑x`$. Note that, in general, both the Dirichlet integral *and* the volume form depend on the parameter $`\lambda `$. The fundamental difference to the case of divergence type operators lies in the dependence of the inner product on $`\lambda `$. In particular, for $`m3`$, there is no monotonicity.

Floquet-decomposition implies that it is enough to study $`\vartheta `$-periodic b.v.p.’s on $`Q`$. Furthermore, Dirichlet-Neumann-bracketing yields enclosures of bands by Neumann and Dirichlet eigenvalues $`E_k^{N/D}(\lambda )`$ (ordered by min-max and repeated according to multiplicities). Similarly, denote by $`E_k^N`$ the Neumann eigenvalues of $`M_0=MQ`$.

Since the norm of the Hilbert space depends on $`\lambda `$, we cannot apply Theorem 1.1. Instead, the convergence of the eigenvalues on $`Q`$ to the Neumann eigenvalues on $`M_0`$ is obtained directly through an application of the min-max-principle and a careful comparison of the corresponding Raleigh quotients.

In particular, we show that $`E_k^N`$ is (approximately) a lower bound for $`E_k^N(\lambda )`$ (using the restriction $`u|_{M_0}`$ as transition operator from $`^1(Q)`$ to $`^1(M_0)`$). In the same way, we prove that $`E_k^N`$ is also an approximate upper bound for $`E_k^D(\lambda )`$ provided $`m3`$ (using an extension operator from $`^1(M_0)`$ to $`_0^1(Q)`$, the $`^1`$-closure of $`C_c^{\mathrm{}}(Q)`$, as transition operator). Therefore, we have

###### 5.1 Theorem.

(\[P-1,2\]) Let $`m3`$. Then $`\sigma (B_\lambda )`$ converges to the spectrum of the Neumann Laplacian on $`M`$, as $`\lambda \mathrm{}`$. In particular, gaps open up, as $`\lambda \mathrm{}`$.

###### 5.2 Remark.

The case $`m=2`$ is special. Here, only the inner product depends on $`\lambda `$. Inverting the rôle of the norm and the quadratic form yields a limiting operators in each fiber which, however, depend on $`\theta `$. In the special case of a $`2`$-dimensional cylinder one can nevertheless prove the existence of spectral gaps by direct calculations (separation of variables), cf. \[P-1,2\].

###### 5.3 Remark.

In contrast to the Euclidean case, we may choose the dimension $`d`$ of the manifold $``$ independently of the dimension of periodicity $`m`$: think of the surface of a cylinder ($`d=2,m=1`$) or the surface of a periodic jungle gym ($`d=2,m=3`$). Furthermore, one can allow any Abelian finitely generated group $`\mathrm{\Gamma }`$ of infinite order, i.e., products of $``$ and cyclic groups.

###### 5.4 Remark.

In some sense, the case of Laplace-Beltrami operators is the opposite of the case of divergence type operators: while the divergence type operator needs the factor $`(1+\lambda \chi _\mathrm{\Omega })`$ to be *large* on $`\mathrm{\Omega }`$ to produce spectral gaps, the Laplace-Beltrami operator needs the factor $`\varrho _\lambda ^{m2}`$ to be *small* on $`\mathrm{\Omega }`$.

###### 5.5 Remark.

Davies and Harrell \[DH\] proved the existence of at least one gap in the periodic conformally flat case (transforming the conformal Laplacian into a corresponding Schrödinger operator). This result is a special case of our result. Using similar methods as \[DH\], Green \[G\] obtained examples with a finite number of gaps in the $`2`$-dimensional conformally flat case. Recently, Yoshitomi \[Y\] proved the existence of spectral gaps for the (Dirichlet) Laplacian on periodically curved quantum wave guides. In \[DH, G, Y\], the existence of gaps is established by analyzing a one-dimensional problem. Here, in contrast, we directly study the multi-dimensional problem.

###### 5.6 Remark.

More “geometric” constructions of periodic manifolds with spectral gaps can be found in \[P-1,2\], think for example of an infinite number of copies of a compact manifold, joined by thin cylinders. If these cylinders are small or thin enough, one again obtains spectral gaps of the corresponding Laplace-Beltrami operator. These results are closely related to earlier works of Chavel and Feldman \[CF\] and Anné \[A\].

## 6. Maxwell operators

Another important elliptic operator is the *Maxwell operator* given as

$$M=\text{curl}\frac{1}{\epsilon (x)}\text{curl},$$
$`(6.1)`$
acting on solenoïdal vector fields in $`^2`$ or in $`^3`$, where $`\epsilon (x)`$ denotes the dielectricity of the medium.

The existence of gaps for periodic structures is of fundamental importance for photonic crystals. In many cases, one is satisfied to find intervals of low spectral density (quasi-gaps). This class of problems has been studied in depth in a series of papers by A. Figotin, P. Kuchment and others (cf., e.g., \[FK-1,2\], \[KK\]) where a complete analysis for high-contrast and/or thin-wall asymptotics is given.


format: text

text:
# Untitled Document

Education Curricula of the Regional Centres for Space Science and Technology Education (Affiliated to the United Nations)

Hans J. Haubold

Office for Outer Space Affairs, United Nations, Vienna International Centre, P.O. Box 500, 1400 Vienna, Austria (hans.haubold@unvienna.org)

ABSTRACT: Since 1988, the United Nations, through the Programme on Space Applications, is supporting the establishment and operation of regional Centres for Space Science and Technology Education in Africa, Asia and the Pacific, Latin America and the Caribbean, and Western Asia. Simultaneously, education curricula have been developed for remote sensing and geographic information system, satellite communications, satellite meteorology and global climate, and space and atmospheric science. The paper reviews briefly these developments and highlights the most recent updated education curricula in the four disciplines that are made available in 2002, in the six official languages of the United Nations, for implementation at the regional Centres and beyond.

Space science and technology education can be pursued at the elementary, secondary and university levels. In space faring nations, National Science Education Standards for space science and technology have been applied to science curricula at those levels (National Research Council, 1996). Such an innovation has not taken place in many developing nations, partly because the benefits of space science and technology have not been appreciated enough and partly because the facilities and resources for teaching science and technology at educational institutions are not yet well developed. Education in space science and technology in developed nations has become highly interactive; the World Wide Web and other information technologies have become useful tools in education programmes at all levels.

The incorporation of elements of space science and technology into university-level science curricula can serve a dual purpose for developed and developing nations (Daniel, 2002). It can enable nations to take advantage of the benefits inherent in the new technologies, which, in many cases, are spin-offs from space science and technology. It can also revitalize the educational system, introduce the concepts of high technology in a non-esoteric fashion and help create national capacities in science and technology in general. In that regard, Pyenson and Sheets-Pyenson (1999) emphasized in their recent work entitled ’Servants of Nature’ that:

“Both geographical decentralization and interdisciplinary innovation have become watchwords in academic science. Electronic information processing to some extent obviates the necessity for a scientist or scholar to reside at an ancient college of learning. Universities everywhere have adapted to new socioeconomic conditions by expanding curricula. They have always responded in this way, although never as quickly as their critics would like. Measured and deliberate innovation is one of academia’s heavy burdens. It is also a great strength. Emerging fields of knowledge become new scientific disciplines only after they have found a secure place in universities. We look to universities for an authoritative word about the latest innovations. New scientific ideas emerge in a variety of settings, but they become the common heritage of humanity only when processed by an institution for advanced instruction like the modern university.”

There are many challenges in the teaching of science at university level, both in developing and developed nations, but the challenges are of a higher magnitude in developing nations. The general problem confronting science education is the inability of students to see or experience the phenomena being taught, which often leads to an inability to learn basic principles and to see the relationship between two or more concepts and their practical relevance to problems in real life. Added to those problems is a lack of skills in the relevant aspects of mathematics and physics and in problem-solving strategies. There are also language problems in nations in which science is not taught in the national language(s). Over the years, developed nations have overcome most of the basic problems, except perhaps a psychological problem, namely that students may consider science to be a difficult subject. In developing nations, however, basic problems linger, exacerbated by the fact that there are not enough academically and professionally well-trained teachers.

The General Assembly of the United Nations, in its resolution 45/72 of 11 December 1990, endorsed the recommendation of the Working Group of the Whole of the Scientific and Technical Subcommittee, as endorsed by the Committee on the Peaceful Uses of Outer Space, that the United Nations should lead, with the active support of its specialized agencies and other international organizations, an international effort to establish regional centres for space science and technology education in existing national/regional educational institutions in the developing nations.111UN Document A/AC.105/456, annex II, para. 4 (n), Report of the Scientific and Technical Subcommittee on the Work of its Twenty-Seventh Session, 12 March 1990, United Nations, New York.

The General Assembly of the United Nations, in its resolution 50/27 of 6 December 1995, also endorsed the recommendation of the Committee on the Peaceful Uses of Outer Space that those centres be established on the basis of affiliation to the United Nations as early as possible and that such affiliation would provide the centres with the necessary recognition and would strengthen the possibilities of attracting donors and of establishing academic relationships with national and international space-related institutions.222UN Document A/RES/50/27, paragraph 30, International Cooperation in the Peaceful Uses of Outer Space, 5 February 1996, United Nations, New York.

Regional centres have been established in India for Asia and the Pacific333See web site of the regional Centre for Asia and the Pacific at http://www.cssteap.org/., in Morocco and Nigeria for Africa, in Brazil and Mexico for Latin America and the Caribbean, and in Jordan for Western Asia, under the auspices of the Programme on Space Applications, implemented by the Office for Outer Space Affairs.444UN Document A/AC.105/749, Regional Centres for Space Science and Technology Education (affiliated to the United Nations), 27 December 2000, United Nations, Vienna. The objective of the centres is to enhance the capabilities of Member States, at the regional and international levels, in various disciplines of space science and technology that can advance their scientific, economic, and social development. Each of the centres provides postgraduate education, research, and application programmes with emphasis on (i) remote sensing and geographic information system, (ii) satellite communications, (iii) satellite meteorology and global climate, and (iv) space and atmospheric science for university educators and research and application scientists. Since 1996, all centres are implementing nine-month postgraduate courses in the four areas based on model curricula that emanated from the United Nations Meeting of Experts on the Development of Education Curricula for the Regional Centres for Space Science and Technology Education, held in Granada, Spain, in 1995. Since 1995, the curricula555UN Document A/AC.105/649, Centres for Space Science and Technology Education: Education Curricula, United Nations, Vienna; the document is also available in English and French at http://www.oosa.unvienna.org/SAP/centres/centres.htm. have been reviewed at regional and international educational meetings.

The Third United Nations Conference on the Exploration and Peaceful Uses of Outer Space (UNISPACE III, 1999), recommended that collaboration should be established between the regional centres and other national, regional and international organizations to strengthen components of their education curricula.666Report of the Third United Nations Conference on the Exploration and Peaceful Uses of Outer Space, Vienna, 9-30 July 1999 (United Nations publication, Sales No. E.00.I.3), chap. II, sect. G, para. 220; the document is also available at http://www.oosa.unvienna.org/. In its resolution 54/68 of 6 December 1999, the General Assembly of the United Nations endorsed the resolution of UNISPACE III entitled ”The Space Millennium: Vienna Declaration on Space and Human Development”, in which action was recommended to ensure sustainable funding mechanisms for the regional centres.777Ibid., chap. I, resolution 1, para. 1 (e) (ii). The Declaration is also available on the home page of the Office for Outer Space Affairs (http://www.oosa.unvienna.org).

The Office for Outer Space Affairs of the Secretariat of the United Nations organized, in cooperation with the European Space Agency, a second United Nations Expert Meeting on the Regional Centres for Space Science and Technology Education: Status and Future Development in Frascati, Italy, in 2001.888UN Document A/AC.105/782, Regional Centres for Space Science and Technology Education (affiliated to the United Nations), 14 March 2002, United Nations, Vienna. The Meeting reviewed the status of establishment and operation of the regional centres with a view to enhancing cooperation between the centres. However, the main objective of the Meeting was to review and update the curricula at the university level and across cultures in the four areas. The Meeting considered that education varied significantly between nations and even between institutions within the same country which led to differences in space science and technology education curricula in terms of course content and modes of presentation. The Meeting noted that the model curricula999UN Document A/AC.105/649, Centres for Space Science and Technology Education: Education Curricula, United Nations, Vienna; the document is also available in English and French at http://www.oosa.unvienna.org/SAP/centres/centres.htm. had contributed to resolving such problems.

The Meeting updated the four education curricula and drew up course syllabuses that differ from most of those available in literature and on the World Wide Web. They are based on physics, mathematics and engineering as taught in many universities around the world. They are not tailored to any specific space-related project or mission that may have been or will be executed by any specific institution. The curricula are available, in the six official languages (Arabic, Chinese, English, French, Russian, Spanish) of the United Nations, on request, as United Nations documents.101010UN Documents A/AC.105/L.238 (satellite meteorology and global climate), A/AC.105/L.239 (satellite communications), A/AC.105/L.240 (space and atmospheric science), and A/AC.105/L.241 (remote sensing and geographic information system).

LIST OF REFERENCES

Daniel, R.R. (ed.) (2002). Concepts in Space Science, Universities Press, Orient Longman, Hyderabad, India.

National Research Council (1996). National Science Education Standards, National Academy Press, Washington, D.C.

Pyenson, L. and Sheets-Pyenson, S. (1999). Servants of Nature: A History of Scientific Institution, Enterprises, and Sensibilities, W. W. Norton and Company, New York.


text:
# Roughness of Interfacial Crack Front: Correlated Percolation in the Damage Zone

## Abstract

We show that the roughness exponent $`\zeta `$ of an in-plane crack front slowly propagating along a heterogeneous interface embeded in a elastic body, is in full agreement with a correlated percolation problem in a linear gradient. We obtain $`\zeta =\nu /(1+\nu )`$ where $`\nu `$ is the correlation length critical exponent. We develop an elastic brittle model based on both the 3D Green function in an elastic half-space and a discrete interface of brittle fibers and find numerically that $`\nu =1.5`$. We conjecture it to be 3/2. This yields $`\zeta =3/5`$. We also obtain by direct numerical simulations $`\zeta =0.6`$ in excellent agreement with our prediction. This modelling is for the first time in close agreement with experimental observations.

An important motivation for studying interfacial crack pinning srvm95 ; sm97  is to simplify the study of the origin of the scaling properties of brittle crack surfaces mpp84 ; bs85 . These scaling properties are seen, for example, in the height-height correlation of the fracture roughness (i.e. out-plane roughness), which shows self-affinity. That is, the conditional probability density $`p(x,h)`$, i.e. the probability that the crack surface passes within $`dh`$ of the height $`h`$ at position $`x`$ when it had height zero at $`x=0`$, shows the invariance

$$\lambda ^\zeta p(\lambda x,\lambda ^\zeta h)=p(x,h),$$
(1)
where $`\zeta `$ is the roughness exponent. It is now generally believed that the roughness exponent shows a universal value of about 0.80 at larger scales blp90-mhhr92-sgr93-cw93-sss95 , while a lower value of about 0.5 might be seen at smaller scales dhbc96-dnbc97 .

Direct observations of the interfacial crack-pinning have been performed recently. The problem consists of following the roughness of a crack front moving along the flat interface between two elastically connected blocks. The experimental study of constrained crack propagation between two sintered plexiglass plates presented in sm97  resulted in the estimate of the in-plane roughness exponent: $`\zeta =0.55\pm 0.05`$. This work was followed by a longer study leading to the estimate $`\zeta =0.63\pm 0.03`$ dsm99 . A recent study of the motion of a helium-4 meniscus along a disordered substrate — a problem closely related to the motion of a crack line — gave $`\zeta =0.56\pm 0.03`$ prg02 .

Numerous models for interfacial crack propogation in heterogeneous material have been proposed. The numerical simulation presented in srvm95  is based on a perturbative Green function approach following the quasistatic evolution of the interfacial crack front position $`a(x,t)`$ gr89 , which is treated as a single-valued function of position $`x`$ along the orthogonal direction to the crack advancement direction, and time $`t`$. The linearized Green function used binds together points only along the fracture front. The stress intensity factor at a point $`x`$ along the fracture front is then found to be

$$\frac{K(x,t)}{K_0}=1+\frac{1}{2\pi }\mathrm{p}.\mathrm{v}._{\mathrm{}}^+\mathrm{}\frac{a(x^{},t)a(x,t)}{(xx^{})^2}𝑑x^{},$$
(2)
where $`K_0`$ is the stress intensity factor that would result if the crack were straight gr89 . The fracture is advanced by identifying the most stressed point along the fracture line and advancing this by a small step. The roughness exponent of the crack front was estimated numerically to be $`\zeta =0.35`$, while a direct dynamical renormalization group calculation gave $`\zeta =1/3`$ to lowest order ek94 . Higher-order corrections to this quasistatic analysis increases the value of $`\zeta `$ to 0.48 cdw01 , while a different quasistatic analytical technique suggests $`0.390\pm 0.002`$ rk02 . Dynamic effects have beeen largely studied numerically and analytically pr94-rf97-rf98-mr98-mr00 , in particular in the form of crack front waves. They lead to an increased roughness exponent compared to the initial 1/3-value up to $`\zeta =0.5`$ but still smaller than experimental observations.

Hence, the situation today is that there is a large and significant discrepancy between theoretical and experimental estimates: theoretical estimates being considerably smaller than the experimental ones.

In this letter we present a numerical calculation of the roughness exponent of a crack front that propagates quasi-statically along a heterogeneous interface, based on a Green function technique which differs from that previously used. Indeed it does not reduce the crack tip to a single tortuous line but describes the tip as a region of interactions between microcracks (see Fig. 1). Our technique is based on the static solution of the elastic equations for the deformation of the surface of an elastic half space ll58 . The local deformation, $`u`$, at position $`(x,y)`$ along the plane is related to the normal stress field $`\sigma `$ by the expression

$$u(x,y)=G(xx^{},yy^{})\sigma (x^{},y^{})𝑑x^{}𝑑y^{},$$
(3)
where the Green function is j85 :

$$G(x,y)=\frac{1s^2}{\pi e}\frac{1}{|(x,y)|},$$
(4)
with $`s`$ is the Poisson ratio and $`e`$ the elastic constant. For large $`r=|(xx^{},yy^{})|`$, $`G1/r`$. For the sake of simplicity, as generally done for the study of contact between two elastic bodies j85 , we substitute one of the elastic plates by an infinitely rigid one. The other plate is modeled as an elastic block for which Eq. (4) is valid. We discretize the model, so that the forces and deformations are described by the discrete version of Eq. (3),

$$u_i=\underset{j}{}\overline{G}_{i,j}f_j.$$
(5)
$`\overline{G}_{i,j}`$ is the Green function Eq. (4) averaged over an area $`b^2`$,

$$\overline{G}_{i,j}=\frac{1}{b^2}_{b/2}^{+b/2}_{b/2}^{+b/2}𝑑x^{}𝑑y^{}G(ix^{},jy^{}),$$
(6)
where $`b`$ is the lattice spacing, $`u_i`$ is the deformation of the elastic body at site $`i`$, and $`f_i`$ the force acting at that point. The indices $`i`$ and $`j`$ in Eq. 5 run over all $`L^2`$ sites.

The elastic block is connected to the infinitely stiff plate by a discrete interface made of an array of elastic harmonic springs. Each spring is brittle and has a breaking threshold randomly drawn from a uniform distribution between zero and one. The spacing between the springs is $`b`$ in both the $`x`$ and $`y`$ directions. The force $`f_i`$ that an unbroken spring $`i`$ is carrying, is transferred over an area of size $`b^2`$ to the soft surface and given by Hooke’s law:

$$f_i=k(u_iD),$$
(7)
where $`k`$ is the spring constant ($`k=1`$ for all springs). $`D`$ is the displacement of the infinitely stiff medium, and is a function of $`y`$, i.e. solid rotation of the stiff medium. The quantity $`(u_iD)`$ is, therefore, the length that spring $`i`$ is stretched. We assume periodic boundary conditions both in the $`x`$ and $`y`$ directions. In order to model the mode I fracturing of the interface between the two media (Double Cantilever Beam (DCB) test) in a way compatible with the biperiodic boundary conditions, we let $`D(y)`$ be wedgeshaped with a maximum at $`y=0`$ and $`L`$, and zero for $`y=L/2`$ (i.e. mirroring conditions). The load the system, $`D(y)`$ is increased and the springs break one by one. The numerical technique to solve the equations that ensue is presented in detail in bhs02 .

We show in Fig. 1 a typical damage front in a $`128\times 128`$ system. We define the fracture front $`a(x,t)`$ in this model as the set of nodes that form a continuous path separating the infinite cluster of broken springs from the infinite cluster of unbroken springs. This definition is similar to that of Ref. sm97 . From Fig. 1, we clearly see that even if a front can be defined, an extended damage zone exists. Accordingly the front does not capture all the active tip of the fracture. This observation supports the relevance of a non line description of the crack front.

In Fig. 2, we show the width of the fracture front $`w=\sqrt{a^2a^2}`$ as a function of its average position $`a`$ — which acts as a time parameter in this quasi-static model, for various system sizes. By collapsing the width evolution for the different system sizes, we see that the crack front follows a Family-Vicsek scaling fv85 . Two important exponents are, thus, obtained: the roughness exponent $`\zeta =0.6`$ and dynamical exponent $`z=0.9`$. The roughness exponent is in close agreement with the experimentally obtained value sm97 ; dsm99 ; ms01 ; prg02 , while the dynamical exponent $`z`$ was found in Ref. ms01  to be $`1.2`$.

Our model distinguishes itself from earlier numerical models in three major ways: (1) Most previous models are based on a small perturbation approach and include linearizations gr89  that are not used here. Furthermore, it was assumed in the earlier studies that (2) the fracture line was a single-valued function, hence ruling out overhangs. This assumption prevents islands of unbroken bonds from forming in the wake of the advancing fracture line. (3) The assumption of a single advancing fracture line also prevents the formation of a damage zone in front of the fracture line. None of these three assumptions are necessary in the present model. In order to test whether assumptions (2) or (3) are responsible for the difference in roughness exponent, we imposed both on the present model. No change in $`\zeta `$ was observed. Hence, we conclude that it is the linearization assumption that is responsable for the difference.

In the following, we demonstrate that the roughness exponent $`\zeta `$ is a result of a correlated percolation process in a gradient imposed by the mode I loading of the system. We start by charactarizing the correlated percolation. For this goal, we consider a similar problem but where the loading is obtained without any gradient, i.e. horizontal and parallel displacement of the rigid medium bhs02 . In this case, when the homogeneous displacement $`D`$ of the rigid medium reaches a maximum value, $`D_c`$ the system goes unstable, and unless $`D`$ is decreased again, catastrophic failure sets in. In Ref. bhs02 , the size distribution of clusters of broken springs was studied, and a power law was found with an exponent $`\tau =1.6`$. This value is different from ordinary percolation where $`\tau 2.05`$ sa92  and shows that correlated percolation takes place. In Fig. 3, we show the fluctuations of density of broken springs at $`D=D_c`$. If there is a diverging correlation length in the problem, these fluctuations scale as $`L^{1/\nu }`$. We find that $`1/\nu =0.65`$, leading to $`\nu =1.54`$. We conjecture that the exact value is $`\nu =3/2`$. Hence, the fracture process in this system is in a universality class which is different from standard percolation where $`\nu =4/3`$.

When $`D`$ is no longer constant, but is given by the wedge shape (DCB load), the form of interactions between the springs does not change. Hence, the critical properties of the model when run in the “horizontal” mode are still present under DCB loading. The wedge shape of $`D`$ leads to a gradient in the loading of the springs, going from very high loading where the damage is large to very low loading well into the still intact part of the interface. In Fig. 4, we show the damage profile $`p(y)`$, i.e. the density of broken springs averaged in the $`x`$ direction, for systems of different sizes. The profile is clearly linear. Somewhere along this damage profile, there is a line in the $`x`$-direction at $`y=y_c`$ where the damage density is critical, $`p(y_c)=p_c`$. In the vicinity of this line, there is a critical region which is characterized by being on the edge between stability and instability and corresponds to the crack front. Following the arguments of Sapoval et al. for percolation in a gradient srg85 , if $`p(y)`$ follows a linear law,

$$p(y)=1\frac{yy_0}{l_y},$$
(8)
on the interval $`y_0yy_0+l_y`$, where $`l_y`$ is the length scale characterizing the width of the damage zone, then $`y_c=y_0+l_y(1p_c)`$ is the position of the fracture front (when ignoring small corrections for finite-size systems). The critical region has a width $`\xi =|y_wy_c|`$ which is related to the damage $`p(y)`$ as $`\xi |p(y_w)p_c|^\nu `$. Eliminating $`|y_cy_w|`$ between these two expressions for $`\xi `$ gives

$$\xi l_y^{\frac{\nu }{1+\nu }}.$$
(9)
Since the gradient in $`D`$ is inversely proportional to $`L`$ and the average strength of the springs does not change with $`L`$, we have that

$$l_yL.$$
(10)
In Fig. 4, the reduced variable $`y/L`$ was used resulting in data collapse for different system sizes, thus validating Eq. (10). Furthermore, the width of the fracture front, $`w`$ is proportional to the width of the critical region $`\xi w`$. Hence, we find

$$wL^{\frac{\nu }{1+\nu }}=L^{3/5}.$$
(11)
where we have used our estimate $`\nu =3/2`$. This result is in excellent agreement with our numerical simulations and with the experimental results.

A similar idea has been proposed for the outplane roughness of fracture surfaces where the gradient and the $`\nu `$ exponent are different and leads to $`\zeta =4/5`$ in excellent agreement with the experimental value 0.80 hs02 . We finally note that Zapperi et al. zhr00  have pointed out the connection between gradient percolation and the interface fracture problem in the limit of infinitely stiff plates. We have shown in this letter that this mechanism stays intact also when the plates do respond elastically. However, the universality class of the corresponding correlated percolation problem is different from ordinary percolation, leading to the observed roughness exponent $`\zeta =3/5`$ rather than 4/7 for ordinary percolation in a gradient.

This work was partially funded by the CNRS PICS contract $`\mathrm{\#}753`$ and the Norwegian Research Council, NFR. J.S. and A.H. thank Fernando A. Oliveira and the ICCMP for financial support during our stay in Brasília. Discussions with K. J. Måløy are gratefully acknowledged.


text:
# Imaging the Near Field

## I Introduction

Conventional optics is a highly developed subject, but has limitations of resolution due to the finite wavelength of light. It has been thought impossible to obtain images with details finer than this limit. Recently it has been shown that a ‘perfect lens’ is in principle possible and that arbitrarily fine details can be resolved in an image provided that the lens was constructed with sufficient precision. The prescription is simple: take a slab of material, thickness $`d`$, and with electrical permittivity and magnetic permeability given by,

$$\epsilon =1,\mu =1$$
(1)
Given that these conditions are realised, the slab will produce an image of any object with perfect resolution. The key to this remarkable behaviour is that the refractive index of the slab is,

$$n=\sqrt{\epsilon \mu }=1$$
(2)
It was Veselago in 1968  who first realised that negative values for $`\epsilon ,\mu `$ would result in a negative refractive index and he also pointed out that such a negative refractive material (NRM) would act as a lens but it took more than 30 years to realise the concept of negative refractive index at microwave frequencies .

It was only in recent times  that the lens’s remarkable property of perfect resolution was noted. For the first time there is the possibility of manipulating the near field to form an image. The physics of negative refractive index has caught the imagination of the physics community as evidenced by the publications in the past two years

Although the conditions for a perfect lens are simple enough to specify, realising them is in practice rather difficult. There are two main obstacles. First the condition of negative values for $`\epsilon ,\mu `$ also implies that these quantities depend very sensitively on frequency so that the ideal condition can only be realised at a single carefully selected frequency. Second it is very important that absorption, which shows up as a positive imaginary component of $`\epsilon `$ or $`\mu `$, is kept to a very small value. Resolution of the lens degrades rapidly with increasing absorption. It is the objective of this paper to explore how the effects of absorption can be minimised.

Let us probe a little deeper into the operation of the perfect lens. Any object is visible because it emits or scatters electromagnetic radiation. The problem of imaging is concerned with reproducing the electro-magnetic field distribution of objects in a two dimensional (2-D) plane in the 2-D image plane. The electromagnetic field in free space emitted or scattered by a 2-D object (x-y plane) can be conveniently decomposed into the Fourier components $`k_x`$ and $`k_y`$ and polarization defined by $`\sigma `$:

$$E(x,y,z;t)=\underset{k_x,k_y,\sigma }{}E(k_x,k_y,k_z)\mathrm{exp}[i(k_xx+k_yy+k_zz\omega t)],$$
(3)
where the source is assumed to be monochromatic at frequency $`\omega `$, $`k_x^2+k_y^2+k_z^2=\omega ^2/c^2`$ and $`c`$ is the speed of light in free space. Obviously when we move out of the object plane the amplitude of each Fourier component changes (note the z-dependence) and the image becomes blurred. The electromagnetic field consists of a radiative component of propagating modes with real $`k_z,`$ and a near-field component of non-propagating modes with imaginary $`k_z`$ whose amplitudes decay exponentially with distance from the source. Provided that $`k_z`$ is real, $`\omega ^2/c^2>k_x^2+k_y^2`$, it is only the phase that changes with z and a conventional lens is designed to correct for this phase change. The evanescent near-field modes are the high-frequency Fourier components describing the finest details in the object and to restore their amplitudes in the image plane requires amplification, which is of course beyond the power of a conventional lens and hence the limitations to resolution.

Thus the Perfect lens performs the dual function of correcting the phase of the radiative components as well as amplifying the near-field components bringing them both together to make a perfect image and thereby eliminating the diffraction limit on the image resolution. In general the conditions under which this perfect imaging occurs are :

$$\epsilon _{}=\epsilon _+,\mu _{}=\mu _+,$$
(4)
where $`\epsilon _{}`$ and $`\mu _{}`$ are the dielectric permittivity and magnetic permeability of the NRM slab, and $`\epsilon _+`$ and $`\mu _+`$ are the dielectric permittivity and magnetic permeability of the surrounding medium respectively.

An important simplification of these conditions can be had in the case that all length scales are much less than the wavelength of light. Under these circumstances electric and magnetic fields decouple: the P-polarised component of light becomes mainly electric in nature, and the S-polarised component mainly magnetic. Therefore in the case of P-polarised light we need only require that $`\epsilon =1`$, and the value of $`\mu `$ is almost irrelevant. This is a welcome relaxation of the requirements especially at optical frequencies where many materials have a negative values for $`\epsilon `$, but show no magnetic activity. We shall concentrate our investigations on these extreme near field conditions and confine our attentions to P-polarised light.

In Section-2, we investigate the properties of a layered structure comprising extremely thin slabs of silver and show that layered structures are less susceptible to the degrading effects of absorption, than are single element lenses. In section-3, we present some detailed calculations of how the multilayer lens transmits the individual Fourier components of the image.

## II The layered perfect lens - an unusual effective medium

Reference to figure 2 shows that extremely large amplitudes of the electric field occur within the lens when the near field is being amplified. This is especially true for the high frequency Fourier components which give the highest resolution to the image. Unless the lens is very close to the ideal lossless structure, these large fields will result in dissipation which will kill the amplifying effect. However there is a way to restructure the lens to ameliorate the effects of dissipation. We observe that in the ideal lossless case we can perfectly well divide the lens into separate layers each one making its contribution to the amplification process (Shamomina et al have made a similar observation  and Zhang et al. have considered a similar system ). Provided that the total length of vacuum between the object and image is equal to the total length of lens material, the lens will still work and produce a perfect image. However this subdivision of the lens makes a big difference to how the lens performs when it is less than ideal and absorption is present. The point is that by distribution the amplification, the fields never grow to the extreme values that they do when the lens is a single slab and therefore the dissipation will be much less. Figure 3 illustrates this point.

First let us estimate the resolution of a lens constituted as a single slab. According to our original calculations  in the near field limit the transmission coefficient through the lens for each Fourier component is,

$$\frac{\mathrm{exp}\left(2\sqrt{k_x^2+k_y^2}d\right)}{\frac{1}{4}\left(\epsilon _{}^{^{\prime \prime }}\right)^2+\mathrm{exp}\left(2\sqrt{k_x^2+k_y^2}d\right)},$$
(5)
where,

$`\epsilon _{}=\epsilon _{}^{^{}}+i\epsilon _{}^{^{\prime \prime }}`$

Obviously when,

$$\sqrt{k_x^2+k_y^2}d<\mathrm{ln}\left[\frac{\epsilon _{}^{^{\prime \prime }}}{2}\right]$$
(6)
the lens’ power to amplify begins to fall away. Fourier components of higher spatial frequency do not contribute and hence the resolution is limited to,

$$\mathrm{\Delta }=\frac{2\pi d}{\mathrm{ln}\left[\frac{\epsilon _{}^{^{\prime \prime }}}{2}\right]}$$
(7)
The easiest way to investigate the properties of a layered system is to recognise that, provided that the slices are thin enough, it will behave as an effective anisotropic medium whose properties we calculate as follows. Applying a uniform displacement field, $`D`$, perpendicular to the slices gives electric fields of $`\epsilon _0^1\epsilon _+^1D`$ and $`\epsilon _0^1\epsilon _{}^1D`$ in the positive dielectric medium and in the negative material of the lens respectively. Therefore the average electric field is given by,

$$E=\frac{1}{2}\left(\epsilon _0^1\epsilon _+^1D+\epsilon _0^1\epsilon _{}^1D\right)=\epsilon _0^1\epsilon _z^1D,$$
(8)
where

$$\epsilon _z^1=\frac{1}{2}\left(\epsilon _+^1+\epsilon _{}^1\right),$$
(9)
is the effective dielectric function for fields acting along the $`z`$-axis. By considering an electric field along the $`x`$-axis we arrive at

$$\epsilon _x=\frac{1}{2}\left(\epsilon _++\epsilon _{}\right),$$
(10)
where $`\epsilon _x`$ is the effective dielectric function for fields acting along the $`x`$-axis. We have assumed for simplicity that the thickness of each material component is the same, but it is also possible to have unequal thicknesses. Now under the perfect lens conditions, $`\epsilon _{}=\epsilon _+`$, we have

$$\epsilon _z\mathrm{},\epsilon _x0$$
(11)
Thus the stack of alternating extremely thin layers of negative and positive refractive media in the limiting case of layer thickness going to zero behaves as a highly anisotropic medium.

Radiation propagates in an anisotropic medium with the following dispersion,

$$\frac{k_x^2+k_y^2}{\epsilon _z}+\frac{k_z^2}{\epsilon _x}=\frac{\omega ^2}{c^2}$$
(12)
and hence for the perfect lens conditions it is always true that,

$$k_z=0$$
(13)
Each Fourier component of the image passes through this unusual medium without change of phase or attenuation. It is as if the front and back surfaces of the medium were in immediate contact.

Here we have a close analogy with an optical fibre bundle where each fibre corresponds to a pixel and copies the amplitude of the object pixels to the image pixels without attenuation and with the same phase change for each pixel, preserving optical coherence. Our layered system performs exactly the same function with the refinement that in principle the pixels are infinitely small, and the phase change is zero. In figure 5 we illustrate this point with an equivalent system: an array of infinitely conducting wires embedded in a medium where $`\epsilon =0`$. In the latter case it is more obvious that an image propagates through the system without distortion. Indeed in the trivial zero frequency limit the system simply connects object to image point by point.

Coming back to our point that the layered system reduces the effect of absorption, we estimate the transmission for P-polarized light through such a system in the near field limit as

$$T_P\frac{1}{\mathrm{cos}\left(\frac{i}{2}\epsilon _{}^{^{\prime \prime }}k_x2d\right)+\frac{1}{2}\left(\epsilon _++\epsilon _+^1\right)\mathrm{sin}\left(\frac{i}{2}\epsilon _{}^{^{\prime \prime }}k_x2d\right)}.$$
(14)
Evidently for small values of $`k_x`$ the transmission coefficient is unity and these Fourier components contribute perfectly to the image, but for large values of $`k_x`$ transmission is reduced. We estimate the resolution limit to be,

$$\mathrm{\Delta }\frac{2\pi }{k_{\mathrm{max}}}=2\pi \left(\frac{1}{2}\epsilon _{}^{^{\prime \prime }}2d\right)=2\pi \epsilon _{}^{^{\prime \prime }}d.$$
(15)
Therefore the smallest detail resolved by the lens decreases linearly with decreasing absorption ($`\epsilon _{}^{^{\prime \prime }}`$). In contrast the original single slab of lens had a much slower improvement of resolution, being inversely as $`\mathrm{ln}\left[\epsilon _{}^{^{\prime \prime }}\right]`$. Thus it appears to be a case of two lenses are better than one but many lenses are the best of all.

## III Image Simulations for a multilayer stack

In the previous section we gave some qualitative arguments as to the properties of metal-dielectric multilayer stacks and is clear that for P-polarized light in the quasi-static limit this structure would behave as a near-perfect ‘fibre optic bundle’. In the electrostatic (magnetostatic) limit of large $`k_xk_zq_z`$, there is no effect of changing $`\mu `$($`\epsilon `$) for the P(S)-polarization. The deviation from the quasi-static limit caused by the non-zero frequency of the electromagnetic wave would, however, not allow this decoupling. When the effects of retardation are included, a mismatch in the $`\epsilon `$ and $`\mu `$ from the perfect-lens conditions would always limit the image resolution and also leads to large transmission resonances associated with the excitation of coupled surface modes that could introduce artifacts into the image. For the negative dielectric (silver) lens, the magnetic permeability $`\mu =1`$ everywhere, and this is a large deviation from the perfect lens conditions. The dispersion of these coupled slab plasmon polaritons and their effects on the image transfer has been extensively studied in Ref. .

Essentially, for a single slab of negative dielectric material which satisfies the conditions for the existence of a surface plasmon on both the interfaces, the two surface plasmon states hybridise to give an antisymmetric and a symmetric state, whose frequencies are detuned away from that of a single uncoupled surface state. The transmission as a function of the transverse wave-vector remains reasonably close to unity up to the resonant wave-vector for the coupled plasmon state, after which it decays exponentially with larger wavevectors. The secret for better image resolution is to obtain a flat transmission coefficient for as large a range of wave vectors as possible. This is possible by using a thinner slab in which case the transmission resonance corresponding to a coupled slab mode occurs at a much larger $`k_x`$. For the transfer of the image over useful distances, we would then have to resort to a layered system of very thin slabs of alternating positive and negative media.

Let us now consider a layered system consisting of thin slabs of silver (negative dielectric constant $`\epsilon _{}`$) and any other positive dielectric medium ($`\epsilon _+`$). Since the dielectric constant of silver is dispersive***An empirical form for the dielectric constant of silver in the visible region of the spectrum is $`\epsilon _{}(\omega )=5.79.0^2(\mathrm{}\omega )^2+i0.4`$, ($`\mathrm{}\omega `$ in eV). The imaginary part can be taken to be reasonably constant in this frequency range., we can choose the frequency ($`\omega `$) of the electromagnetic radiation so as to satisfy the perfect lens condition at the interfaces between the media ($`\epsilon _{}(\omega )=\epsilon _+(\omega )`$). We use the transfer matrix method to compute the transmission through the layered medium as a function of the transverse wave-vector at a frequency at which the perfect lens condition is satisfied. We will denote by $`N`$ the number of slabs with negative dielectric constant in the alternating structure, each period consisting of a negative and positive slab as shown in figure 4. Now the total length of the system is $`2d=N\delta `$ where $`\delta `$ is the period of the multilayer stack (the negative and positive slabs being of equal thickness of $`\delta /2`$). Note that the total thicknesses of positive and negative dielectric media between the object plane to the image plane are also equal.

The transmission across the multilayer system is shown in figure 6, where the thickness of the individual slabs is kept constant, but the number of layers is increased, thereby increasing the total length of the system. We get divergences in the transmission at wave-vectors corresponding to the coupled plasmon resonances. The number of the resonances increases with the number of layers, corresponding to the number of surface modes at the interfaces. For the system with the (hypothetical) lossless negative media, one notes that as we increase the number of layers, the transmission coefficient is almost constant and close to unity with increasing $`k_x`$, until it passes through the set of resonances and decays exponentially beyond. The range of $`k_x`$ for which the transfer function is constant is independent of the total number of layers and depends only on the thickness of the individual layers which sets the coupling strength for the plasmon states at the interfaces. In the presence of absorption in the negative medium, however, the decay is extremely fast for the system with larger $`N`$ simply as a consequence of the larger amount of absorptive medium present. Also note that the absorption removes all the divergences in the transmission. As noted by us in earlier publications, the absorption is actually vital in this system to prevent the resonant divergences which would otherwise create artifacts that dominate the image.

Next we keep the total length of the stack fixed and change the number of layers. In the lossless case, the range of $`k_x`$ for which there is effective amplification of the evanescent waves, simply increases with reducing layer thickness as can be seen in figure 7. Of course, the number of transmission resonances which depend on the number of surface states increases with the number of layers. With absorptive material, however, the transmission decays faster with $`k_x`$ for larger $`k_x`$ in the case of the thicker slabs (10nm) than in the case of the thinner slabs (5nm). This reconfirms our analytical result that the effects of absorption would be less deleterious for the image resolution in the case of thinner layers. Note that the total amount of absorptive material in this case is the same in both the cases.

In any case, the absorption in the negative dielectric (metal) appears to set the ultimate limit on the image resolution in this case of the layered medium. We have noted earlier in Ref.  that the effects of absorption could be minimised by using a large dielectric constant, GaAs say ($`\epsilon _+=12`$), for the positive medium and tuning to the appropriate frequency where the perfect lens condition $`\epsilon _{}^{^{}}=\epsilon _+`$ is satisfied for the real part of the dielectric constant $`\epsilon _{}^{^{}}`$ of the metal. In the case of silver, the imaginary part of the permittivity or the absorption is reasonably constant ($``$0.4) over the frequency range of interest. Hence, it is immediately seen that the fractional deviation from the perfect lens condition in the imaginary part is smaller when the real part of the permittivity is large and hence the amplification of the evanescent waves becomes more effective. Now we show the transmission obtained across a multilayer stack where $`\epsilon _+=12`$ and $`\epsilon _{}=12+i0.4`$, corresponding to alternating slabs of silver and GaAs, in figure 8. We must first note that the wavelength of light at which the perfect lens condition for the permittivity of silver is satisfied is different in the two cases. Using the empirical formula for the dispersion of silver, we obtain $`\epsilon _{}^{^{}}=1`$ at 356 nm and $`\epsilon _{}^{^{}}=12`$ at 578 nm. In figure 8 for the lossless system, the transmission resonances appear to occur at higher values of $`k_x/k_0`$ for the high index system, but it must realised that $`k_0=2\pi /\lambda `$ is smaller in this case and the corresponding image resolution would actually be lower. However, when we compare the transmission with absorption included, the beneficial effects of using the larger value of the dielectric constant become obvious. The transmission coefficient indeed decays much more slowly with $`k_x`$ in this case. Also note that we have taken the source to be in air and the image to be formed inside the high-index dielectric medium.

Finally, we show in figure 9 the images of two slits of 15nm width and a peak-to-peak separation of 45 nm obtained by using a single slab of silver as the lens and a layered medium of alternating layers of silver and a positive dielectric medium as the lens. The total distance from the object plane to the image plane in both cases is $`2d=`$ 80 nm. The images of the slits in the case of the single slab lens are hardly resolved, whereas the images of the slits are well separated and clearly resolved in the case of the layered lens. The enhancement in the image resolution for the layered lens is obvious from the figure. The bump seen in between the slits is an artifact due to the fact that the transmission function is not exactly a constant for all wave-vectors.

## IV Conclusions

We have elaborated the design of the perfect lens by considering a multilayer stack and shown that this has advantages over the original configuration of a single slab of material. In particular the effects of absorption are much reduced by the division into mutilayers. The limiting case of infinitesimal multilayers was also considered and shown to be equivalent to an effective medium through which the image propagates without distortion as if it were conveyed by an array of very fine infinitely conducting wires. We went on to make a detailed analysis of how imperfections in the lens affects the image quality. The effects of retardation and the coupled slab plasmon resonances can be minimized by considering very thin layers of 5 to 10 nm thickness. The effects of absorption then dominate the image transfer, but are less deleterious when the individual layer thicknesses are smaller. The effects of absorption can also be minimized by using materials with higher dielectric constants, and tuning the frequency of the radiation to meet the perfect lens conditions.

## Acknowledgments

SAR would like to acknowledge the support from DoD/ONR MURI grant N00014-01-1-0803.


text:
# References

UA/NPPS-4-2002

CONSERVATION OF PARTICLE MULTIPLICITIES BETWEEN CHEMICAL AND THERMAL FREEZE-OUT
A. S. Kapoyannis
University of Athens, Division of Nuclear and Particle Physics,
GR-15771 Athens, Greece

## Abstract

The evolution of a hadronic system after its chemical decomposition is described through a model that conserves the hadronic multiplicities to their values at chemical freeze-out. In the partition function describing the model all known hadronic resonances with masses up to 2400 MeV have been included. The state of the system is found as function of temperature and the corresponding baryon density is evaluated. The baryon density at thermal decoupling is also computed.

PACS numbers: 25.75.-q, 12.40.Ee, 05.70.Ce, 12.38.Mh

Keywords: chemical, thermal, freeze-out, hadron gas, baryon density

1. Introduction

Thermal approaches have extensively been used to describe the particle multiplicities which emerge from high energy collisions \[1-16\]. The results of such approaches are satisfactory since they are able to predict quite accurately a large number of different experimentally measured hadronic multiplicities as function of a few thermodynamic variables, such as temperature, volume and chemical potentials.

The extracted parameters from such approaches define the “chemical freeze-out”, i.e. the point where the chemical composition of the system that produces the particle multiplicities is fixed. Another point called “thermal freeze-out” can also be defined. This second point is associated with the particle momentum distribution which is measured experimentally. After this point this distribution remains fixed and the particles no longer interact among themselves.

There is evidence that in a lot of circumstances the two points corresponding to the same system do not occur at the same temperature. Generally the freeze-out temperature is lower than the chemical freeze-out one. Since the particles are measured once they have reached the experimental apparatus, after any kind of interaction among themselves has ceased and since the thermodynamic parameters of the chemical freeze-out point predict quite well all those particles, one has to infer that all these abundances will have to remain fixed through the whole process between chemical and thermal freeze-out.

In this paper the main focus will be to construct a model of relativistic hadronic particles formulated in the grand canonical ensemble that will be able to conserve the particle multiplicities between points which correspond to different temperatures. Of course the two points with greatest interest are the chemical and the freeze-out point. The question that arises then is why there is need for another model and why not use one of the existing thermal models. In these models the particles are described as thermally equilibrated entities. In some of these models the hadrons are non-interacting particles \[2-7\] and in others a kind of interaction among them has been included \[1,8-16\]. But usually the free thermodynamic parameters are the volume $`V`$, the temperature $`T`$ and a set of a few chemical potentials each of which is associated with the conservation of quantum numbers like baryon number $`B`$, charge $`Q`$, strangeness $`S`$, etc. During the evolution of the system among states with different temperatures its content can alter. This is done by adjusting the chemical potentials appropriately so as to keep the relevant quantum numbers fixed. The only parameter that remains unfixed after the conservation of quantum numbers is the volume. If there is need to fix the particle numbers as well it is clear that this cannot be accomplished with the existing parameters. Even if it is assumed that conservation of all the particle numbers automatically conserves the quantum numbers, meaning that only the conservation of the particles is enough, this limits the application of the above models to only situations where the number of the particle entities is equal to the number of the chemical potentials plus one111 In such a case the quantum number potentials and the volume would have to be adjusted to keep the particles fixed between points with different temperatures.. But generally the hadrons which have to be considered are more numerous than the quantum numbers potentials.

As it is evident these thermal models cannot evolve the system from chemical to thermal freeze-out. In this work the necessity to have fixed particle numbers will be used to construct a model which will determine the evolution of the hadronic system after its chemical freeze-out. This newly constructed model will coincide at chemical freeze-out with “Ideal Hadron Gas” model (IHG) \[2-5\], one of the aforementioned thermal models, which is formulated in the grand canonical ensemble and describes hadrons as relativistic non-interacting particles.

2. Formulation of the model

Before discussing the new model, the IHG model will be presented. In the context of IHG the grand canonical partition function, formulated in the Boltzmann approximation, has the form

$$\mathrm{ln}Z(V,T,\{\lambda \})_{IHG}=V\underset{i}{}\lambda _{QN_i}\underset{j}{}Z_{H_{ij}}(T)V\underset{i}{}\lambda _{QN_i}\underset{j}{}\frac{T}{2\pi ^2}g_{ij}m_{ij}^2K_2(\frac{m_{ij}}{T}),$$
(1)
where $`i`$ runs over all hadronic families such as mesons, $`N`$ Baryons, $`\mathrm{\Lambda }`$ Baryons, etc. and $`j`$ represents the specific member of the family with degeneracy factor $`g_{ij}`$ and mass $`m_{ij}`$. $`\lambda _{QN_i}`$ stands for the product of all the fugacities associated with the particular family. These fugacities can either be quantum numbers fugacities related to Baryon number, Strangeness, etc. or to quark flavour222For example, for $`\mathrm{\Xi }^{}`$ Baryons, $`\lambda _{QN}`$ would read $`\lambda _B\lambda _Q^1\lambda _S^2\gamma _s^2`$ or $`\lambda _d\lambda _s^2\gamma _s^2`$. One can look in , eq. 14, to find out how the two sets of fugacities are related..

One can evaluate particle abundances if the above partition function is extended by the introduction of a fugacity $`\lambda _{ij}`$ for every particle. After calculating the particle number one has to set in IHG $`\lambda _{ij}=1`$ , so again the particle number is only expressed as function of the quantum numbers fugacities.

Now, if someone wishes to keep the particle numbers fixed it is only natural to extend (1) by the use of fugacities $`\lambda _{ij}`$ (corresponding to particle numbers), but with the difference that the constraint $`\lambda _{ij}=1`$ will not be imposed. This model will be called Fixed Particle Numbers (FPN) model and accordingly its partition function will depend on $`\lambda _{ij}`$’s

$$\mathrm{ln}Z(V,T,\{\lambda \})_{FPN}=V\underset{ij}{}\lambda _{H_{ij}}Z_{H_{ij}}(T)V\underset{ij}{}\lambda _{H_{ij}}\frac{T}{2\pi ^2}g_{ij}m_{ij}^2K_2(\frac{m_{ij}}{T}),$$
(2)
where $`\lambda _{H_{ij}}`$ is product of quantum numbers as well as particle number fugacities333For example, for $`\mathrm{\Xi }(1530)^{}`$ Baryon, $`\lambda _H`$ would read $`\lambda _B\lambda _Q^1\lambda _S^2\gamma _s^2\lambda _{\mathrm{\Xi }(1530)^{}}`$.. Since the experimentally measured multiplicities usually contain feeding from the decay of resonances, all known hadrons with masses up to 2400 MeV have been included in the FPN partition function. The same hadrons have also been included to IHG partition function (1). The mean particle number can be evaluated through the relation

$$<N_{ij}>_{FPN}=\lambda _{ij}\frac{\mathrm{ln}Z(V,T,\{\lambda \})_{FPN}}{\lambda _{ij}}|_{\{\lambda \}\lambda _{ij}},$$
(3)
where $`\{\lambda \}\lambda _{ij}`$ means that for the evaluation of the partial derivative all fugacities except $`\lambda _{ij}`$ are considered as constants.

The next point that has to be elucidated is at what values the particle numbers will stay fixed. Since after chemical freeze-out these values do not alter it is useful to formulate FPN so as to keep the particle numbers fixed at their chemical freeze-out values. One has to remember that the chemical freeze-out values are extracted from a thermal model, like IHG, by a successful fit to the experimentally measured values. The particle numbers can then be fixed at the values calculated through IHG for the chemical freeze-out thermodynamic variables. Thus the results of IHG and FPN should coincide at chemical freeze-out point. So it has to be required444The primed variables in this paper will generally be related to subsequent points of the chemical freeze-out point.

$$<N_{ij}>_{IHG}=<N_{ij}^{}>_{FPN}V\lambda _{QN_i}Z_{H_{ij}}(T)=V^{}\lambda _{H_{ij}}Z_{H_{ij}}(T^{})$$

$$\lambda _{H_{ij}}=\frac{V\lambda _{QN_i}Z_{H_{ij}}(T)}{V^{}Z_{H_{ij}}(T^{})}.$$
(4)
The above equation can be used to calculate the total product of fugacities $`\lambda _{H_{ij}}`$ consisting of quantum numbers fugacities and of hadron fugacities at temperature $`T^{}`$. It has to be pointed out that it is not possible to evaluate each quantum number fugacity separately, but this is irrelevant since the full product of fugacities can be calculated. Another focal point is that all quantum numbers are automatically conserved as linear combination of the particle numbers.

In the right hand side of eq. (4) the only term which is left undetermined after the imposition of the conservation of particle numbers is the multiplicand factor $`\frac{V}{V^{}}`$. So an additional constraint has to be applied. For example conservation of entropy can be assumed555A lot of authors assume isentropic evolution of the system, e.g. see ..

The entropy of the system can in general be calculated from666The symbol of entropy is tilded in order not to be confused with the symbol of Strangeness. $`K`$ can be set equal to one.

$$\stackrel{~}{S}=\left(\frac{[T\mathrm{ln}Z(V,T,\{\mu \})]}{T}\right)_{V,\{\mu \}},$$
(5)
where $`\mu `$ represents the chemical potential associated with fugacity $`\lambda =\mathrm{exp}(\mu /T)`$. Applying (5) to the partition function (2)777Since the IHG and FPN partition functions coincide at chemical freeze-out point, the IHG partition function (1), where $`\lambda _{ij}=1`$, can be used for the evaluation of the entropy at this point. the constraint of fixed entropy will read

$$\stackrel{~}{S}=\stackrel{~}{S}^{}$$

$$\mathrm{ln}Z_{IHG}(V,T,\{\mu \})+VT\underset{ij}{}\lambda _{QN_i}\frac{Z_{H_{ij}}(T)}{T}VT\underset{ij}{}\lambda _{QN_i}\frac{\mu _{QN_i}}{T^2}Z_{H_{ij}}(T)=$$

$$=\mathrm{ln}Z_{FPN}(V^{},T^{},\{\mu ^{}\})+V^{}T^{}\underset{ij}{}\lambda _{H_{ij}}\frac{Z_{H_{ij}}(T^{})}{T^{}}V^{}T^{}\underset{ij}{}\lambda _{H_{ij}}\frac{\mu _{H_{ij}}}{T^2}Z_{H_{ij}}(T^{}).$$

With the use of (4) the last equation becomes

$$VT\underset{ij}{}\lambda _{QN_i}\frac{Z_{H_{ij}}(T)}{T}V\underset{ij}{}\lambda _{QN_i}\mathrm{ln}(\lambda _{QN_i})Z_{H_{ij}}(T)=$$

$$=V^{}T^{}\underset{ij}{}\frac{V\lambda _{QN_i}Z_{H_{ij}}(T)}{V^{}Z_{H_{ij}}(T^{})}\frac{Z_{H_{ij}}(T^{})}{T^{}}V^{}\underset{ij}{}Z_{H_{ij}}(T^{})\frac{V\lambda _{QN_i}Z_{H_{ij}}(T)}{V^{}Z_{H_{ij}}(T^{})}\mathrm{ln}(\frac{V\lambda _{QN_i}Z_{H_{ij}}(T)}{V^{}Z_{H_{ij}}(T^{})})$$

$$T\underset{ij}{}\lambda _{QN_i}\frac{Z_{H_{ij}}(T)}{T}\underset{ij}{}\lambda _{QN_i}\mathrm{ln}(\lambda _{QN_i})Z_{H_{ij}}(T)=$$

$$=T^{}\underset{ij}{}\lambda _{QN_i}\frac{Z_{H_{ij}}(T)}{Z_{H_{ij}}(T^{})}\frac{Z_{H_{ij}}(T^{})}{T^{}}\underset{ij}{}\lambda _{QN_i}Z_{H_{ij}}(T)\mathrm{ln}(\frac{V\lambda _{QN_i}Z_{H_{ij}}(T)}{V^{}Z_{H_{ij}}(T^{})}).$$
(6)
Setting $`x\frac{V^{}}{V}`$, (6) can be solved for $`x`$ to give

$$x_{FPN}=$$

$$\mathrm{exp}\left[\frac{_{ij}\lambda _{QN_i}Z_{H_{ij}}(T)\mathrm{ln}(\frac{Z_{H_{ij}}(T)}{Z_{H_{ij}}(T^{})})+T_{ij}\lambda _{QN_i}\frac{Z_{H_{ij}}(T)}{T}T^{}_{ij}\lambda _{QN_i}\frac{Z_{H_{ij}}(T)}{Z_{H_{ij}}(T^{})}\frac{Z_{H_{ij}}(T^{})}{T^{}}}{_{ij}\lambda _{QN_i}Z_{H_{ij}}(T)}\right].$$
(7)
Equation (7) can be used to evaluate the volume expansion ratio as the system has cooled to a temperature $`T^{}`$ less than the chemical freeze-out temperature $`T`$. With the use of the same equation, quantities like the baryon density of the system can be calculated at $`T^{}`$. The baryon chemical potential at which the system is found at temperature $`T^{}`$ cannot be calculated separately from the rest of chemical potentials in the context of FPN. But the baryon density has no problem to be evaluated. One has to remember that baryon number is also fixed with the imposition of the constraints (4). So

$$n_{BFPN}=\frac{}{V^{}}=\frac{^{ch}}{V^{}}=\frac{V}{V^{}}\frac{^{ch}}{V}=\frac{n_B^{ch}}{x_{FPN}},$$
(8)
where $`n_B^{ch}`$ is the baryon density calculated at chemical freeze-out.

3. Application

The newly constructed model, FPN, can describe thermally equilibrated hadronic systems with fixed particle numbers when their temperature is known. As an example the systems formed at different interactions at SPS will be considered. For this reason the chemical freeze-out parameters obtained for these systems through fits to their experimentally measured values will be used. These parameters are listed in Table 1 along with the references where they can be found. From a variety of thermal analyses performed by different authors the particular ones have been chosen because they allow for partial strangeness equilibrium ($`\gamma _s1`$) and they use most recent available values for the experimentally measured hadronic multiplicities. The values of Table 1 are then taken, for each interaction separately, as input to the equations $`<S>=0`$ and $`\frac{<B>}{2<Q>}=\beta `$,888$`\beta `$ is fixed from the baryon number and charge of the participant nucleons, e.g. see . to determine the rest of the fugacities. Thus the whole set of chemical freeze-out parameters $`(T,\mu _B,\mu _Q,\mu _S,\gamma _s)`$ are calculated and also the products of fugacities $`\lambda _{QN_i}`$ in (1) are also set.

Giving different values to temperature $`T`$, equation (8) can be used to calculate the corresponding baryon density. The resulting paths for FPN for $`S+S`$, $`S+Ag`$ and $`Pb+Pb`$ interactions are shown in Figure 1 with solid curves. For the $`Pb+Pb`$ interaction the thermal freeze-out temperature is calculated in Refs.  and . For these values baryon density at thermal freeze-out $`n_B^{ther}`$ can be evaluated. The results are listed in the last column of Table 2. The path for $`Pb+Pb`$ is followed until the lower temperature (of the two given in Refs. ) is reached. The points which correspond to the thermal freeze-out temperatures of these references are depicted with squares on the FPN curve.

FPN has the unique attribute to conserve each particle species separately. So it is not possible to compare its results directly with another thermal model. In order to have a general view we shall depict on the graphs with the FPN results IHG states for different temperatures. The IHG model to be used only conserves the quantum numbers $`<B>`$, $`<Q>`$ and $`<S>`$ and the entropy $`<\stackrel{~}{S}>`$. The IHG states are represented by a dotted curve in Figure 1. It has to be pointed out that IHG does not take the system from chemical to thermal freeze-out. The system remains in the context of IHG all the time at chemical equilibrium where the particles can transform into one another.

The usefulness of Figure 1 is that, as can be observed, a “soup” of fixed particles cools more efficiently than a chemically equilibrated IHG state. This means that at the same baryon density (which through the conservation of the baryon number is equivalent to equal volumes) smaller temperature corresponds to the FPN state than to the IHG state. The reason is that IHG only conserves a few quantum numbers. As temperature drops the number of particles diminishes but without affecting the preservation of quantum numbers. For example an equal reduction to the number of protons and antiprotons will not affect the conservation of baryon number. In the FPN case, on the contrary, no particle number is allowed to diminish. So in the IHG state the available energy has to be distributed among less particles than FPN, and their mean kinetic energy has to be greater, leading to greater temperature.

In Figure 2 the FPN ratios $`x=V^{}/V`$, where $`V`$ is the chemical freeze-out volume, are plotted as function of temperature for the three SPS interactions. For comparison the ratios $`x=V_{IHG}^{}/V_{IHG}`$ for the particular IHG model discussed above are plotted with dotted curves. The two volumes in the IHG ratios correspond to chemically equilibrated states at different temperatures. $`V_{IHG}^{}`$ and $`V_{IHG}`$ are calculated for the same temperatures as $`V^{}`$ and $`V`$ respectively.

Also, in order to show the effect of hadrons with large masses, FPN calculations with hadrons only up to the Delta mass have been included in figures 1 and 2. The chemical freeze-out parameters for $`Pb+Pb`$ of Table 1 have been used for these calculations. From figure 1 it is evident that the baryon density calculated through the truncated version of FPN is considerably less than the baryon density calculated at the same temperature with the FPN model using all the hadrons. The expansion ratio $`x`$ for the truncated FPN, though, is close to the calculated ratio with FPN using all the hadrons, as it can be seen from figure 2.

4. Conclusion

After chemical freeze-out the collisions among hadrons that compose the hadronic gas can no longer change its chemical composition. Following this requirement a non-interacting hadron gas model (FPN) has been presented that keeps the multiplicity of every particle fixed to the value dictated by the chemical freeze-out conditions. In the context of FPN the constraints of conservation of quantum numbers are broken up to a larger number of constraints, these of conservation of particle numbers. The chemical potentials of quantum numbers are no longer “good” variables to describe the evolution of the system. Of course the fugacities of particle numbers used as variables in FPN are not “free” parameters. Their values are fixed from the given set of the quantum numbers fugacities at chemical freeze-out. So the evolution of a hadronic system is described as function of temperature (after imposing conservation of entropy). This is done for three SPS interactions.

Following this evolution and using values of thermal freeze-out temperature extracted for the $`Pb+Pb`$ interaction the baryon density at freeze-out is evaluated. As the temperature at thermal decoupling for various interactions can be calculated using transverse mass spectra or HBT analysis  the same procedure can be applied to evaluate the corresponding baryon density before free streaming for these interactions.

In this paper the IHG has been used as the thermal model which would coincide with FPN at chemical freeze-out point. Any other thermal model, interacting or non-interacting, can also be used in place of IHG and, with the use of particle fugacities, a model that conserves particle multiplicities can also be formed. Recently in model II of 999 appeared during the refereeing period of the present paper. the particle fugacities have been used to conserve particle multiplicities. This a three dimensional hydrodynamic model describing radial and elliptic flow but it includes hadrons with masses only up to Delta mass and it is restricted to zero baryon chemical potential.

Acknowledgement I would like to thank Professor N. G. Antoniou and Professor C. N. Ktorides for reviewing the manuscript and for useful remarks.

Figure Captions

###### Figure 1

Contours (solid lines) that follow hadronic systems after chemical freeze-out on ($`T,n_B`$) plane for 3 interactions at SPS, calculated through FPN (model of Fixed Particle Numbers). On the same graph points at chemical equilibrium (dotted lines) calculated through an IHG model that conserves baryon number, charge, strangeness and entropy are depicted. The slashed curve represents calculations with FPN including hadrons with masses only up to the Delta mass (FPN (a)) and using as chemical freeze-out parameters the ones that correspond to $`Pb+Pb`$ interaction of Table 1.

###### Figure 2

The ratios of the volumes $`V^{}`$ of the SPS hadronic systems at a certain temperature to their volumes $`V^{ch}`$ at chemical freeze-out calculated for FPN (solid lines). The dotted lines correspond to the ratios of volumes of chemically equilibrated states of IHG model of Figure 1 calculated at the same temperatures as the volumes of the FPN ratios. The slashed curve represents calculations with model FPN (a) of Figure 1.

Table Captions

###### Table 1

Chemical freeze-out parameters calculated for different interactions at SPS and the corresponding references.

###### Table 2

Thermal freeze-out temperature calculated in two different references for the $`Pb+Pb`$ interaction and the corresponding computation of baryon density through FPN. The upper errors of baryon density correspond to the upper errors of temperature. The same is true for the lower errors.

| Experiment | $`T^{ch}(MeV)`$ | $`\mu _B^{ch}(MeV)`$ | $`\gamma _s^{ch}`$ | Reference |
| --- | --- | --- | --- | --- |
| S+S 200 $`AGeV`$ | $`180.5\pm 10.9`$ | $`220.2\pm 18.0`$ | $`0.747\pm 0.048`$ |  |
| S+Ag 200 $`AGeV`$ | $`178.9\pm 8.1`$ | $`241.5\pm 14.5`$ | $`0.711\pm 0.063`$ |  |
| Pb+Pb 158 $`AGeV`$ | $`174.7\pm 6.7`$ | $`240\pm 14`$ | $`0.900\pm 0.049`$ |  |

Table 1.

| Experiment | $`T^{ther}(MeV)`$ | Reference | $`n_B^{ther}`$ ($`fm^3`$) |
| --- | --- | --- | --- |
| Pb+Pb 158 $`AGeV`$ | $`120\pm 12`$ |  | $`0.099_{0.019}^{+0.022}`$ |
| Pb+Pb 158 $`AGeV`$ | $`95.8\pm 3.5`$ |  | $`0.0627_{0.0045}^{+0.0047}`$ |

Table 2.