text:
# References

A set of orthogonal polynomials, dual to

alternative $`q`$-Charlier polynomials

Natig M. Atakishiyev

Instituto de Matemáticas, UNAM, CP 62210 Cuernavaca, Morelos, México

and

Anatoliy U. Klimyk

Bogolyubov Institute for Theoretical Physics, 03143 Kiev, Ukraine

E-mail: natig@matcuer.unam.mx and aklimyk@bitp.kiev.ua

## Abstract

The aim of this paper is to derive (by using two operators, representable by a Jacobi matrix) a family of $`q`$-orthogonal polynomials, which turn to be dual to alternative $`q`$-Charlier polynomials. A discrete orthogonality relation and a three-term recurrence relation for these dual polynomials are explicitly obtained. The completeness property of dual alternative $`q`$-Charlier polynomials is also established.

Mathematical Subject Classification (2000): 33D45, 47B36, 81Q10

1. Introduction

It is well known that orthogonal polynomials are closely connected with spectral properties of symmetric operators, which can be represented in some basis by a Jacobi matrix, and with the classical moment problem. Namely, a spectrum of an operator, represented by a Jacobi matrix, is determined by an orthogonality measure for corresponding orthogonal polynomials. If orthogonal polynomials admit many orthogonality relations, then the corresponding symmetric operator is not self-adjoint and it leaves room for infinitely many self-adjoint extensions. These extensions are determined by orthogonality measures for the appropriate orthogonal polynomials.

Contrary to orthogonal polynomials of the hypergeometric type (Wilson, Jacobi, Laguerre and so on), basic hypergeometric polynomials (or $`q`$-orthogonal polynomials) are not so deeply understood yet. These polynomials are collected in the $`q`$-analogue of the Askey–scheme of orthogonal polynomials (see, for example, ), which starts with Askey–Wilson polynomials and $`q`$-Racah polynomials, introduced in  and . Importance of these polynomials is magnified by the fact that they are closely related to the theory of quantum groups. As an instance of such connection we refer to a paper , in which Al-Salam–Chihara $`q`$-orthogonal polynomials have been employed to construct locally compact quantum group $`SU_q(1,1)`$. Another application of $`q`$-orthogonal polynomials is related to the theory of $`q`$-difference equations, which often surface in contemporary theoretical and mathematical physics.

The purpose of present paper is to study completeness and duality properties of alternative $`q`$-Charlier polynomials. We shall show below that this originates a novel type of $`q`$-orthogonal polynomials and a discrete orthogonality relation for them. To achieve this, we essentially use two operators $`I_1`$ and $`q^{J_0}`$, which are certain representation operators for the quantum algebra $`U_q(\mathrm{su}_{1,1})`$ with a lowest weight (however, we do not use explicitly the theory of representations in what follows). The operator $`I_1`$ is related to the three-term recurrence relation for alternative $`q`$-Charlier polynomials. We diagonalize the trace class operator $`I_1`$ and obtain two bases in the Hilbert space: an initial basis and a basis of normalized eigenvectors of $`I_1`$. These bases are connected by an orthogonal matrix. The orthogonality relations for rows and columns of this matrix lead to orthogonality relations for alternative $`q`$-Charlier polynomials and for the functions, which are dual to these polynomials. We extract from the latter functions a dual set of polynomials and obtain a discrete orthogonality relation for them. As a result, one is led to the completeness property of dual alternative $`q`$-Charlier polynomials.

Observe that the present paper is a continuation of our research, initiated in  and .

Throughout the sequel we always assume that $`q`$ is a fixed positive number such that $`q<1`$. We use (without additional explanation) notations of the theory of special functions and the standard $`q`$-analysis (see, for example,  or ).

2. Pair of operators $`(I_1,J)`$

Let $``$ be a separable complex Hilbert space with an orthonormal basis $`|n`$, $`n=0,1,2,\mathrm{}`$. We define on $``$ two operators. The first one, denoted as $`q^{J_0}`$, acts on the basis elements as $`q^{J_0}|n=q^n|n`$, and the second one, denoted as $`I_1`$, is given by the formula

$$I_1|n=a_n|n+1+a_{n1}|n1+b_n|n,$$
$`(1)`$
$$a_n=(aq^{3n+1})^{1/2}\frac{\sqrt{(1q^{n+1})(1+aq^n)}}{(1+aq^{2n+1})\sqrt{(1+aq^{2n})(1+aq^{2n+2})}},$$

$$b_n=q^n\left(\frac{1+aq^n}{(1+aq^{2n})(1+aq^{2n+1})}+aq^{n1}\frac{1q^n}{(1+aq^{2n1})(1+aq^{2n})}\right),$$

where $`a`$ is any fixed positive number. Clearly, $`I_1`$ is a symmetric operator.

Since $`a_n0`$ and $`b_n0`$ when $`n\mathrm{}`$, the operator $`I_1`$ is bounded. Therefore, we assume that it is defined on the whole Hilbert space $``$. For this reason, $`I_1`$ is a self-adjoint operator. Let us show that $`I_1`$ is a trace class operator. For the coefficients $`a_n`$ and $`b_n`$ from (1), we have $`a_{n+1}/a_nq^{3/2}`$ and $`b_{n+1}/b_nq`$ when $`n\mathrm{}`$. Since $`0<q<1`$, for the sum of all matrix elements of the operator $`I_1`$ in the basis $`|n`$, $`n=0,1,2,\mathrm{}`$, we have $`_n(2a_n+b_n)<\mathrm{}`$. This means that $`I_1`$ is a trace class operator. Thus, a spectrum of $`I_1`$ is discrete and has a single accumulation point at 0. Moreover, a spectrum of $`I_1`$ is simple, since $`I_1`$ is representable by a Jacobi matrix with $`a_n0`$ (see , Chapter VII).

To find eigenfunctions $`\xi _\lambda `$ of the operator $`I_1`$, $`I_1\xi _\lambda =\lambda \xi _\lambda `$, we set $`\xi _\lambda =_n\beta _n(\lambda )|n`$, where $`\beta _n(\lambda )`$ are appropriate numerical coefficients. Acting by the operator $`I_1`$ upon both sides of this relation, one derives that $`_{n=0}^{\mathrm{}}\beta _n(\lambda )(a_n|n+1+a_{n1}|n1+b_n|n)=\lambda _{n=0}^{\mathrm{}}\beta _n(\lambda )|n`$, where $`a_n`$ and $`b_n`$ are the same as in (1). Collecting in this identity all factors, which multiply $`|n`$ with fixed $`n`$, one derives the recurrence relation for the coefficients $`\beta _n(\lambda )`$:

$$\beta _{n+1}(\lambda )a_n+\beta _{n1}(\lambda )a_{n1}+\beta _n(\lambda )b_n=\lambda \beta _n(\lambda ).$$

The substitution

$$\beta _n(\lambda )=\left(\frac{(a;q)_n(1+aq^{2n})}{(q;q)_n(1+a)(a/q)^n}\right)^{1/2}q^{n(n+3)/4}\beta _n^{}(\lambda )$$

reduces this relation to the following one

$$A_n\beta _{n+1}^{}(\lambda )C_n\beta _{n1}^{}(\lambda )+(A_n+C_n)\beta _n^{}(\lambda )=\lambda \beta _n^{}(\lambda ),$$

$$A_n=q^n\frac{1+aq^n}{(1+aq^{2n})(1+aq^{2n+1})},C_n=aq^{2n1}\frac{1q^n}{(1+aq^{2n1})(1+aq^{2n})}.$$

This is the recurrence relation for the alternative $`q`$-Charlier polynomials

$$K_n(\lambda ;a;q):={}_{2}{}^{}\varphi _{1}^{}(q^n,aq^n;\mathrm{\hspace{0.33em}0};q,q\lambda )$$

(see, formulas (3.22.1) and (3.22.2) in ). Therefore, $`\beta _n^{}(\lambda )=K_n(\lambda ;a;q)`$ and

$$\beta _n(\lambda )=\left(\frac{(a;q)_n(1+aq^{2n})}{(q;q)_n(1+a)a^n}\right)^{1/2}q^{n(n+1)/4}K_n(\lambda ;a;q).$$
$`(2)`$
For the eigenvectors $`\xi _\lambda `$ we thus have the expression

$$\xi _\lambda =\underset{n=0}{\overset{\mathrm{}}{}}\left(\frac{(a;q)_n(1+aq^{2n})}{(q;q)_n(1+a)a^n}\right)^{1/2}q^{n(n+1)/4}K_n(\lambda ;a;q)|n.$$
$`(3)`$
Since the spectrum of the operator $`I_1`$ is discrete, only for a discrete set of values of $`\lambda `$ these vectors belong to the Hilbert space $``$. This discrete set of eigenvectors determines a spectrum of $`I_1`$.

Now we look for a spectrum of the operator $`I_1`$ and for a set of polynomials, dual to alternative $`q`$-Charlier polynomials. To this end we use the action of the operator

$$J:=q^{J_0}aq^{J_0}$$

upon the eigenvectors $`\xi _\lambda `$, which belong to the Hilbert space $``$. In order to find how this operator acts upon these vectors, one can use the $`q`$-difference equation

$$(q^naq^n)K_n(\lambda )=aK_n(q\lambda )+\lambda ^1K_n(\lambda )\lambda ^1(1\lambda )K_n(q^1\lambda )$$
$`(4)`$
for the alternative $`q`$-Charlier polynomials $`K_n(\lambda )K_n(\lambda ;a;q)`$ (see formula (3.22.5) in ). Multiply both sides of (4) by $`d_n|n`$ and sum up over $`n`$, where $`d_n`$ are the coefficients of $`K_n(\lambda ;a;q)`$ in the expression (2) for $`\beta _n(\lambda )`$. Taking into account the formula (3) and the fact that $`J|n=(q^naq^n)|n`$, one obtains the relation

$$J\xi _\lambda =a\xi _{q\lambda }+\lambda ^1\xi _\lambda \lambda ^1(1\lambda )\xi _{q^1\lambda }.$$
$`(5)`$
We shall see in the next section that the spectrum of the operator $`I_1`$ consists of the points $`q^n`$, $`n=0,1,2,\mathrm{}`$. This means that $`J`$ has the form of a Jacobi matrix in the basis of eigenvectors of $`I_1`$; that is, the pair of the operators $`I_1`$ and $`J`$ form a Leonard pair (see  for the corresponding definition).

3. Spectrum of $`I_1`$ and orthogonality of alternative $`q`$-Charlier polynomials

The aim of this section is to find, by using the operators $`I_1`$ and $`J`$, a basis in the Hilbert space $``$, which consists of eigenvectors of the operator $`I_1`$ in a normalized form, and to derive explicitly the unitary matrix $`U`$, connecting this basis with the basis $`|n`$, $`n=0,1,2,\mathrm{}`$, in $``$. This matrix leads directly to the orthogonality relation for alternative $`q`$-Charlier polynomials. For this purpose we first find a spectrum of $`I_1`$.

Let us analyze a form of the spectrum of $`I_1`$ from the point of view of the spectral theory of trace class operators. If $`\lambda `$ is a spectral point of the operator $`I_1`$, then (as it is easy to see from (5)) a successive action by the operator $`J`$ upon the vector (eigenvector of $`I_1`$) $`\xi _\lambda `$ leads to the eigenvectors $`\xi _{q^m\lambda }`$, $`m=0,\pm 1,\pm 2,\mathrm{}`$. However, since $`I_1`$ is a trace class operator, not all of these points may belong to the spectrum of $`I_1`$, since $`q^m\lambda \mathrm{}`$ when $`m+\mathrm{}`$ if $`\lambda 0`$. This means that the coefficient $`1\lambda ^{}`$ of $`\xi _{q^1\lambda ^{}}`$ in (5) must vanish for some eigenvalue $`\lambda ^{}`$. Clearly, it vanishes when $`\lambda ^{}=1`$. Moreover, this is the only possibility for the coefficient of $`\xi _{q^1\lambda ^{}}`$ in (5) to vanish, that is, the point $`\lambda =1`$ is a spectral point for the operator $`I_1`$. Let us show that the corresponding eigenfunction $`\xi _1\xi _{q^0}`$ belongs to the Hilbert space $``$.

By formula (II.6) of Appendix II in , one has $`K_n(1;a;q)={}_{2}{}^{}\varphi _{1}^{}(q^n,aq^n;\mathrm{\hspace{0.33em}0};q,q)=(a)^nq^{n^2}`$. Therefore,

$$\xi _1,\xi _1=\underset{n=0}{\overset{\mathrm{}}{}}\frac{(a;q)_n(1+aq^{2n})}{(1+a)(q;q)_na^nq^{n(n+1)/2}}K_n^2(1;a;q)=\underset{n=0}{\overset{\mathrm{}}{}}\frac{(a;q)_n(1+aq^{2n})a^n}{(1+q)(q;q)_nq^{n(3n1)/2}}.$$
$`(6)`$
In order to calculate this sum, we take the limit $`d,e\mathrm{}`$ in the equality

$$\underset{n=0}{\overset{\mathrm{}}{}}\frac{(1+aq^{2n})(a;q)_n(d;q)_n(e;q)_n}{(1+a)(aq/d;q)_n(aq/e;q)_n(q;q)_n}\left(\frac{aq}{de}\right)^nq^{n(n1)/2}=\frac{(aq;q)_{\mathrm{}}(aq/de;q)_{\mathrm{}}}{(aq/d;q)_{\mathrm{}}(aq/e;q)_{\mathrm{}}}$$

(see formula in Exercise 2.12, Chapter 2 of ). Since

$$\underset{d,e\mathrm{}}{lim}(d;q)_n(e;q)_n(aq/de)^n=q^{n(n1)}(aq)^n,$$

we obtain from here that the sum in (6) is equal to $`(aq;q)_{\mathrm{}}`$, that is, $`\xi _1,\xi _1<\mathrm{}`$ and $`\xi _1`$ belongs to the Hilbert space $``$. Thus, the point $`\lambda =1`$ does belong to the spectrum of the operator $`I_1`$.

Let us find other spectral points of the operator $`I_1`$ (recall that a spectrum of $`I_1`$ is discrete). Setting $`\lambda =1`$ in (5), we see that the operator $`J`$ transforms $`\xi _{q^0}`$ into a linear combination of the vectors $`\xi _q`$ and $`\xi _{q^0}`$. Moreover, $`\xi _q`$ belongs to the Hilbert space $``$, since the series

$$\xi _q,\xi _q=\underset{n=0}{\overset{\mathrm{}}{}}\frac{(a;q)_n(1+aq^{2n})}{(1+a)(q;q)_na^n}q^{n(n+1)/2}K_n^2(q;a;q)$$

is majorized by the corresponding series (6) for $`\xi _{q^0}`$. Therefore, $`\xi _q`$ belongs to the Hilbert space $``$ and the point $`q`$ is an eigenvalue of the operator $`I_1`$. Similarly, setting $`\lambda =q`$ in (5), we find that $`\xi _{q^2}`$ is an eigenvector of $`I_1`$ and the point $`q^2`$ belongs to the spectrum of $`I_1`$. Repeating this procedure, we find that all $`\xi _{q^n}`$, $`n=0,1,2,\mathrm{}`$, are eigenvectors of $`I_1`$ and the set $`q^n`$, $`n=0,1,2,\mathrm{}`$, belongs to the spectrum of $`I_1`$. So far, we do not know yet whether other spectral points exist or not.

The vectors $`\xi _{q^n}`$, $`n=0,1,2,\mathrm{}`$, are linearly independent elements of the Hilbert space $``$ (since they correspond to different eigenvalues of the self-adjoint operator $`I_1`$). Suppose that values $`q^n`$, $`n=0,1,2,\mathrm{}`$, constitute a whole spectrum of $`I_1`$. Then the set of vectors $`\xi _{q^n}`$, $`n=0,1,2,\mathrm{}`$, is a basis in the Hilbert space $``$. Introducing the notation $`\mathrm{\Xi }_k:=\xi _{q^k}`$, $`k=0,1,2,\mathrm{}`$, we find from (5) that

$$J\mathrm{\Xi }_k=a\mathrm{\Xi }_{k+1}+q^k\mathrm{\Xi }_kq^k(1q^k)\mathrm{\Xi }_{k1}.$$

As we see, the matrix of the operator $`J`$ in the basis $`\mathrm{\Xi }_k`$, $`k=0,1,2,\mathrm{}`$, is not symmetric, although in the initial basis $`|n`$, $`n=0,1,2,\mathrm{}`$, it was symmetric. The reason is that the matrix $`(a_{mn})`$ with entries $`a_{mn}:=\beta _m(q^n)`$, $`m,n=0,1,2,\mathrm{}`$, where $`\beta _m(q^n)`$ are the coefficients (2) in the expansion $`\xi _{q^n}=_m\beta _m(q^n)|n`$, is not unitary. This fact is equivalent to the statement that the basis $`\mathrm{\Xi }_n=\xi _{q^n}`$, $`n=0,1,2,\mathrm{}`$, is not normalized. To normalize it, one has to multiply $`\mathrm{\Xi }_n`$ by corresponding numbers $`c_n`$ (which are not known at this moment). Let $`\widehat{\mathrm{\Xi }}_n=c_n\mathrm{\Xi }_n`$, $`n=0,1,2,\mathrm{}`$, be a normalized basis. Then the matrix of the operator $`J`$ is symmetric in this basis. Since $`J`$ has in the basis $`\{\widehat{\mathrm{\Xi }}_n\}`$ the form

$$J\widehat{\mathrm{\Xi }}_n=c_{n+1}^1c_na\widehat{\mathrm{\Xi }}_{n+1}+q^n\widehat{\mathrm{\Xi }}_nc_{n1}^1c_nq^n(1q^n)\widehat{\mathrm{\Xi }}_{n1},$$

then its symmetricity means that $`c_{n+1}^1c_na=c_n^1c_{n+1}q^{n1}(1q^{n+1})`$, that is, $`c_n/c_{n1}=\sqrt{aq^n/(1q^n)}`$. Therefore,

$$c_n=c(a^nq^{n(n+1)/2}/(q;q)_n)^{1/2},$$

where $`c`$ is a constant.

The expansions

$$\widehat{\xi }_{q^n}(x)\widehat{\mathrm{\Xi }}_n(x)=\underset{m}{}c_n\beta _m(q^n)|m\underset{m}{}\widehat{a}_{mn}|m$$
$`(7)`$
connect two orthonormal bases in the Hilbert space $`_l`$. This means that the matrix $`(\widehat{a}_{mn})`$, $`m,n=0,1,2,\mathrm{}`$, with entries

$$\widehat{a}_{mn}=c_n\beta _m(q^n)=c\left(\frac{a^nq^{n(n+1)/2}}{(q;q)_n}\frac{(a;q)_m(1+aq^{2m})}{(1+a)(q;q)_ma^mq^{m(m+1)/2}}\right)^{1/2}K_m(q^n;a;q)$$
$`(8)`$
is unitary, provided that the constant $`c`$ is appropriately chosen. In order to calculate this constant, we use the relation $`_{m=0}^{\mathrm{}}|\widehat{a}_{mn}|^2=1`$ for $`n=0`$. Then this sum is a multiple of the sum in (6) and, consequently, $`c=(aq;q)_{\mathrm{}}^{1/2}`$.

The matrix $`(\widehat{a}_{mn})`$ is real and orthogonal, that is,

$$_n\widehat{a}_{mn}\widehat{a}_{m^{}n}=\delta _{mm^{}},_m\widehat{a}_{mn}\widehat{a}_{mn^{}}=\delta _{nn^{}}.$$
$`(9)`$
Substituting into the first sum over $`n`$ in (9) the expressions for $`\widehat{a}_{mn}`$, we obtain the identity

$$\underset{n=0}{\overset{\mathrm{}}{}}\frac{a^nq^{n(n+1)/2}}{(q;q)_n}K_m(q^n;a;q)K_m^{}(q^n;a;q)=\frac{(aq^m;q)_{\mathrm{}}a^m(q;q)_m}{(1+aq^{2m})}q^{m(m+1)/2}\delta _{mm^{}},$$
$`(10)`$
which must yield the orthogonality relation for alternative $`q`$-Charlier polynomials. An only gap, which remains to be clarified, is the following. We have assumed that the points $`q^n`$, $`n=0,1,2,\mathrm{}`$, exhaust the whole spectrum of $`I_1`$. Let us show that this is the case.

Recall that the self-adjoint operator $`I_1`$ is represented by a Jacobi matrix in the basis $`|n`$, $`n=0,1,2,\mathrm{}`$. According to the theory of operators of such type (see, for example, , Chapter VII), eigenvectors $`\xi _\lambda `$ of $`I_1`$ are expanded into series in the basis $`|n`$, $`n=0,1,2,\mathrm{}`$, with coefficients, which are polynomials in $`\lambda `$. These polynomials are orthogonal with respect to some positive measure $`d\mu (\lambda )`$ (moreover, for self-adjoint operators this measure is unique). The set (a subset of $``$), on which these polynomials are orthogonal, coincides with the spectrum of the operator under consideration and the spectrum is simple.

We have found that the spectrum of $`I_1`$ contains the points $`q^n`$, $`n=0,1,2,\mathrm{}`$. If the operator $`I_1`$ would have other spectral points $`x`$, then on the left-hand side of (10) there would be other summands $`\mu _{x_k}K_m(x_k;a;q)K_m^{}(x_k;a;q)`$, corresponding to these additional points. Let us show that these additional summands do not appear. To this end we set $`m=m^{}=0`$ in the relation (10) with the additional summands. Since $`K_0(x;a;q)=1`$, we have the equality

$$\underset{n=0}{\overset{\mathrm{}}{}}\frac{a^nq^{n(n+1)/2}}{(q;q)_n}+\underset{k}{}\mu _{x_k}=(aq;q)_{\mathrm{}}.$$

According to the formula for the $`q`$-exponential function $`E_q(a)`$ (see formula (II.2) of Appendix II in ), we have $`_{n=0}^{\mathrm{}}a^nq^{n(n+1)/2}(q;q)_n^1=(aq;q)_{\mathrm{}}`$. Hence, $`_k\mu _{x_k}=0`$ and all $`\mu _{x_k}`$ disappear. This means that additional summands do not appear in (10) and it does represent the orthogonality relation for alternative $`q`$-Charlier polynomials.

As we have shown, the orthogonality relation for the alternative $`q`$-Charlier polynomials is given by formula (10). Due to this orthogonality, we arrive at the following statement:

Proposition. The spectrum of the operator $`I_1`$ coincides with the set of points $`q^n`$, $`n=0,1,2,\mathrm{}`$. The spectrum is simple and has one accumulation point at 0.

4. Dual alternative $`q`$-Charlier polynomials

Now we consider the second identity in (9), which gives the orthogonality relation for the matrix elements $`\widehat{a}_{mn}`$, considered as functions of $`m`$. Up to multiplicative factors these functions coincide with the functions

$$F_n(x;a|q)={}_{2}{}^{}\varphi _{1}^{}(x,a/x;\mathrm{\hspace{0.33em}0};q,q^{n+1}),$$
$`(11)`$
considered on the set $`x\{q^m|m=0,1,2,\mathrm{}\}`$. Consequently,

$$\widehat{a}_{mn}=\left(\frac{a^nq^{n(n+1)/2}}{(q;q)_n}\frac{(1+aq^{2m})}{(aq^m;q)_{\mathrm{}}(q;q)_ma^mq^{m(m+1)/2}}\right)^{1/2}F_n(q^m;a|q)$$

and the second identity in (9) gives the orthogonality relation for $`F_n(q^m;a|q)`$:

$$\underset{m=0}{\overset{\mathrm{}}{}}\frac{(1+aq^{2m})}{a^m(aq^m;q)_{\mathrm{}}(q;q)_mq^{m(m+1)/2}}F_n(q^m;a|q)F_n^{}(q^m;a|q)=\frac{(q;q)_n}{a^nq^{n(n+1)/2}}\delta _{nn^{}}.$$
$`(12)`$
The functions $`F_n(x;a,b|q)`$ can be represented in another form. Indeed, taking in the relation (III.8) of Appendix III in  the limit $`c\mathrm{}`$, one derives the relation

$${}_{2}{}^{}\varphi _{1}^{}(q^m,aq^m;\mathrm{\hspace{0.33em}0};q,q^{n+1})=(a)^mq^{m^2}{}_{3}{}^{}\varphi _{0}^{}(q^m,aq^m,q^n;q,q^n/a).$$

Therefore, we have

$$F_n(q^m;a|q)=(a)^mq^{m^2}{}_{3}{}^{}\varphi _{0}^{}(q^m,aq^m,q^n;;q,q^n/a).$$
$`(13)`$
The basic hypergeometric function $`{}_{3}{}^{}\varphi _{0}^{}`$ in (13) is a polynomial of degree $`n`$ in the variable $`\mu (m):=q^maq^m`$, which represents a $`q`$-quadratic lattice; we denote it by

$$d_n(\mu (m);a;q):={}_{3}{}^{}\varphi _{0}^{}(q^m,aq^m,q^n;;q,q^n/a).$$
$`(14)`$
Then formula (12) yields the orthogonality relation

$$\underset{m=0}{\overset{\mathrm{}}{}}\frac{(1+aq^{2m})a^m}{(aq^m;q)_{\mathrm{}}(q;q)_m}q^{m(3m1)/2}d_n(\mu (m))d_n^{}(\mu (m))=\frac{(q;q)_n}{a^nq^{n(n+1)/2}}\delta _{nn^{}}$$
$`(15)`$
for the polynomials (14) when $`a>0`$. As far as we know this orthogonality relation is new. We call the polynomials $`d_n(\mu (m);a;q)`$ dual alternative $`q`$-Charlier polynomials. Thus, we proved the following theorem.

Theorem. The polynomials $`d_n(\mu (m);a;q)`$, given by formula (14), are orthogonal on the set of points $`\mu (m):=q^maq^m`$, $`m=0,1,2,\mathrm{}`$, and the orthogonality relation is given by formula (15).

Remark. The duality of polynomials is the well-known notion (see  and ) and, in particular, in the case of polynomials, orthogonal with respect to a finite number of discrete points, it reflects the simple fact that a finite-dimensional matrix, orthogonal by rows, is also orthogonal by its columns (cf (9)). There is also an analytical way of deriving a dual orthogonality for polynomials, whose weight functions are supported on an infinite number of discrete points (see, for example, -). But this derivation is given in terms of a dual set of functions (for instance, their explicit form (11) for the case of alternative $`q`$-Charlier polynomials is given above) and one needs to make one step further in order to extract an appropriate family of dual polynomials from these functions. Observe that in our approach to the duality of $`q`$-polynomials it is not assumed that an initial family is orthogonal, this orthogonality property is straightforwardly derived. Besides, one naturally extracts from a dual set of orthogonal functions an appropriate dual family of $`q`$-polynomials.

Let $`𝔩^2`$ be the Hilbert space of functions on the set $`m=0,1,2,\mathrm{}`$ with the scalar product

$$f_1,f_2=\underset{m=0}{\overset{\mathrm{}}{}}\frac{(1+aq^{2m})a^m}{(aq^m;q)_{\mathrm{}}(q;q)_m}q^{m(3m1)/2}f_1(m)\overline{f_2(m)},$$
$`(16)`$
where the weight function is taken from (15). The polynomials (14) are in one-to-one correspondence with the columns of the orthogonal matrix $`(\widehat{a}_{mn})`$ and the orthogonality relation (15) is equivalent to the orthogonality of these columns. Due to (9) the columns of the matrix $`(\widehat{a}_{mn})`$ form an orthonormal basis in the Hilbert space of sequences $`𝐚=\{a_n|n=0,1,2,\mathrm{}\}`$ with the scalar product $`𝐚,𝐚^{}=_na_n\overline{a_n^{}}`$. This scalar product is equivalent to the scalar product (16) for the polynomials $`d_n(\mu (m);a;q)`$. For this reason, the set of polynomials $`d_n(\mu (m);a;q)`$, $`n=0,1,2,\mathrm{}`$, form an orthogonal basis in the Hilbert space $`𝔩^2`$. This means that the point measure in (15) is extremal for the dual alternative $`q`$-Charlier polynomials $`d_n(\mu (m);a;q)`$ (for the definition of an extremal orthogonality measure see, for example, in ). This means that the dual alternative $`q`$-Charlier polynomials (14) form a complete system in the $`L^2`$-space with respect to the point measure in (15). Observe that the completeness of alternative $`q`$-Charlier polynomials in the $`L^2`$-space with respect to the point measure in (10) is a consequence of the fact that the operator (1) is bounded.

A recurrence relation for the polynomials $`d_n(\mu (m);a;q)`$ is derived from (4). It has the form

$$(q^maq^m)d_n(\mu (m))=ad_{n+1}(\mu (m))+q^nd_n(\mu (m))q^n(1q^n)d_{n1}(\mu (m)),$$
$`(17)`$
where $`d_n(\mu (m))d_n(\mu (m);a;q)`$. A $`q`$-difference equation for $`d_n(\mu (m);a;q)`$ can be obtained from the three-term recurrence relation for alternative $`q`$-Charlier polynomials.

Note that for the polynomials $`d_n(\mu (m);a;q^1)`$ with $`q<1`$ we have the expression

$$d_n(\mu (m);a;q^1)={}_{3}{}^{}\varphi _{2}^{}(q^m,aq^m,q^n;\mathrm{\hspace{0.33em}0},0;q,q).$$
$`(18)`$
However, the recurrence relation for these polynomials (which can be obtained from the relation (17), does not satisfy the positivity condition $`A_nC_{n+1}>0`$, that is, they are not orthogonal polynomials for $`a>0`$ (as it is the case for alternative $`q`$-Charlier polynomials). This positivity condition holds only if we demand $`a<0`$. In this case, the polynomials (18) are the continuous big $`q`$-Hermite polynomials $`H_n(x;a|q)`$ (for an explicit form of these polynomials see, for example, , formula (3.18.1)), which are orthogonal on a certain continuous set of points.

Acknowledgments

This research has been supported in part by the SEP-CONACYT project 41051-F and the DGAPA-UNAM project IN112300 ”Óptica Matemática”. A. U. Klimyk acknowledges the Consejo Nacional de Ciencia y Technología (México) for a Cátedra Patrimonial Nivel II.


format: text

text:
# Nonperturbative scalar-meson resonances with open charm and beauty

## 1 Introduction

On studying the vast amount of literature  dealing with the phenomena observed in $`S`$-wave meson-meson scattering, one may wonder whether the question “What is a resonance?” has any answer at all. Different theoretical considerations even come to opposite conclusions  for the same phenomenon.

One reason for the difficulty to unambiguously explain meson-meson-scattering phenomena, may stem from the following. In atomic physics, the perturbative approach towards electromagnetic (EM) transitions is justified, since the EM binding forces are apparently more decisive for the energy-level splittings than photon radiation. We may thus predict the radiation frequencies for atomic transitions to a good accuracy, by just determining the atomic levels of the binding forces. Atomic states have lifetimes which are large as compared with the internal EM frequencies. Hence, electrons complete many revolutions about the nucleus before the system decays. Such systems can be considered stable, like planetary ones. We may, therefore, determine the masses of excited atoms, without making a drastic approximation, just from the EM bound-state spectrum, which procedure is only correct for stable systems. The wave functions for very unstable systems have to be determined by different methods. Most probably, they can be expressed as linear combinations of the bound-state solutions, in which case not just one, but all excitations contribute to the transition amplitude.

In Table 1 we express the stability of quark-antiquark systems by estimating the number of complete cycles of the system before decay. For the internal frequency of $`q\overline{q}`$ systems, we take the universal frequency $`\omega `$ of the Nijmegen model (see Table 2).

We read from Table 1 that pions and Kaons are very stable, whereas the eta mesons are reasonably stable. However, strongly decaying mesons are very unstable.

For example, we are convinced that the structure observed at energies below 1 GeV in $`P`$-wave $`\pi \pi `$ scattering, is caused by a quark-antiquark state, the $`\rho `$ meson, with a mass of about 770 MeV. But if the $`\rho `$ resonance stems from a state described by a quark confined to an antiquark, then that state is far from stable. With a size of approximately 1 fm, its internal frequency must be some 200 MeV. However, the lifetime of the $`q\overline{q}`$ state is of the same order of magnitude, indicated by its width of 150 MeV. Consequently, the $`q\overline{q}`$ system associated with the $`\rho `$ resonance barely completes ONE complete revolution during its existence. We may thus wonder whether it makes sense to refer to such a phenomenon as a state, or even a particle, while we just concluded that it does not really exist.

Resonances in scattering are somehow supposed to stem from an underlying spectrum of bound states which are broken up by a weaker interaction. Such a philosophy could in principle work very well for the stable mesons $`\pi `$ and $`K`$, since their formation is associated with strong interactions, whereas their decays are due to electro-weak processes. Nevertheless, strong decay can have dramatic effects even for the masses of pions and Kaons . Once these effects are fully included in the internal dynamics of pions and Kaons, we may indeed handle electro-weak processes perturbatively .

In atomic systems the mechanism for binding is the same as the mechanism for transitions, namely electromagnetism. Nevertheless, since the fine-structure constant $`\alpha `$ is small, and transitions are of higher order in $`\alpha `$ than the binding forces, one may still approximate EM transitions of atoms by the Born term in a perturbative approach. However, the strong binding force of quarks is of comparable order of magnitude as the strong forces that cause hadronic decay. Consequently, in order to analyse hadronic scattering data, one must treat both aspects of strong interactions at least on the same footing, maybe even attributing more importance to strong decay.

The consequences of such a combined approach are, for the time being, still to be studied in effective models, which, however primitive they are and no matter what details need be addressed in the future, have already indicated two maior properties of hadronic resonances:

1. Upon assuming an underlying spectrum of bare mesonic states with $`q\overline{q}`$ degrees of freedom only, resonances in meson-meson scattering come out at energies which are very different from the energy levels of the underlying confinement spectrum. This has been confirmed by he Cornell group, focusing on heavy mesonic resonances of the charmonium and beautonium systems , the Nijmegen group, describing, in principle, all mesonic resonances , and the Lisbon group, in the framework of a chiral quark model . All report mass shifts, due to coupled channels, of up to several hundred MeV, which is comparable in magnitude to the level splittings of hadronic resonances.

2. In experiment, more resonances may be observed than just those stemming from the assumed underlying confinement spectrum. The Nijmegen Model predicts several mesonic resonances  and even rather stable states  which do have a different origin than the majority of mesonic resonances.

## 2 Resonances

Although resonances have not yet been well defined by us, we turn in this section to the issue of their description. Let us first discuss atomic transitions. When a photon passes through a container filled with a gas of non-excited atoms, then it can be absorbed once its frequency matches some allowed transition from the ground state to an excited state of the atom’s specrum. Ideally, the absorption of photons is measured as a bell-shaped fluctuation of the light beam’s intensity, when the light frequency is continuously varied. The reason why a Breit-Wigner (BW) form rather accurately fits the shape of this intensity fluctuation is that the Born term of the interaction dominates the process, and moreover, the excited state is relatively stable. Measured over a larger frequency interval, the absorption curve for photons can be described by a sum of BW-type contributions, each one corresponding to a specific excited state of the atom and its respective life time.

Were the fine-structure constant $`\alpha `$ not small, then not only higher-order terms in $`\alpha `$ would significantly contribute, but also would the photon noticeably couple to all higher excitations. Hence, the form of the scattering amplitude would then certainly not coincide with BW approximations about each peak of the absorption curve. As a consequence, summing up BW-type amplitudes is not justified for strong interactions. This can easily be demonstrated with a simple model.

In Fig. 1, we display the cross section for harmonic-oscillator confinement coupled to meson-meson scattering. The parameters, i.e., the nonstrange and strange quark masses and the oscillator frequency, are chosen (see Table 2) such that the

radial excitations come out at 1.39, 1.77, 2.15, 2.53, … GeV. In the (elastic) meson-meson scattering channel, we have chosen the meson masses at 0.14 GeV (pion) and 0.5 GeV (Kaon). The intensity of the coupling, parametrised by $`\lambda `$, is varied. For $`\lambda =0.1`$ we are in the small-coupling limit (upper-left picture of Fig. 1), for which a sum of BW amplitudes works reasonably well. We also notice that, for $`\lambda =0.1`$, all peaks are near the energy levels of the pure harmonic oscillator.

But as we increase the coupling, we observe that the cross section has structures which are not of a BW type. For $`\lambda =0.4`$, we see a dominant tructure at 0.8 GeV, which in a more realistic approach  can be identified with the $`K_0^{}(800)`$ resonance. From the results of Fig. 1, we conclude that, for the description of strong interactions, models are needed treating the binding forces and the transition interactions on an equal footing. Our proposal for elastic meson-meson scattering is to use the following expression for the partial-wave phase shift:

$$\text{cotg}\left(\delta _{\mathrm{}}(p)\right)=\frac{\pi \lambda ^2\mu p{\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{𝒥_n\mathrm{}^{}(p)𝒩_n\mathrm{}(p)}{E(p)E_{n\mathrm{}_c}}}\mathrm{\hspace{0.33em}1}}{\pi \lambda ^2\mu p{\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{𝒥_n\mathrm{}^{}(p)𝒥_n\mathrm{}(p)}{E(p)E_{n\mathrm{}_c}}}}.$$
(1)
The details of this formula and its derivation can be found in Refs. . It expresses the fact that in strong hadronic decay all possible confinement states contribute to the scattering amplitude, through the factors $`𝒥_n\mathrm{}(p)`$ and $`𝒩_n\mathrm{}(p)`$, which are convolution integrals over all possible confinement excitations, the transition potential, and spherical Bessel resp. Neumann functions. In the limit of small coupling $`\lambda `$, Eq. (1) reduces to a BW form, with a dominant contribution from the nearest confinement state.

At this stage we may formulate our definition of resonances. Structures like those depicted in Fig. 1 are related to complex-energy singularities in Eq. (1). Those for which a pole in the complex-energy plane (in the nearest Riemann sheet) can be identified are given by us the status of resonance. For small coupling (see Fig. 1 for $`\lambda =0.1`$), Eq. (1) has poles near the confinement spectrum, close to the real axis and with small negative imaginary parts. However, poles can be found, too, with very large (up to several GeVs) imaginary parts, though possibly having no measurable effect. The latter poles can be attributed the status of mathematical resonance, but they may be unobservable. Hence, whether or not we include them in our resonance list depends on the accuracy of measurement, not on their actual existence. For stronger coupling we see them showing up as genuine structures in the $`S`$-wave cross section. Suddenly our spectrum has become richer than the original confinement spectrum. This is the realm of strong interactions!

Expression (1), which is the exact solution of a simple dynamical model, is very compact and elegant. But the summations over the radial excitations of the confined system do not converge fast. Hence, for the numerical calculation of scattering quantities it is not suitable. A more practical method has been described in Ref. . Nevertheless, because of its close contact with the physical situation, Eq. (1) is extremely transparent, and useful as a starting point for making approximations.

For example, when we are studying a limited range of energies, we may restrict the summation to only the few terms near by in energy, characterised by $`E_{n\mathrm{}_c}`$ in Eq. (1), and approximate the remaining terms by a simple function of $`p`$, the relative meson-meson linear momentum. Notice, however, that the resulting expression is only of the type of a BW expansion if the strength of the transition, characterised by $`\lambda `$, is small. For larger values of $`\lambda `$, formula (1) behaves very differently from BW expansions.

Another in practice very useful approximation is inspired by the radial dependence of the transition potential, which was derived in Ref.  for meson-meson scattering in the cases $`J^P=0^{}`$, $`J^P=1^{}`$, and $`J^P=0^+`$. We find there that it is peaked at short distances. Moreover, if we also take into account that confinement wave functions must be of short range, we may just define a transition radius $`a`$ and approximate the spatial integrals in Eq. (1) by choosing a spherical delta shell $`\delta (ra)`$ for the transition potential. Then we end up with the expression

$$\text{cotg}\left(\delta _{\mathrm{}}(p)\right)\frac{2a^4\lambda ^2\mu pj_{\mathrm{}}(pa)n_{\mathrm{}}(pa){\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\left|_{n\mathrm{}_c}(a)\right|^2}{E(p)E_{n\mathrm{}_c}}}\mathrm{\hspace{0.33em}1}}{2a^4\lambda ^2\mu pj_{\mathrm{}}^2(pa){\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\left|_{n\mathrm{}_c}(a)\right|^2}{E(p)E_{n\mathrm{}_c}}}}.$$
(2)
This way we have pulled the $`p`$ dependence outside the infinite summation over the radial confinement spectrum, which makes it much easier to handle truncations, as now the rest term does not depend on $`p`$ and may thus be chosen constant.

## 3 Open charm

Our recent detailed analyses of scalar mesons have actually been performed in the approximation of formula (2). For $`DK`$ elastic $`S`$-wave scattering, we adopt the same expression that has been used for $`K\pi `$ as described in Ref. . We only have to substitute the confinement masses of $`d\overline{s}`$ by those for $`c\overline{s}`$, thereby using the parameters of Table 2. For the ground state of the $`J^P=0^+`$ $`c\overline{s}`$ confinement spectrum, the mass comes out at 2.545 GeV. The first radial excitation is then 380 MeV (the model’s level spacing) higher in mass. The summation over the higher radial excitations in Eq. (2) we approximate by a constant, normalised to 1. The relative couplings of the ground state and the first radial excitation to the $`DK`$ channel are 1.0 and 0.2, respectively . Thus we obtain

$$\underset{n=0}{\overset{\mathrm{}}{}}\frac{\left|_{n,1}(a)\right|^2}{E(p)E_{n,1}}\frac{1.0}{E(p)E_0}+\frac{0.2}{E(p)E_1}\mathrm{\hspace{0.33em}1}\text{GeV}^{\mathrm{\hspace{0.25em}2}},$$
(3)
where $`E_0=2.545`$ and $`E_1=2.925`$. For Eq. (2) with the substitution (3), and furthermore for the physical value of $`\lambda `$ (meaning the value that fits $`K\pi `$ elastic $`S`$-wave scattering from threshold up to 1.6 GeV ), we show the resulting elastic $`S`$-wave $`DK`$ cross section in Ref. . The corresponding pole positions in the complex-energy plane for the $`D_{s0}^{}`$ ground state and first radial excitation are shown in Ref. . These are found at 2.28 GeV and ($`2.78i\mathrm{\hspace{0.17em}0.093}`$) GeV, respectively . We find no structure in the cross section in the region around 2.55 GeV, where typical quark models predict $`J^P=0^+`$ $`c\overline{s}`$ states , nor poles in the scattering amplitude. On the contrary, employing confinement plus unitarisation, we find structures in the cross section some 250 MeV below and 300–400 MeV above that energy, which explains the observation of the (likely) $`D_{s0}^{}`$ ground state below the $`DK`$ threshold, and predicts for the first radial excitation a resonance somewhere between 2.8 and 3.0 GeV, with a width of roughly 100–200 MeV. We may also predict the scattering length for $`DK`$ elastic $`S`$-wave scattering. Keeping in mind that the $`J^P=0^+`$ $`c\overline{s}`$ ground state comes out a little bit too low in mass with our value of $`\lambda `$, we estimate the theoretical error by varying $`\lambda `$. For a conservative choice, we obtain for the scattering length in $`DK`$ $`S`$ wave $`a_{DK}^0=5\pm 1`$ GeV-1.

## 4 Open beauty

For $`J^P=0^+`$ systems containing one $`b`$ quark, we may repeat the above-described procedure. We determine the ground-state masses of the uncoupled systems from the model parameters of Table 2. This results in 5.605, 5.707, and 6.761 GeV for $`\overline{b}u/d`$, $`b\overline{s}`$, and $`\overline{b}c`$, respectively. The first radial excitations lie 0.38 GeV higher. The remaining parameters $`a`$ and $`\lambda `$ are kept the same as before, but now scaled with the reduced constituent $`q\overline{q}`$ mass $`\mu _{qq}`$, in order to guarantee flavour invariance of the strong interactions , i.e.,

$$a_{xy}\sqrt{\mu _{xy}}=\text{constant}\text{and}\lambda _{xy}\sqrt{\mu _{xy}}=\text{constant},$$
(4)
where $`x`$ and $`y`$ represent the two flavours involved . The constants in Eq. (4) are fixed by  $`a_{us}=3.2`$ GeV-1 and $`\lambda _{us}=0.75`$ GeV-3/2 for $`S`$-wave $`K\pi `$ scattering. Henceforth, we shall quote for the parameters $`\lambda `$ and $`a`$ their values scaled down to elastic $`K\pi `$ scattering. The actual values of $`\lambda `$ and $`a`$ can, in each case, be obtained through Eq. (4) and Table 2.

In Fig. 2 we show how in our model the scattering poles move through the complex-energy plane for $`B\pi `$ elastic scattering, when the coupling is varied. For the lowest pole we obtain the very interesting result that the physical pole (for $`\lambda =0.75`$) ends up on top of the threshold. Within the accuracy of our model, this means that one has to expect in experiment either a virtual, or a real bound state, just below threshold. Only through an accurate measurement of the scattering length in the $`B\pi `$ $`S`$ wave, the precise position can be obtained experimentally.

For higher energies, we only find poles at energies far above 5.6 GeV. This implies that we do not expect any resonance in the $`S`$ wave below roughly 5.9 GeV, where the confinement spectrum of $`\overline{b}n`$ has its first radial excitation. In our model we obtain two resonances, each some 200 MeV wide, and about 150 MeV apart.

In Fig. 3 we depict similar pole trajectories, but now for $`BK`$ and $`BD`$ elastic scattering.

We find a $`B_{s0}^{}`$(5570) bound state below the $`BK`$ threshold, which experimentally should show up as a narrow $`\overline{b}u/d`$ resonance. In the $`BD`$ case we find a $`B_{c0}^{}`$(6490) bound state below threshold, which represents a $`b\overline{s}`$ bound state. Also here, nothing is found at 5.707 or 6.761 GeV, where the confinement spectrum has its ground states. Possible resonances above threshold and their widths can be read off from Fig. 3. We summarise these results in Table 3.

## 5 Who is who?

On applying the flavour-independence relations (4) also to $`S`$-wave $`DK`$ scattering, we obtained an interesting result, which we shall report here.

In our initial analysis , we identified the $`D_{sJ}^{}(2317)^+`$ meson with the lowest-lying extra pole, which stems from the background, whereas the higher resonance at about 2.9 GeV, with a width of some 150 MeV, was linked to the confinement ground state. However, if we choose for the delta-shell radius $`a`$ its corrected flavour-independent value, then we find that the roles of these two poles get interchanged.

This phenomenon is depicted in Fig. 4, where, upon a very small variation in $`a`$, the pole trajectories interchange their physically relevant limiting points, which nevertheless change only very little themselves.

As a consequence, we must conclude that the “identity” of a pole cannot always be inferred with certainty. While in the small-coupling situation (see e.g. the $`\lambda =0.1`$ case in Fig. 1) it is unquestionable which poles are to be associated with the confinement spectrum, for intermediate and large couplings one can always find other parameters in Eq. (1) that, when varied simultaneously with the coupling, lead to an identity change. Since intermediate-size couplings are typical for strong hadronic decay, one must conclude that such processes are inherently nonperturbative, and thus must be treated accordingly.

We should mention, moreover, that the flavour-independence relation (4) mildly alters our previous results. The $`D_{sJ}^{}(2317)^+`$ pole now comes out at 2.32 GeV (actually very close to the experimental mass), instead of at the published  value of 2.28 GeV. Furthermore, we find that the higher resonance shifts even more considerably. Here we obtain a width of some 400 MeV instead of 150 MeV, and a central mass around 2.8 GeV, roughly 100 MeV below the initially proposed value.

## 6 Summary and conclusions

The main issue of this paper is schematically depicted in Fig. 5, where we show how the various sectors of meson physics contribute to resonances in meson-meson scattering.

Glue is probably dominant for the dynamics of meson-meson scattering. It provides confinement, contributes dominantly to the effective masses of the light quarks, and is also involved, though in a complicated and not entirely understood fashion, in the creation of quark-antiquark pairs with vacuum quantum numbers. The degrees of freedom of mesons are given by their $`q\overline{q}`$ content. Higher-order configurations of $`q\overline{q}`$ pairs mediate the transitions of a mesonic state to multi-meson states. Figure 5 also indicates schematically the origin of deformation of the confinement sector by the multi-meson sector, which results in mass shifts, resonance widths, extra resonances, and bound states. For all that, processes involving the creation of one or more $`q\overline{q}`$ pairs are generally strong, and should not be handled perturbatively.

As to the spectroscopy of open-charm and open-beauty mesons, we predict several bound states and resonances in $`DK`$, $`B\pi `$, $`BK`$, and $`BD`$ $`S`$ waves.

Acknowledgements

We wish to thank the organisers of the Workshop, in particular Lídia Ferreira, for the warm hospitality. We are also indebted to David Bugg and Frieder Kleefeld for very useful discussions. This work was partly supported by the Fundação para a Ciência e a Tecnologia of the Ministério da Ciência e do Ensino Superior of Portugal, under contract numbers CERN/FIS/43697/2001 and POCTI/FNU/49555/2002.


text:
# Persistence of mean-field features in the energy spectrum of small arrays of Bose-Einstein condensates

## 1 Introduction

In a recent publication , the present authors reported an exhaustive study of the mean-field (MF) dynamics of a system of three interacting Bose-Einstein condensates. Such a system, henceforth referred to as trimer, deserves attention since it can be considered a paradigm for larger arrays. Indeed, similar to longer chains and more in general to larger arrays, its MF dynamics is non-integrable and its phase space is characterized by extended regions of instability. On the other hand, the system is still sufficiently small to allow a thorough analytical study, which provides a valuable operational tool in both interpreting experimental observations and identifying significant experimental setups triggering macroscopic effects . Furthermore, such exhaustive investigation makes the trimer a candidate benchmark for verifying the effectiveness of the MF description, which, being essentially equivalent to the Bogolubov approximation, is expected to provide a reliable picture of the dynamics of an array of interacting BECs for sufficiently large boson populations.

As to the experimental realization of interacting BEC arrays, optical fragmentation allowed to attain 1D , 2D  and even 3D  lattices consisting of a rather large number of sites. Though Bose-Einstein condensation has been already achieved both in an optomagnetic  and, more recently, in an entirely magnetic  double-well trap, microtrap technology  is probably the most promising approach for realizing the small systems addressed by the present article, as well as any complex array characterized by a trapping potential lacking the typical periodicity of optical lattices.

We will focus on two small systems, the asymmetric dimer, i.e. a condensate trapped into an asymmetric double well trap, and the open trimer of reference , where an energy offset between the central and the lateral wells of the trapping potential is present (a multiwell trapping potential characterized by uneven well depths may possibly be achieved within a supercell of a suitable optical superlattice ). Our aim is to compare the spectrum of the quantum Hamiltonians describing these systems with the characteristic energies of the relevant mean-field picture, in the spirit of reference . There Aubry and co-workers illustrated a remarkable relation between the energy of the mean-field saddle point (separatrix) of the symmetric dimer and some features of the spectrum of the quantum counterpart of the same system.

In the following we will extend these results to the more general case of the asymmetric dimer. Subsequently, exploiting the exact mapping between the mean-field picture of the dimer and that of the integrable sub regime of the trimer, which we show to apply at a qualitative level in the purely quantum treatment of the systems, we will point out some more interesting manifestations of MF features in the spectrum of the trimer.

Before proceeding with our analysis, we briefly recall some generalities. The essential physics of an array of coupled BECs is captured by the Bose-Hubbard model , describing a gas of interacting bosons hopping across the sites of an ambient lattice. The relevant Hamiltonian reads $`H=_j[Un_j^2v_jn_jT/2_{kj}(a_j^+a_k+a_ka_j^+)]`$, where $`a_j^+`$ ($`a_j`$) is the creation (annihilation) operator on site $`j`$, $`n_j=a_j^+a_j`$ counts the number of bosons at the same on site, $`U`$ is the on-site interaction parameter, $`v_j`$ is a the energy offset of site $`j`$, $`T`$ is the hopping parameter and the subscript $`kj`$ restricts the sum to the nearest neighbours of site $`j`$.

A mean-field treatment of the Bose-Hubbard model can be achieved by means of a time dependent variational principle based on Glauber’s coherent states . According to such picture each site in the array is described by a macroscopic complex variable $`z_j=\sqrt{\rho _j}\mathrm{}^{\mathrm{}\varphi _j}`$ whose square modulus $`\rho _j`$ and argument $`\varphi _j`$ represent the population and phase associated to the condensate sitting at that site, respectively. Barring constant terms, the effective MF Hamiltonian has the form $`=_j[U|z_j|^4v_j|z_j|^2T/2_{kj}(z_j^{}z_k+z_jz_k^{})]`$ and, equipped with the Poisson brackets among the (canonical) macroscopic variables, $`\{z_j,z_k^{}\}=\delta _{jk}`$ determines the dynamics of the system. 111The set of equations governing the dynamics of the system are basically the same discrete nonlinear Schr̈odinger equations (DNLEs) deriving from the application of the tight-binding approximation to the Gross-Pitaevski equation describing the evolution of a fragmented BEC . Reference  reports on some analytical solutions of similar DNLSEs on small arrays such as the ones presently under concern.. The essential features of the dynamics are determined by the phase space structure, which can be analyzed working out the location and the stability character of the corresponding fixed points. Notice that since $`H`$ features only two-boson operators, it commutes with the the total number operator $`n=_jn_j`$. Likewise, the total population $`\rho =_j\rho _j`$ is a constant of the motion mirroring the global phase invariance of $``$. Hence the fixed point equations have the form $`\{z_j,\mu \rho \}=0`$, where $`\mu `$ is a Lagrange multiplier constraining the total number of bosons. Though the solutions of the fixed point equations are in principle complex vectors, they can be safely restricted to real values, with the only exception of the topological solutions  which are to be found in cyclically invariant structures (i.e. homogeneous closed chains and their higher dimensional generalizations). Since the structures we will consider are neither cyclic nor homogeneous, in solving the fixed point equations we will set $`z_j=\sqrt{\rho }x_j`$, where $`x_j`$ and $`_jx_j^2=1`$.

## 2 The asymmetric dimer

The Bose-Hubbard Hamiltonian for an asymmetric dimer has the form

$$H_\mathrm{d}=U_\mathrm{d}\left(n_{\mathrm{}}^2+n_r^2\right)\frac{w_\mathrm{d}}{2}(n_{\mathrm{}}n_r)\frac{T_\mathrm{d}}{2}\left(a_{\mathrm{}}a_r^++a_{\mathrm{}}^+a_r\right)$$
(1)
where the subscripts $`\mathrm{}`$ and $`r`$ denote the left and the right sites of the structure, respectively, and $`w_\mathrm{d}`$ is a local energy offset deriving from the asymmetry of the double-well trap confining the BEC. As we briefly recalled above, the MF counterpart of Hamiltonian (1) is $`_\mathrm{d}=U_\mathrm{d}(|z_{\mathrm{}}|^4+|z_r|^4)\frac{1}{2}w_\mathrm{d}(|z_{\mathrm{}}|^2|z_r|^2)\frac{1}{2}T_\mathrm{d}(z_{\mathrm{}}z_r^{}+z_{\mathrm{}}^{}z_r)`$, and determines the evolution of the macroscopic semiclassical variables through the (integrable) equations of motion

$$\mathrm{}\mathrm{}\dot{z}_{\mathrm{}}=\left(2U_\mathrm{d}|z_{\mathrm{}}|^2\frac{w_\mathrm{d}}{2}\right)z_{\mathrm{}}\frac{T_\mathrm{d}}{2}z_r,\mathrm{}\mathrm{}\dot{z}_r=\left(2U_\mathrm{d}|z_r|^2+\frac{w_\mathrm{d}}{2}\right)z_r\frac{T_\mathrm{d}}{2}z_{\mathrm{}}$$
(2)
According to the discussion in the previous section, and after introducing the conveniently rescaled parameters

$$\nu _\mathrm{d}=\frac{w_\mathrm{d}}{U_\mathrm{d}\rho _\mathrm{d}},\tau _\mathrm{d}=\frac{T_\mathrm{d}}{U_\mathrm{d}\rho _\mathrm{d}},\mu _\mathrm{d}=\frac{\chi _\mathrm{d}}{U_\mathrm{d}\rho _\mathrm{d}}$$
(3)
where $`\rho _\mathrm{d}=|z_{\mathrm{}}|^2+|z_r|^2=\rho _\mathrm{d}(x_{\mathrm{}}^2+x_r^2)`$, the consequent fixed point equations read

$$\left(2x_{\mathrm{}}^2\frac{\nu _\mathrm{d}}{2}\mu _\mathrm{d}\right)x_{\mathrm{}}\frac{\tau _\mathrm{d}}{2}x_r=0,\left(2x_r^2+\frac{\nu _\mathrm{d}}{2}\mu _\mathrm{d}\right)x_r\frac{\tau _\mathrm{d}}{2}x_{\mathrm{}}=0$$
(4)
Notice that $`\tau _\mathrm{d}`$ and $`\nu _\mathrm{d}`$ are the only significant parameters, whereas $`\mu _\mathrm{d}`$ has an auxiliary function and its value is different for different fixed points. In the following we will consider $`|\nu _\mathrm{d}|`$’s up to the order of $`2`$, corresponding to relative energy offsets of the order of $`2U\rho _\mathrm{d}`$. If large total populations are involved, the limit situation may correspond to an overly asymmetric double well potential. Notice however that comparable relative energies are easily attained in the MF integrable subregime of the trimer, which can be described in terms of an asymmetric dimer, as we will discuss in the next section.

Depending on whether $`|\nu _\mathrm{d}|`$ is smaller or larger than $`(2^{2/3}\tau _\mathrm{d}^{2/3})^{3/2}`$, (4) have either four or two solutions. In the former case there are two maxima, a saddle point and a minimum. As $`|\nu _\mathrm{d}|`$ grows, the local maximum and the saddle get closer and closer, until they collapse and disappear when $`|\nu _\mathrm{d}|`$ crosses the ($`\tau _\mathrm{d}`$-dependent) critical value 222Notice that the coalescence of the (possibly local) maximum and the saddle can be likewise achieved also through a variation of $`\tau _\mathrm{d}`$, and occurs at a $`|\nu _\mathrm{d}|`$-dependent critical value . A similar coupling-driven change in the semiclassical phase space of a 2D BEC array has been recently proposed as a mechanism to trigger an abrupt delocalization in the system .. We mention that the presence of two maxima is an aspect of the mirror-symmetry breaking induced by nonlinearity, leading to interesting phenomena such as self-trapping and $`\pi `$-phase oscillations, that are investigated in reference  for the symmetric dimer.

Aubry and co-workers  focused on the symmetric dimer ($`w_\mathrm{d}=0\nu _\mathrm{d}=0`$), and found that when the MF saddle point (separatrix) is present, i.e. for $`\tau _\mathrm{d}<2`$, its energy is strictly related to some interesting features of the spectrum of $`H_\mathrm{d}`$. They noticed that the density of states of $`H_\mathrm{d}`$ features a sharp peak centered around the saddle point energy. Furthermore this MF energy marks the boundary between two regions of the spectrum characterized by significantly different magnitudes of the splitting energy of the doublets arising from the mirror symmetry of the symmetric dimer. More precisely the energy splitting experiences a drop of several orders of magnitude as the average energy of the doublet becomes larger than the threshold value marked by the saddle point energy. Further analytical considerations about the splitting pattern of the quantum dimer doublets and its relation with the saddle point energy of the corresponding mean-field model can be found in reference . Notice that the magnitude of the splitting between even-odd doublets of the symmetric dimer is also related to the characteristic time-scale for the symmetry-restoring quantum tunneling between semiclassically self-trapped states .

Here we show that the results by Aubry and co-workers can be extended and generalized to the asymmetric dimer, where $`w_\mathrm{d}`$ (and of course $`\nu _\mathrm{d}`$) has a finite value. To this end it proves useful to plot the energies of the MF fixed points of $`_\mathrm{d}`$ versus the energy offset among the two sites of the dimer, $`w_\mathrm{d}`$, which gives a measure of its asymmetry. This can be accomplished recasting the dynamics of the trimer in terms of the alternative set of canonical variables (See e.g. reference  for further details about this alternative canonical representation). $`\rho _\mathrm{d}=|z_1|^2+|z_2|^2`$, $`D=|z_1|^2|z_2|^2`$, $`\theta =(\varphi _1\varphi _2)/2`$, $`\psi =(\varphi _1+\varphi _2)/2`$, whose significant Poisson brackets are $`\{D,\theta \}=\{\rho _d,\psi \}=\mathrm{}^1`$. This way, the relevant dynamic equations can be used to give a parametric representation allowing to plot the fixed point energies versus $`w_\mathrm{d}`$ for a given value of the hopping amplitude. Since it is interesting to compare situations characterized by different total populations, we provide this parametric representation in terms of the rescaled parameters $`\tau _\mathrm{d}`$, $`\nu _\mathrm{d}`$,

$$\nu _\mathrm{d}(d;\tau _\mathrm{d})=2d\left(1+\frac{\tau _\mathrm{d}\mathrm{cos}2\theta }{2\sqrt{1d^2}}\right),E_\mathrm{d}(d;\tau _\mathrm{d})=\frac{(1d^2)^{\frac{3}{2}}\tau _\mathrm{d}\mathrm{cos}2\theta }{2\sqrt{1d^2}},$$
(5)
where $`E_\mathrm{d}`$ is the energy in units of $`U\rho _\mathrm{d}^2`$, $`d=D/\rho _\mathrm{d}[1,1]`$ is the rescaled population imbalance , and $`\theta `$ is either zero or $`\pi /2`$. The former choice pertains to the minimum fixed point no matter for the value of the hopping amplitude. As we discussed above, when $`\tau _\mathrm{d}>2`$ there is only one more point, a maximum, independent of the value of the energy offset $`|\nu _\mathrm{d}|`$. In this case (5) the remaining choice for the angular variable, $`\theta =\pi /2`$, describes the energy of the maximum on the entire domain, $`d[1,1]`$. If, conversely, $`\tau _\mathrm{d}<2`$, the choice $`\theta =\pi /2`$ pertains to different fixed points depending on the domain of evaluation. More in detail, the energy of the saddle, local maximum and absolute maximum are obtained evaluating (5) with $`|d|`$ belonging to the interval $`[0,d_1]`$, $`[d_1,d_2]`$ and $`[d_2,1]`$, respectively, where $`d_1=\sqrt{1(\tau _\mathrm{d}/2)^{2/3}}`$ and $`d_2=\sqrt{1(\tau _\mathrm{d}/2)^2}`$ .

The circles, the thick solid line and the thick dashed lines appearing in both panels of figure 1 were obtained as described above, and show the behaviour of the energy of the MF saddle, minimum and maxima, respectively, for $`\tau _\mathrm{d}=0.1`$. In the left panel of the same figure we also plotted the energy levels of $`H_\mathrm{d}`$ for a total population of $`\rho _\mathrm{d}=16`$ bosons (thin solid lines), using alternating colors for emphasizing the avoided level crossing feature. Notice indeed that at a discrete set of values of the energy offset ($`w_\mathrm{d}=2kU_\mathrm{d}`$, $`k=1,\mathrm{},\rho _\mathrm{d}`$), signalled by the vertical dotted lines, some pairs of adjacent energy levels become so close that, if plotted with the same color, they would give the illusion of crossing as they actually do when $`\tau _\mathrm{d}=0`$. If conversely $`\tau _\mathrm{d}>0`$, the adjacent energy levels never cross333We mention that the avoided level crossing is discussed to some extent in reference . However there the varying parameter is the hopping amplitude between the sites of a dimer and a third site lacking of the repulsive term., but rather they “bounce” on each other, so that their behaviour with varying $`\nu _\mathrm{d}`$ resembles the trajectory of particles repelling each other (see e.g. ). Hence, whenever the energy offset between the wells is an even multiple of $`U_\mathrm{d}`$ ($`\nu _\mathrm{d}=2k/\rho _\mathrm{d}`$), the situation becomes similar to that of the symmetric dimer: the couples of degenerate levels possibly present for $`\tau _\mathrm{d}=0`$ split for finite values of the hopping parameter. Once again, the extent of the splitting of a given doublet can vary of order of magnitudes depending on whether its average energy is smaller or larger than the MF saddle point energy. Hence also in the asymmetric case such energy marks a clear boundary between two regions of the spectrum with quite different features. This is rather evident in the inset of the right panel of figure 1, where we plotted the energy levels for a quantum dimer with a total population of $`\rho _\mathrm{d}=200`$ bosons (black and gray solid lines) in a small region around the saddle point energy (white circles). The main plot of the right panel of figure 1 shows the (Lorentzian broadened) density of states of Hamiltonian (1) for the same total population. It appears evident that also in the asymmetric case the sharp peak of the density of states (dark shading) follows very closely the saddle point energy (white filled circles). These features are observed whenever the hopping amplitude is sufficiently small ($`\tau _\mathrm{d}<2,T_\mathrm{d}<2U_\mathrm{d}\rho _\mathrm{d}`$) to allow the appearance of a saddle point for some values of $`\nu _\mathrm{d}`$.

## 3 The trimer

In reference  the present authors analyzed the mean-field dynamics of a trimer described by the following Bose-Hubbard Hamiltonian

$$H_\mathrm{t}=U_\mathrm{t}\underset{j=1}{\overset{3}{}}n_j^2+\frac{w_\mathrm{t}}{2}\left(n_1n_2+n_3\right)\frac{T_\mathrm{t}}{2}[(a_1+a_3)a_2^++(a_1^++a_3^+)a_2]$$
(6)
where both the relative offset of the central well, $`w_\mathrm{t}`$, and the hopping amplitude $`T_\mathrm{t}`$ were regarded as adjustable parameters. The results of this study were gathered into two stability diagrams summarizing, for any choice of the significant parameters, the character and linear stability of the entire set of MF fixed points, and hence describing the essential features of the phase space structure. This knowledge would allow to produce the trimeric counterparts of figure 1 with no further effort other than the diagonalization of Hamiltonian (6). Notice however that the size of this operator scales as $`(\rho _\mathrm{t}+1)(\rho _\mathrm{t}+2)/2`$, which could be so large to make analytical diagonalization impracticable even for a relatively small total boson population $`\rho _\mathrm{t}`$. Thus the only readily available analytical information to provide a systematic description of the quite intricate spectrum of the trimer (see e.g. figure 3) is the $`Z_2`$ symmetry of (6), allowing to classify its energy levels according to their parity.

This is why here we prefer to proceed discussing a close relation between the fixed point of the trimer and those of the dimer which will turn out to persist also at the quantum level, allowing us to use the results reported in section 2 as a guide to gain deeper insight in the spectrum of the trimer.

### 3.1 The mean-field picture of the trimer and its analogy to the symmetric dimer

The effective Hamiltonian ensuing from (6), $`_\mathrm{t}=U_\mathrm{t}_j|z_j|^4+\frac{1}{2}w_\mathrm{t}(|z_1|^2|z_2|^2+|z_3|^2)\frac{1}{2}T_\mathrm{t}[z_2^{}(z_1+z_3)+z_2(z_1^{}+z_3^{})]`$, yields the following dynamic equations

$$\mathrm{}\mathrm{}\dot{z}_j=\left(2U|z_j|^2+\frac{w_\mathrm{t}}{2}\right)z_j\frac{T_\mathrm{t}}{2}z_2,\mathrm{}\mathrm{}\dot{z}_2=\left(2U|z_2|^2\frac{w_\mathrm{t}}{2}\right)z_2\frac{T_\mathrm{t}}{2}\left(z_1+z_3\right)$$
(7)
where the first equation applies to both of the lateral sites, $`j=1,3`$. This feature is a consequence of the mirror symmetry that $`_\mathrm{t}`$ inherited from $`H_\mathrm{t}`$, and allows room for an integrable sub regime in the generally non-integrable dynamics of the trimer . Notice indeed that the initial condition $`z_1=z_3`$ is propagated in time by (7). As we will discuss in some detail in the following, this sub regime of the trimer can be mapped exactly onto the asymmetric dimer considered in the previous section, which is why it is also referred to as dimeric sub regime of the trimer. However, since the analogy with the dimer applies more in general to the energies of most of the fixed points of the trimer, either exactly or approximately, we will refrain from assigning the adjective “dimeric” exclusively to the integrable sub regime of the trimer and the relevant fixed points. To avoid confusion, the latter will be referred to as symmetric fixed points (SFPs) alluding to their feature $`z_1=z_3`$. Likewise we will term asymmetric fixed point (AFP) any fixed point such that $`z_1z_3`$, with the only exception of the so-called central depleted well  ($`z_1=z_3`$, $`z_2=0`$) which, having no counterpart in the asymmetric dimer, will not be considered in the following discussion.

The close relation between the integrable sub regime of the trimer and the asymmetric dimer can be emphasized by recasting equations (7), subject to the (dynamically propagated) constraint $`z_1=z_3`$, in terms of the alternative variables 444Likewise, if the identification $`(z_{\mathrm{}},z_r)=(z_1/\sqrt{2},z_2)`$ is made, the resulting equations have the same form as those describing a dimer with site-dependent on-site interaction . $`d=(|z_1|^2/2|z_2|^2)/\rho _\mathrm{t}`$ and $`\varphi =(\varphi _1\varphi _2)/2`$.The resulting equations are exactly in the same form as those ensuing describing the dynamics of an asymmetric dimer with the following parameters

$$U_\mathrm{d}=\frac{3}{4}U_\mathrm{t},w_\mathrm{d}=w_\mathrm{t}\frac{1}{2}\rho _\mathrm{t}U_\mathrm{t},T_\mathrm{d}=\sqrt{2}T_\mathrm{t},\rho _\mathrm{d}=\rho _\mathrm{t}$$
(8)
Hence, as we mentioned in section 2, the integrable sub regime of a trimer with no energy offset ($`w_\mathrm{t}=0`$) could provide a realization of a rather asymmetric dimer, the corresponding relative offset being proportional to the (possibly very large) total population, $`w_\mathrm{d}=\frac{1}{2}\rho _\mathrm{t}U_\mathrm{t}`$.

Here we will not give further analytic details about the study of the trimer dynamics (solution of the “asymmetric” fixed point equations, stability analysis), which can be found elsewhere , our present concern being the analogy between the trimer and the dimer. In this respect we remark that the above described exact mapping applies also to the stability character of the SFP only within the “symmetric” phase subspace $`z_1=z_3`$ pertaining to the integrable sub regime. Indeed, when considering the entire phase space, the SFP corresponding to one of the (stable) maxima of the asymmetric dimer may turn out to be unstable. This is clearly shown in the left panel of figure 2, where we plotted the energies the SFPs, as obtained with the straightforward modifications to equation (5) ensuing from the mapping (8): $`\nu _\mathrm{t}(d;\tau _\mathrm{t})=\frac{1}{2}+\frac{3}{4}\nu _\mathrm{d}(d;\tau _\mathrm{d})`$, $`E_\mathrm{t}(d;\tau _\mathrm{t})=\frac{3}{4}E_\mathrm{d}(d;\tau _\mathrm{d})`$, $`\tau _\mathrm{d}=\frac{3}{4\sqrt{2}}\tau _\mathrm{t}`$.

As a further, somewhat unexpected, aspect of the dimer-trimer analogy, the parametric representation (5) allows to describe also the energies of most of the asymmetric fixed points of the trimer, though approximately and for not exceedingly large values of $`\tau _\mathrm{t}`$. An example of this situation is portrayed in the right panel of figure 2, where we plotted the energy of the AFPs (gray markers) along with those of the fixed points of an asymmetric dimer (black lines, circles) with exactly the same parameters as the trimer:

$$U_\mathrm{d}=U_\mathrm{t},w_\mathrm{d}=w_\mathrm{t},T_\mathrm{d}=T_\mathrm{t},\rho _\mathrm{d}=\rho _\mathrm{t}$$
(9)
As we mentioned above, the larger the rescaled hopping amplitude $`\tau _\mathrm{t}`$ the less satisfactory the overlap between the fixed points of dimer and the AFPs of the trimer. However, the comparison among situations relevant to different hopping amplitudes seems to suggest that a more complex parameter mapping than (9) could produce better overlaps. Due to the limited space, here we do not develop this issue, but we rather limit our qualitative discussion to sufficiently small values of the hopping amplitude (up to $`\tau _\mathrm{d}0.5`$) for the simple mapping in (9) to give satisfactory results. In this respect we mention that a closer analysis of the relevant fixed point equations  shows that the dimer analogy is more effective for the AFPs characterized by a very low population of one of the lateral sites. A further feature worth noticing is the correspondence between the dimeric minimum (thick solid line) and the stable saddle points 555Recall that, according to the (linear) stability theory, the maxima and the minima are always stable, whereas the saddle point can be either stable or unstable. In the case of the trimer, the analytically predicted stability of the saddle points has been widely confirmed by numerical simulations . typical of the non-integrable regime of the trimer (gray squares). We conclude this somewhat qualitative analysis of the dimer-trimer correspondence by emphasizing that one of the dimeric maxima (black dashed lines) has no counterpart among the AFPs. As we will briefly discuss in the following, a similar feature appears in the in the purely quantum treatment of the problem, in terms of parity of the energy levels.

### 3.2 Quantum persistence of the dimer-trimer analogy

As we discussed at the beginning of the present section, the Hilbert space for a trimer is so large to require a numerical diagonalization of the (number conserving) Hamiltonian (6) even for relatively small boson populations $`\rho _\mathrm{t}`$. Likewise, the pattern produced by the energy levels of the trimer as the relative energy offset $`w_\mathrm{t}`$ is varied is much more complex than its dimeric counterpart (see e.g. figure 3). Some help in the classification of these energy levels comes from the observation that parity is a good quantum number, which casts some light on the rather complex avoided crossing pattern occurring as $`w_\mathrm{t}`$ varies. The crossing of two levels is either forbidden or allowed depending on whether they have the same or different parity. Indeed levels with different parity belong to independent sectors of the Hilbert space, and there is no hindrance to their crossing.

The discussion presented in section 3.1 provides the very interesting suggestion of interpreting the energy level pattern of the trimer in terms of the spectrum of a dimer with the same population and suitably chosen parameters. Actually, the dimer-trimer analogy persists also at the quantum level, as it is clearly shown in figure 3, where we plotted the 36 levels of a trimer containing $`\rho _\mathrm{t}=7`$ bosons along with the spectrum of two different dimers with the same population. The dimeric spectra appearing in the left and right panels of the figure were obtained from the diagonalization of Hamiltonian (1) with the parameter choices in (8) and (9), respectively. Notice that in both cases each dimeric level follows quite closely a level of the trimer. In this respect we emphasize that the recognition of dimeric structures in the spectrum of (6) is quite difficult even after the introduction of an alternative canonical scheme incorporating the parity symmetry of the trimer , and involves non-trivial manipulations hardly conceivable except as an afterthought following the MF analysis.

Notice that any definite-parity quantum state yields $`n_1n_3=a_1^+a_3a_3^+a_1=0`$, which is a possible quantum counterpart of the MF integrability condition $`z_1=z_3`$666$`a_\mathrm{e}^+a_\mathrm{e}=0`$, where $`a_\mathrm{e}=(a_1+a_3)/\sqrt{2}`$ makes probably a better candidate. Due to limited space, we refer to  for a discussion about this issue.. This is basically the reason why each dimeric level in the right panel of figure 3 overlaps with an odd-even doublet of the trimer: the latter is the minimal set of eigenstates producing a state lacking a definite parity, which is the quantum counterpart of a MF asymmetric configuration. Recall indeed that the dimer under concern is the quantum counterpart of that qualitatively describing the asymmetric fixed points of the trimer. Conversely, though intrinsically asymmetric, the dimer considered in the left panel of figure 3 pertains to the symmetric configurations of the trimer and hence there is no surprise that in some instances a dimeric level corresponds to a solitary level of the trimer rather than to an even-odd doublet. If proved to occur in general, the fact that such solitary levels are always even, as in figure 3, remains to be understood.

Notice that the (even) absolute maximum of the trimer for $`\nu _\mathrm{t}<0`$ seems to have a dimeric counterpart both in the left and in the right panel of figure 3. Actually, as it becomes evident considering a larger value of $`\tau _\mathrm{t}`$, the overlap in the right panel is accidental, and hence it is not an exception to our discussion above. This means that, as we mentioned in section 3.1, one of the dimeric levels in the right panel (the absolute maximum for $`\nu _\mathrm{t}<0`$) has no correspondent among the levels of the trimer.

Owing to the quantum persistence of the dimer-trimer correspondence, most of the results about the asymmetric trimer presented in section 2 can be extended in a straightforward manner to the case of the trimer.

As we mentioned at the beginning of the present section, figure 3 provides a quite useful key for interpreting the otherwise quite obscure trimeric version of figure 1. Actually, owing to the quantum persistence of the dimer-trimer correspondence, most of the results about the asymmetric trimer presented in section 2 are so straightforwardly extended to the case of the trimer, that there is no need of including such figure in the present work.

## 4 Conclusions and perspectives

We focused on the Bose-Hubbard Hamiltonians describing two small arrays of interacting BECs, the asymmetric dimer and the trimer. In section 2 we related some interesting features of the energy spectrum of the former Hamiltonian, whose asymmetry ensues from an imbalance between the on-site energy terms, to the energies of the fixed point of the corresponding mean-field effective Hamiltonian. These results are inspired by and generalize those concerning the symmetric dimer reported in reference . Section 3 is devoted to a similar investigation on the spectrum of the trimer. In order to find our bearings in the rather complex level pattern, in section 3.1 we exploited an interesting mapping among the mean-field dynamics of the trimer and that of the asymmetric dimer discussed in the previous section. In particular we showed that such mapping, which is exact only for the integrable sub regime of the trimer, applies approximately also to the energies of most of its fixed points. Based on these mean-field results, in section 3.2 we observed that a similar approximate mapping exists also at the quantum level. We emphasize that the (quasi) dimeric character of the spectrum of the trimer is by no means evident from the structure of the relevant Bose-Hubbard Hamiltonian, and it strongly relies on the observations reported in section 3.1. The results gathered in sections 3.1 and 3.2, along with the underlying systematic study of references  and , allow to extend straightforwardly the observations reported in section 2 to the case of the trimer.

The results presented in this work form the basis for a systematic study of the purely quantum trimer. Among the issues worth developing we mention a possible characterization of the quantum counterpart of chaos. Recall indeed that, owing to the non-integrability of the mean-field dynamics, the trajectories based in the vicinity of unstable fixed points are likely to display strong instabilities. Hence the mapping between fixed points and energy spectrum may provide a guide in recognizing possible quantum manifestations of chaos.

Furthermore, if the system under concern is about to be realized experimentally, as promised by the recent advances in trapping techniques, the systematic study here proposed could provide a benchmark for verifying the effectiveness of the model Hamiltonians used to describe its physics. In this respect we recall that small bosonic lattices have been proposed for the implementation of quantum computation .

## References


text:
# Anomalous density of states in a metallic film in proximity with a superconductor

## Abstract

We investigated the local electronic density of states in superconductor-normal metal (Nb-Au) bilayers using a very low temperature (60 mK) STM. High resolution tunneling spectra measured on the normal metal (Au) surface show a clear proximity effect with an energy gap of reduced amplitude compared to the bulk superconductor (Nb) gap. Within this mini-gap, the density of states does not reach zero and shows clear sub-gap features. We show that the experimental spectra cannot be described with the well-established Usadel equations from the quasi-classical theory.

At the contact with a superconductor (S), the Andreev reflection of the electrons from a Normal metal (N) locally modifies the N metal electronic properties, including the local density of states (LDOS) Revue . The precise LDOS spectra depend on the N-S structure geometry, in particular the N metal length, and on the electron transport regime. In a diffusive N-S junction with a N metal shorter than the phase coherence length, one expects a fully opened mini-gap, which remains smaller than the superconductor’s energy gap $`\mathrm{\Delta }`$ Belzig-PRB . In a larger N metal, the LDOS shows, within a pseudo-gap, a linear evolution with the energy, reaching zero precisely at the Fermi level. In the ballistic regime, the same distinction holds between a chaotic (mini-gap) and an integrable (pseudo-gap) cavity Beenakker-EPL . In general terms, a mini-gap shows up if every electronic state at the Fermi level can couple to the S interface while maintaining quantum phase coherence. The order of magnitude of the energy gap $`E_g`$ is $`\mathrm{}/\tau _{AR}`$, where $`\tau _{AR}`$ is the characteristic diffusion time before an electron experiences an Andreev reflection. Here, we shall consider the case of diffusive N and S metals brought in contact through a highly transparent interface. Then $`\tau _{AR}L_n^2/D_n`$ and the predicted gap $`E_g`$ is about the Thouless energy $`\mathrm{}D_n/L_{n}^{}{}_{}{}^{2}`$, where $`L_n`$ and $`D_n`$ are the length and the diffusion constant of the normal metal.

Recently, the LDOS of lateral S-N (Nb-Au) structures was probed with solid tunnel junctions Gueron-PRL  and (very) low temperature STM Moussy-EPL ; Vinet-PRB . These studies focused on the pseudo-gap regime in long N metals. Also, a mini-gap was observed in a very thin (20 nm) Au layer on top of a Nb dot Moussy-EPL . A good agreement with the quasi-classical theory based on the Usadel equations Belzig-PRB  and with the Bogoliubov-de Gennes equations Halterman-PRB  was obtained. In the mini-gap regime, a NbSe2 crystal covered with a varying thickness of Au was studied by STM at 2.5 K Dynes-PRL  but the temperature did not allow the observation of a fully open mini-gap. Therefore, the mini-gap regime remains to be investigated, in particular its evolution with the normal metal size and the possible presence of states within the gap. A clear distinction between a fully open gap and a pseudo-gap requires a high energy resolution $`k_BT\mathrm{\Delta }`$, which can be achieved only at very low temperatures ($`T1K`$) and with a large gap $`\mathrm{\Delta }`$.

In this paper, we present measurements of the local density of states at the N metal surface of S-N bilayers with a varying N metal length. We used a very low temperature STM as the local tunneling probe, which permitted us to avoid the averaging over the surface inhomogeneities. The achieved high energy resolution enabled us to unveil a new phenomenon, namely a non-zero density of states appearing in the vicinity of the Fermi level together with clear sub-gap structures, and an anomalous mini-gap amplitude dependence with the N metal length. We demonstrate that these new features cannot be explained within the well-established quasi-classical theory.

We fabricated simultaneously a series of Nb-Au bilayers, with a fixed Nb thickness and a varying Au thickness, on a single (6$`\times `$40 mm2) Si substrate. The Nb layer thickness was chosen to be significantly larger than the superconducting coherence length in order to avoid any effects due to its finite value. For varying the Au thickness, another Si wafer was used as a mask and moved in situ above the substrate. After depositing the 120 nm Nb film, the Au film (from 10 to 260 nm) was deposited within 15 minutes at a pressure below 10-8 mbar. These conditions minimize the interface contamination and should preserve the best Nb-Au interface transparency. The full Si wafer with the bilayer films thus obtained (see Fig. 1a) was cleaved in air to separate the different samples. Individual Au (260 nm) and Nb (120 nm) layers were characterized by transport measurements which are summarized in Table 1.

Each bilayer sample was cooled separately in a STM working at 60 mK in an upside-down dilution fridge Moussy-RSI . A fresh-cut Pt-Ir wire was used as the STM tip. A typical STM image (of the 72 nm sample) acquired at very low temperature (100 mK) is shown in Fig. 1b. We observe a polycrystalline structure with a grain size of 30 - 50 nm. The same grain size was observed in every sample, except in the thinnest one (10 nm). Moreover, this value is consistent with the measured elastic mean free path $`l_{e,n}`$ in the 260 nm sample. Our Au films are thus in the diffusive regime, except for the smallest thicknesses where ballistic effects may occur. On every sample, we acquired series of $`I(V)`$ tunnel characteristics at a tunnel resistance of 5 to 12 $`M\mathrm{\Omega }`$ and at numerous places in order to check reproducibility. The spectra did not depend on the tunnel resistance value and were very reproducible over a sample surface with minor details (in particular small asymmetries at zero bias) slightly evolving, and also punctual defects on a few locations. The $`I(V)`$ data were numerically differentiated to obtain the differential conductance $`dI/dV(V)`$, which gives the LDOS at the energy eV with an accuracy limited by the thermal smearing. A selection of the tunneling spectra is shown in Fig. 2. The spectra are flat at large bias voltages ($`\left|V\right|`$ 2.5 meV) and have been normalized as 1 there.

For the bilayer with the smallest Au thickness (10 nm), the spectrum qualitatively resembles a BCS spectrum with a gap amplitude very close to the expected bulk Nb gap : $`\mathrm{\Delta }1.5`$ meV. For larger thicknesses, the spectra show a mini-gap that reduces in width with increasing Au thickness. There seems to be a crossover from one type of spectra to another at a Au thickness between 72 and 100 nm. For small thicknesses, there is a sharp rise to a peak in the LDOS at $`\mathrm{\Delta }`$, and a relatively slow decrease as we go towards zero bias. For large thicknesses, there is a small dip at $`\mathrm{\Delta }`$, the LDOS then rises slowly to a peak and falls off rapidly inside the mini-gap. This dip at $`\mathrm{\Delta }`$ actually implies that the superconducting gap in the bulk of Nb is not significantly influenced by a large thickness of Au.

Inside the mini-gap, the density of states is strikingly different from zero and features a few sub-gap structures. The observed sub-gap features do not scale with the Nb energy gap $`\mathrm{\Delta }`$, as their energy position evolves with the normal metal length $`L_n`$. At the Fermi level (V = 0), the mean LDOS monotonously increases with the Au thickness. This effect is therefore not a ballistic one which would be restricted to the smaller thicknesses.

We checked the energy resolution of our STM by measuring the LDOS at the surface of plain Nb and Al films and fitting the data with a BCS density of states without any depairing parameter but with a thermal smearing. The fit was very good and gave an effective temperature of 210 mK which corresponds to an energy resolution of 36 $`\mu `$eV Moussy-RSI . This proven high resolution rules out the thermal smearing as a cause of the large (up to 50%) zero bias conductance observed here.

We compared our data to the predictions of the quasi-classical theory in the diffusive regime Golubov-JLTP ; Belzig-PRB , by solving the Usadel equations with the help of a numerical code from W. Belzig et al. Belzig-PRB . The input parameters for this calculation are the lengths $`L_{n,s}`$ of N and S metals in units of the related coherence lengths $`\xi _{n,s}=\sqrt{\mathrm{}D_{n,s}/2\mathrm{\Delta }}`$, the specific interface resistance $`r_B`$, the inelastic $`\gamma _{in}`$ and spin-flip $`\gamma _{sf}`$ scattering rates in the N and S metals, the electronic properties mismatch parameter $`\gamma =\rho _s\xi _s/\rho _n\xi _n`$ $`(\sqrt{l_{e,n}/l_{e,s}})`$. Here $`\rho _{n,s}`$ is the normal-state resistivity in N or S. We considered the interface transparency as perfect and we took into account a thermal smearing with an effective temperature of 210 mK.

We adjusted the calculation parameters to fit the tunneling spectra. Fig. 3 displays a comparison of experimental spectra with two sets of calculated curves. We first assumed that the inelastic and spin-flip scattering lengths are much larger than the Au layer thickness, so that we can neglect these processes (Fig. 3 dotted lines). In order to recover the observed LDOS peak amplitude and position, we had to assume a mismatch parameter $`\gamma `$ value of 0.6 instead of 2.6 as estimated from the transport measurements, keeping the other parameters matching precisely the measured values. This discrepancy may be related to the Nb-Au interface roughness which is expected to affect the Andreev reflection rate, and thus the induced superconductivity in the Au layer Pilgram-PRB . The overall spectral shapes are qualitatively reproduced, except for the peak of the 50 nm spectra which is noticeably different from the experimental data. At large $`L_n`$, the predicted mini-gap is clearly smaller than the one observed experimentally. Moreover, we always obtain a fully opened gap with a zero LDOS at the Fermi level. This is in clear disagreement with our experimental results, where we always get a non-zero conductance at zero bias. The dependence of the mean free path with the Au layer thickness $`L_n`$ Fuchs-PCPS  was taken into account through the related variation of the characteristic length $`\xi _n`$ and of the mismatch parameter $`\gamma `$. This modified only slightly the calculated spectra in the regime $`L_n`$$`<`$ 100 nm and did not improve the fit.

With the hope to obtain a better fit, we included a non-zero inelastic scattering rate in the calculation parameters (Fig. 3 dashed lines). We achieved a good agreement of the calculated LDOS with the measured one, but only at the precise Fermi level, by choosing $`\gamma _{in}=0.05\mathrm{\Delta }/\mathrm{}`$. This value corresponds to an inelastic mean free path of about 370 nm, which is significantly smaller than the expected value of about 2 $`\mu m`$. At large thicknesses, the calculated mini-gap remains smaller than in the experiment. Fitting the mini-gap width requires a significant modification of the calculation parameters such that the calculated LDOS close to the Fermi level also diminishes. This makes the disagreement between the experimental data and calculated LDOS even stronger for a fixed value of $`\gamma _{in}`$. We also tried to include a spin-flip scattering Gueron-PRL ; Vinet-PRB  or a small interface resistance in the calculation, but it did not improve the fit.

From the experimental tunneling spectra, we extracted the mini-gap $`E_g`$ as the energy at the inflection point, i.e., with a maximum slope $`d^2I/d^2V`$. The quantity $`(E_g/\mathrm{\Delta })^{1/2}`$ is plotted in Fig. 4 as a function of the N metal length $`L_n`$, so that a Thouless energy scaling should appear as a linear behavior. The calculated mini-gap is close to this scaling : $`E_g=1.3\times 10^4/(76+L_n)^2`$ meV, $`L_n`$ being in nm units. As naively expected, the numerator value is of the order of $`\mathrm{}D_n=1.09\times 10^4`$meV.nm2, based on Table 1 data. The extra length appearing in the numerator expresses the fact that the Andreev reflection happens in the superconductor over this length scale. It also makes $`E_g`$ extrapolate to $`\mathrm{\Delta }`$ as $`L_n`$ goes to zero. The experimental data follows this behavior in the small thickness regime $`L_n`$$`<`$ 100 nm. At large $`L_n`$, the discrepancy between theory and experiment is striking.

Up to now, we described samples with an exceptional interface transparency. We have also fabricated similar Nb-Au samples but with a degraded interface, thanks to an Ar etch of the Nb surface performed just before a 20 nm Au deposition. Different Ar etch times were used, keeping Ar pressure (3$`.10^4`$mbar) and beam voltage (100 eV) fixed. Compared to Fig. 2 data at similar thickness, Fig. 5 tunneling spectra show a drastic change in the LDOS. A mini-gap is observed, with a decreasing amplitude as the etching time is increased. We interpret this behavior as due to a decrease of the interface transparency induced by the surface milling. This enhances the contribution of specular reflection of electrons as compared to the Andreev reflection, which makes the time $`\tau _{AR}`$ increase and the mini-gap $`E_g`$ decrease. Interestingly, these spectra also show a non-zero density of states at the Fermi level, which seems insensitive to the variation of the interface transparency. Although we were able to qualitatively fit the mini-gap width and the peak profile with the Usadel solutions by taking into account a non-zero interface resistance (see Fig. 5 inset), it was again impossible to fit the whole spectra, especially the non-zero LDOS near the Fermi level. A moderate interface transparency thus causes an anomalous LDOS, in a similar way as a large N metal length.

The observed LDOS, especially its amplitude at the Fermi level and the mini-gap width evolution, cannot be described by the quasi-classical theory, even after assuming an extremely strong phase-breaking scattering. In the diffusive regime of relevance here, neither the angular dependence of the Andreev reflection at an imperfectly transparent N-S interface Mortensen-PRB  nor the presence of surface electron states Everson-JVST  should play a role. In terms of electron trajectories, a finite LDOS at the Fermi level means that a significant fraction of the electron states present at the Au surface do not Andreev-reflect at the Nb-Au interface on the time scale of the phase coherence time.

In our opinion, this new effect was not observed before because previous studies either did not feature the relevant thickness range and a confined geometry where the mini-gap enables a clear observation of additional low-energy single particle states Moussy-EPL ; Vinet-PRB , or were carried out at a higher temperature with insufficient energy resolution in the measured LDOS Dynes-PRL . Still, our observations are somewhat reminiscent of some unexplained features reported on the LDOS of presumably quite disordered 20 nm N metal ridges lying over a Nb surface Vinet-PRB . There also, a clear filling of the LDOS at the Fermi level was observed in the absence of any mini-gap width evolution. Let us point out that both this and our systems are far from the strongly disordered regime where a non-zero LDOS for all energies is predicted Feigelman-PRL . Rather, our films have the typical polycrystalline structure of a weakly disordered thin film. The elastic diffusion in the N metal is not homogenous as in an amorphous metal but controlled by the interfaces between ballistic grains. The possible confinement of electron states within a grain can have serious effects on the LDOS in the proximity superconductivity regime of interest here. A significant fraction of electron states may be sufficiently decoupled from the N-S interface to avoid any Andreev reflection. These states would thus remain insensitive to the superconductivity in the S layer and contribute fully to the LDOS.

In conclusion, we investigated with a high energy resolution the LDOS at the N metal surface of diffusive S-N (Nb-Au) bilayers and compared to the theoretical predictions. We uncovered an anomalous mini-gap width evolution and a non-zero LDOS near the Fermi level that increases as the N metal thickness is increased. The latter behavior is also observed at small thickness when the interface transparency is reduced on purpose. This new effect cannot be described within the quasi-classical theory. We suggest that the granular structure of the Au layer enables the decoupling of some electron states from the Andreev reflection at the interface, which results in a non-zero LDOS inside the mini-gap.

We thank W. Belzig for sharing his numerical code, M. Feigelman for discussions and the ACI ”Nanostructures” for financial support.


text:
# A comparison between several correlated stochastic volatility models

, ,

## Abstract

We compare the most common SV models such as the Ornstein-Uhlenbeck (OU), the Heston and the exponential OU (expOU) models. We try to decide which is the most appropriate one by studying their volatility autocorrelation and leverage effect, and thus outline the limitations of each model. We add empirical research on market indices confirming the universality of the leverage and volatility correlations.

The multiplicative diffusion process is considered to be the most popular model in finance. We thus have that the log-price change without average describes the diffusive process $`\dot{X}(t)=\sigma \xi _1(t)`$ (in the Itô sense), where $`\sigma `$ is the volatility (assumed to be constant) and $`\xi _1`$ is Brownian noise with $`\xi _1(t)\xi _1(t^{})=\delta (tt^{})`$. However, this model is unable to capture most of the statistical properties of real markets, those called stylized facts.

The stochastic volatility (SV) models are a natural way out to avoid the inconsistencies of the log-Brownian model by still assuming that $`X`$ follows a diffusion process but now with a random $`\sigma `$. However, the dynamics of $`\sigma `$ has not been definitively attached to any specific SV model. In this sense, we want to discern between the most common SV models. Among others we have:

(a) The Ornstein-Uhlenbeck (OU) SV model

$$\dot{\sigma }(t)=\alpha (\sigma m)+k\xi _2(t),\sigma (t)=m+k_{\mathrm{}}^te^{\alpha (tt^{})}\xi _2(t^{})𝑑t^{}.$$
(1)
(b) The Heston SV model  (with $`V\sigma ^2`$)

$`\dot{V}(t)=\alpha (Vm^2)+k\sqrt{V}\xi _2(t),`$ (2)
$`V(t)=m^2+k{\displaystyle _{\mathrm{}}^t}e^{\alpha (tt^{})}\sqrt{V(t^{})}\xi _2(t^{})𝑑t^{}.`$ (3)
(c) And, finally, the exponential Ornstein-Uhlenbeck (expOU) SV model  (with $`Y\mathrm{ln}(\sigma /m)`$)

$$\dot{Y}(t)=\alpha Y+\frac{k}{m}\xi _2(t),Y(t)=m+\frac{k}{m}_{\mathrm{}}^te^{\alpha (tt^{})}\xi _2(t^{})𝑑t^{}.$$
(4)
Note that the Heston model does not provide a closed equation for the volatility. This can be a serious drawback in case we want to deal with analytical expressions. In addition, we suppose that noises are correlated (i.e., $`\xi _1(t)\xi _2(t^{})=\rho \delta (tt^{})`$ with $`1\rho 1`$). This allows us to explain some stylized facts.

The leverage effect  is measured by: $`L(\tau )dX(t+\tau )^2dX(t)/dX(t)^2`$. For the models above we have

$`L_{\mathrm{OU}}(\tau )=2\rho \left[{\displaystyle \frac{\nu \sqrt{2\alpha }\left(1+\nu ^2e^{\alpha \tau }\right)}{(1+\nu ^2)^2m}}\right]e^{\alpha \tau },L_\mathrm{H}(\tau )4\rho {\displaystyle \frac{\alpha }{k}}e^{\alpha \tau },`$ (5)
$`L_{\mathrm{expOU}}(\tau )={\displaystyle \frac{2\rho k}{m^2}}e^{\nu ^2/2}\mathrm{exp}\left(2\nu ^2e^{\alpha \tau }\alpha \tau \right),`$ (6)
for $`\tau >0`$ and $`L(\tau )=0`$ for negative $`\tau `$. It seems impossible to derive an analytical expression for the Heston leverage because $`V(t)`$ does not have a closed expression (cf. Eq. (3)). The leverage given in Eq. (5) is exact only when $`m^2=k^2/2\alpha `$ but Ref.  shows that this is not true for the Dow Jones. We will have thus to simulate the process to compute the leverage effect in each specific Heston case. Eqs. (5)–(6) are plotted and compared with data in Fig. 1. The data is too noisy to assert which is the most appropriate model. In any case, the three models have the desired exponential decay with a characteristic time of the order of few days. We also note that, without the correlation coefficient $`\rho `$, there is no leverage and data confirms that $`\rho <0`$ .

The volatility autocorrelation is $`C(\tau )[dX(t)^2dX(t+\tau )^2dX(t)^2dX(t+\tau )^2]/\text{Var}[dX(t)^2]`$. In this case, we can exactly derive this correlation for each model. We have

$`C_{\mathrm{OU}}(\tau )={\displaystyle \frac{\nu ^2e^{\alpha \tau }(\nu ^2e^{\alpha \tau }+1)}{4\nu ^2(2+\nu ^2)+1}},C_\mathrm{H}(\tau )=\nu ^2e^{\alpha \tau },`$ (7)
$`C_{\mathrm{expOU}}(\tau )={\displaystyle \frac{\mathrm{exp}\left(4\nu ^2e^{\alpha \tau }\right)1}{3\mathrm{exp}(4\nu ^2)1}}.`$ (8)
The volatility autocorrelation observed in empirical data shows at least two different time scales (see Fig. 2 and Ref. ). The shortest one coincides with the leverage correlation time while the second scale is of the order of years, around 10 times bigger than the leverage time. The OU and Heston models are not capable of reproducing the slowest time decay but the expOU model has a nontrivial decay with the two time scales:

$$C_{\mathrm{expOU}}(\tau )e^{\alpha \tau }(\alpha \tau >1),C_{\mathrm{expOU}}(\tau )e^{(k/m)^2\tau }(\alpha \tau <1),$$
(9)
where the faster decay is given by the diffusive coefficient $`k/m`$ and the slow decay of the volatility is provided by the reversion coefficient $`\alpha `$ (cf. Eq. (4)).

Therefore, we conclude that the expOU model appears to be as more realistic. Further investigations on this model are required to confirm this assertion. However, an alternative choice is to sophisticate the OU and Heston models. A second time scale can be generated by including a third SDE for $`m`$. In a previous work , we have studied this possibility for the OU model with some success by taking: $`\dot{m}(t)=\alpha _m(m\theta )+k_m\xi _3(t)`$. A similar analysis can be performed for the Heston model.

This work has been supported in part by Dirección General de Proyectos de Investigación under contract No. BFM2003-04574, and by Generalitat de Catalunya under contract No. 2001SGR-00061.