format: text

text:
# Introduction

## Introduction

The integrable structure of the Calogero-type systems is of special interest concerning such diverse problems as hierarchies of partial differential equations , Seiberg-Witten solutions of N=2 SUSY , quasi-Hopf deformations of quantum groups and many others. There are many different approaches to describe the Hamiltonian structure of such systems. One of the aims of this paper is to combine the Hitchin description and the language of the universal algebraic-geometric symplectic form. The method of the special inverse spectral problem was successfully applied in the matrix KP case in . This method deals with the auxiliary linear problem coming from the Lax representation of the equation. The Baker-Akhiezer function (FBA) turns to be a function on the spectral algebraic curve, and its specific analytic properties are very restrictive and provide the reconstruction of the potential and the poles of FBA in terms of theta-functions on the Jacobian of the spectral curve. The same procedure was applied to the spin Ruijsenaars-Schneider system in . But for a wide class of integrable systems of the Hitchin type the explicit solutions remain still unknown. In the first section we revisit the construction of the Hitchin system on the degenerate curve and especially the Gaudin elliptic two-point N-particle system with spin. In the second section we investigate the finite gap potentials of the matrix Davey-Stewartson equation. The system of non-linear equations of its poles dynamics is deduced and the algebraic-geometrical correspondence between the phase space of this system and the space of spectral data is constructed. In the third part the FBA is reconstructed in terms of $`\theta `$-functions, the explicit expressions for the potential and for the solution of the poles dynamics are given. The universal symplectic structure is constructed in the fourth section. We also demonstrate that the system of poles dynamics is hamiltonian and identify it with the Gaudin elliptic system with spin. In the last section we discuss the rational and special modular degenerations of the system.

## 1 Gaudin system

The Hitchin system was introduced in  as an integrable system on the cotangent bundle of the moduli space of stable holomorphic bundles on the algebraic curve $`\mathrm{\Sigma }`$. The phase space was obtained by the Hamiltonian reduction from the space of pairs $`d_A^{^{\prime \prime }},\mathrm{\Phi }`$, where $`d_A^{^{\prime \prime }}`$ is the operator defining the holomorphic structure on the bundle $`V`$ and $`\mathrm{\Phi }`$ is the endomorphism of this bundle, more precisely $`\mathrm{\Phi }\mathrm{\Omega }^{0,1}(\mathrm{\Sigma },End(V))`$. The reduction by the action of the loop group $`GL_N(z)`$ was done, i.e. the group of $`GL_N`$-valued functions on $`\mathrm{\Sigma }`$. The invariant symplectic structure on the base space writes

$$\omega =_\mathrm{\Sigma }Tr\delta \mathrm{\Phi }\delta d_A^{^{\prime \prime }}.$$
(1.1)
The zero level of the moment map is described by the condition $`d_A^{^{\prime \prime }}\mathrm{\Phi }=0`$ which means that $`\mathrm{\Phi }`$ is holomorphic with respect to the induced holomorphic structure on the bundle $`End(V)`$. It turns out that the system of quantities $`Tr\mathrm{\Phi }^k`$, treated as vector functions on the phase space, Poisson-commute and their number is exactly half the dimension of the phase space. This system was stated in a formal way and its explicit formulation was obtained for special cases in .

The generalized construction of the Hitchin systems on degenerate curves  provides a more explicit description. One considers the moduli space of semi-stable bundles $`V`$ on an algebraic curve with singularities. The Higgs field in this case is a meromorphic section of the bundle $`End(V)`$. The holomorphic bundle on the singular curve is the holomorphic bundle on its normalization and the set of additional data characterizing the bundle over singularities.

Let us consider the elliptic curve with two simple fixed points $`z_1,z_2`$. The Higgs field is a meromorphic section with two poles on the curve. The moduli space $``$ of semi-stable bundles on the elliptic curve with such type of singularity decomposes into the moduli space of semi-stable bundle $`V`$ on the normalized curve and the finite dimensional space of additional data at the fixed points: $`g_iEnd(V|_{z_i})`$. The cotangent bundle $`𝒯^{}`$ has additional coordinates $`p_i𝒯^{}End(V|_{z_i})`$ where $`i=1,2`$. The modified symplectic structure is

$$\omega =\delta \mathrm{\Phi }\delta d_A^{^{\prime \prime }}+Tr(\delta (g_1^1p_1)\delta g_1+\delta (g_2^1p_2)\delta g_2)$$
(1.2)
and the zero level of the modified moment map is

$$\overline{}\mathrm{\Phi }+[A,\mathrm{\Phi }]+p_1\delta ^2(z_1)+p_2\delta ^2(z_2)=0.$$
(1.3)
We restrict ourselves to the subspace of semi-stable holomorphic bundles of rank $`N`$ which decompose into the sum of line bundles. In the basis associated to the decomposition we can make the matrix $`A`$ diagonal with elements $`x_1,\mathrm{},x_n`$. The equation (1.3) becomes

$$\overline{}\mathrm{\Phi }_{ij}+(x_ix_j)\mathrm{\Phi }_{ij}+(p_1)_{ij}\delta ^2(z_1)+(p_2)_{ij}\delta ^2(z_2)=0.$$
(1.4)
It could be solved as in  in the following way

$$\mathrm{\Phi }_{ij}=\frac{\mathrm{exp}((x_ix_j)\frac{z\overline{z}}{\tau \overline{\tau }})}{2\pi i}((p_1)_{ij}\frac{\sigma (x_ix_j+zz_1)}{\sigma (x_ix_j)\sigma (zz_1)}+(p_2)_{ij}\frac{\sigma (x_ix_j+zz_2)}{\sigma (x_ix_j)\sigma (zz_2)});$$
(1.5)
$$\mathrm{\Phi }_{ii}=\omega _i+(p_1)_{ii}\zeta (zz_1)+(p_2)_{ii}\zeta (zz_2);$$
(1.6)
where $`\tau `$ is the modular parameter of the elliptic curve, $`\sigma `$ is the Weierstraß $`\sigma `$-function and $`\omega _i`$ are additional free parameters. Now $`res_{z=z_1}4\pi ^2Tr\mathrm{\Phi }^2`$ is equal to

$$H_1=\underset{i}{}\omega _i(p_1)_{ii}+\underset{i}{}(p_1)_{ii}(p_2)_{ii}\zeta (z_1z_2)+\underset{ij}{}(p_1)_{ij}(p_2)_{ji}\frac{\sigma (x_ix_j+z_1z_2)}{\sigma (x_ix_j)\sigma (z_1z_2)},$$
(1.7)
and the residue at the point $`z=z_2`$ gives a similar expression for the Hamiltonian $`H_2`$ which differs from $`H_1`$ by the interchanging of $`z_1`$ and $`z_2`$. These functions are called the Hamiltonians of the Gaudin elliptic system with spin. We are interested in the complete integrability of this system. On the moduli space of holomorphic bundles on the regular curve in  the number of coordinates on the set of quantities $`Tr\mathrm{\Phi }^k`$, i. e. the number of Beltrami differentials, was proved to be exactly half the dimension of $`𝒯`$. Our case differs from the classical one and we calculate the dimension directly. The coefficients of the principal part of $`Tr\mathrm{\Phi }^k`$ at the vicinity of the fixed points are the integrals. The mentioned quantities are elliptic functions, so the sum of their residues is equal to zero. The number of integrals is equal to $`N(N+1)`$. The number of coordinates on the phase space $`x_i,\omega _i(p_1)_{ij},(p_2)_{ij},(g_1)_{ij},(g_2)_{ij}`$ is equal to $`2N+4N^2`$. The symplectic structure (1.2) is invariant under the action of the group $`GL_N\times GL_N`$ in the following way: the element $`(g_1^0,g_2^0)`$ acts as

$$g_1g_1g_1^0,g_2g_2g_2^0.$$

The moment map of this action is $`\mu =((g_1)^1p_1g_1,(g_2)^1p_2g_2)`$ and we can fix it to have a diagonal form and this kills $`2N(N1)`$ degrees of freedom. The stabilizer of this choice of moment map consists of the diagonal subgroup $`Diag\times DiagG\times G.`$ It means that the choice of the representative of the stabilizer’s orbit in the moment map subspace kills $`2N`$ additional degrees of freedom. The dimension of the reduced phase space is equal to $`2N+4N^22N(N1)2N=2N(N+1)`$ so that the number of integrals is exactly half the number of variables.

## 2 Davey-Stewartson equation

As was shown in , the Davey-Stewartson equation

$$iu_t+u_{xx}u_{yy}+u\varphi =0\varphi _{xx}+\varphi _{yy}=2(|u|_{xx}^2|u|_{yy}^2)$$
(2.1)
has a Lax representation $`\dot{L}=[L,M]`$ with

$$L=\left(\begin{array}{cc}_z& 0\\ 0& _{\overline{z}}\end{array}\right)\frac{1}{2}\left(\begin{array}{cc}0& u\\ v& 0\end{array}\right)z=x+iy.$$

Let us consider the $`l`$-dimensional generalization of the auxiliary linear problem related to this operator

$$\left(_t+\left(\begin{array}{cc}d_1& 0\\ 0& d_2\end{array}\right)_x+\left(\begin{array}{cc}0& A_1\\ A_2& 0\end{array}\right)\right)\mathrm{\Psi }=D_A\mathrm{\Psi }=0,$$
(2.2)
where $`A_1,A_2`$ are $`l\times l`$-matrices, $`d_1=E,d_2=E`$, and $`E`$ is the identity $`l\times l`$-matrix, and $`\mathrm{\Psi }`$ is a $`2l`$-dimensional vector function. Looking for the quasi-periodic solutions of the matrix Davey-Stewartson equation we put the potential in the special form:

$$A_1=\underset{j=1}{\overset{n}{}}a_1^jb_{1}^{j}{}_{}{}^{+}\frac{\sigma (xx_j+\eta )}{\sigma (\eta )\sigma (xx_j)},A_2=\underset{j=1}{\overset{n}{}}a_2^jb_{2}^{j}{}_{}{}^{+}\frac{\sigma (xx_j\eta )}{\sigma (\eta )\sigma (xx_j)},$$
(2.3)
where $`a_k^j`$ and $`b_k^j`$ are $`l`$-vectors, and $`\sigma `$ is a Weierstraß function. As usual we look for the quasi-periodic solutions of the linear problem. The convenient basis of such functions with simple poles on $`x`$ is provided by

$$\mathrm{\Phi }_i(z,x)=\frac{\sigma (zz_i+x)}{\sigma (x)\sigma (zz_i)}\mathrm{exp}(\zeta (z)x),\text{where}z_2z_1=\eta .$$

We call $`\mathrm{\Psi }_1`$ the first $`l`$ components of the function $`\mathrm{\Psi }`$, and $`\mathrm{\Psi }_2`$ the last $`l`$ components. The blocks of the potential have a monodromy such that the ratio of the monodromies of $`\mathrm{\Psi }_1`$ and $`\mathrm{\Psi }_2`$ must be equal to the monodromy of $`A_1`$. We look for the eigenfunction in the form:

$$\mathrm{\Psi }_i=\underset{j=1}{\overset{n}{}}S_i^j\mathrm{\Phi }_i(z,xx_j)\mathrm{exp}(ky_++\frac{1}{k}y_{}),$$
(2.4)
where $`y_\pm =\frac{1}{2}(t\pm x)`$. We also look for the solution of the dual equation

$$\mathrm{\Psi }^+D_A=0$$
(2.5)
in the form

$$\mathrm{\Psi }_i^+=\underset{j=1}{\overset{n}{}}S_{i}^{j}{}_{}{}^{+}\mathrm{\Phi }_i(z,x_jx)\mathrm{exp}(ky_+\frac{1}{k}y_{}).$$
(2.6)
###### Theorem 1

The equation (2.2) has a solution of the form (2.4) and the equation (2.5) has a solution of the form (2.6) if and only if the poles dynamics is described by the system:

$$\dot{a}_1^j+\lambda _1^ja_1^j+\underset{lj}{}a_1^l(b_{1}^{l}{}_{}{}^{+}a_2^j)\frac{\sigma (x_jx_l+\eta )}{\sigma (\eta )\sigma (x_jx_l)}=0;$$
(2.7)
$$\dot{a}_2^j+\lambda _2^ja_2^j\underset{lj}{}a_2^l(b_{2}^{l}{}_{}{}^{+}a_1^j)\frac{\sigma (x_jx_l\eta )}{\sigma (\eta )\sigma (x_jx_l)}=0;$$
(2.8)
$$\dot{b}_1^j\lambda _2^jb_1^j+\underset{lj}{}b_1^l(b_{2}^{j}{}_{}{}^{+}a_1^l)\frac{\sigma (x_jx_l+\eta )}{\sigma (\eta )\sigma (x_jx_l)}=0;$$
(2.9)
$$\dot{b}_2^j\lambda _1^jb_2^j\underset{lj}{}b_2^l(b_{1}^{j}{}_{}{}^{+}a_2^l)\frac{\sigma (x_jx_l\eta )}{\sigma (\eta )\sigma (x_jx_l)}=0;$$
(2.10)
$$\ddot{x}_j=\underset{lj}{}\left((b_{1}^{l}{}_{}{}^{+}a_2^j)(b_{2}^{j}{}_{}{}^{+}a_1^l)\frac{\sigma (x_jx_l+\eta )}{\sigma (\eta )\sigma (x_jx_l)}(b_{1}^{j}{}_{}{}^{+}a_2^l)(b_{2}^{l}{}_{}{}^{+}a_1^j)\frac{\sigma (x_jx_l\eta )}{\sigma (\eta )\sigma (x_jx_l)}\right);$$
(2.11)
$$\dot{\lambda }_1^j\dot{\lambda }_2^j=\underset{lj}{}(b_{2}^{j}{}_{}{}^{+}a_1^l)(b_{1}^{l}{}_{}{}^{+}a_2^j)\frac{\sigma (x_jx_l+\eta )}{\sigma (\eta )\sigma (x_jx_l)}(\zeta (x_jx_l)\zeta (x_jx_l+\eta )+2\zeta (\eta ))+$$

$$\underset{lj}{}(b_{1}^{j}{}_{}{}^{+}a_2^l)(b_{2}^{l}{}_{}{}^{+}a_1^j)\frac{\sigma (x_jx_l\eta )}{\sigma (\eta )\sigma (x_jx_l)}(\zeta (x_jx_l)\zeta (x_jx_l\eta )2\zeta (\eta )).$$
(2.12)
Proof As was done in , introducing the expression (2.4) into equation (2.2) and eliminating the essential singularity, we calculate the coefficients of the expansion in $`x`$ at the vicinity of the poles $`x=x_j`$. The lowest term is of order $`2`$

$$S_1^j(\dot{x_j}1)+a_1^jb_{1}^{j}{}_{}{}^{+}S_2^j=0;$$

$$S_2^j(\dot{x_j}+1)+a_2^jb_{2}^{j}{}_{}{}^{+}S_1^j=0.$$

The similar condition from the dual equation gives

$$S_{1}^{j}{}_{}{}^{+}(\dot{x_j}1)+S_{2}^{j}{}_{}{}^{+}a_2^jb_{2}^{j}{}_{}{}^{+}=0;$$

$$S_{2}^{j}{}_{}{}^{+}(\dot{x_j}+1)+S_{1}^{j}{}_{}{}^{+}a_1^jb_{1}^{j}{}_{}{}^{+}=0.$$

We can fix the gauge of variables $`a_k^j,b_k^j`$ such that

$$1\dot{x_j}=(b_{1}^{j}{}_{}{}^{+}a_2^j),1\dot{x_j}=(b_{2}^{j}{}_{}{}^{+}a_1^j),$$

and so

$$S_k^j=c^ja_k^j,S_{k}^{j}{}_{}{}^{+}=c_{}^{j}{}_{}{}^{+}b_{\overline{k}}^j\overline{k}=3k,$$
(2.13)
where $`c^j`$ are scalar functions. The terms of order $`1`$ at $`x=x_j`$ implying (2.13) can be represented as (2.7),(2.8) where $`\lambda _k^j`$ is defined from the following equations on $`c^i`$

$$(_t+kI+L_1)C=0;$$
(2.14)
$$(_t+\frac{1}{k}I+L_2)C=0,$$
(2.15)
where

$$(L_1)_{ij}=\delta _{ij}((b_{1}^{i}{}_{}{}^{+}a_2^i)(\zeta (zz_2)\zeta (z)\zeta (\eta ))\lambda _1^i)+(1\delta _{ij})(b_{1}^{i}{}_{}{}^{+}a_2^j)\mathrm{\Phi }_2(x_ix_j,z);$$

$$(L_2)_{ij}=\delta _{ij}((b_{2}^{i}{}_{}{}^{+}a_1^i)(\zeta (zz_1)\zeta (z)+\zeta (\eta ))\lambda _2^i)+(1\delta _{ij})(b_{2}^{i}{}_{}{}^{+}a_1^j)\mathrm{\Phi }_1(x_ix_j,z).$$

Concerning the dual equation, we obtain similar equations for the variables $`c_{}^{i}{}_{}{}^{+}`$ and the equations (2.9), (2.10). Thus for the Lax equation we have

$$_t(L_1L_2)=[L_1,L_2],$$

which is equivalent to the equations (2.11),(2.12). $`\mathrm{}`$

The standard trick in the method of algebraic-geometrical direct and inverse problem is that the spectral curve $`\mathrm{\Gamma }`$, defined by the characteristic equation $`det(k\frac{1}{k}+L_1L_2)=0`$, has specific properties and that $`\mathrm{\Psi }`$-function defines a holomorphic line bundle on the curve, and its essential singularity depending on time is the deformation of the fixed line bundle. The analysis is rather traditional and we venture to outline main propositions on the analytic properties of eigenfunction and properties of the spectral curve. We use the notation $`L(z)=L_1(z)L_2(z)`$ and $`\widehat{k}=k\frac{1}{k}.`$ First, noting that the matrix $`res_{z=z_i}L(z)`$ differs from the unity matrix by a matrix of rank $`l`$ we find the decomposition at the points $`z_1,z_2`$

$$R(\widehat{k},z)=\underset{i=1}{\overset{N}{}}(\widehat{k}+\nu _i^{1,2}z^1+h_i^{1,2}(z))$$
(2.16)
where $`h_i^{1,2}(z)`$ are regular and $`\nu _i^{1,2}=1`$ if $`i>l`$. Using that $`\mathrm{\Gamma }`$ is an $`N`$-sheeted covering of the base curve $`z`$ and the definition of the branch points of the curve, we count the number of branch points; it is equal to $`4Nl2l(l+1)`$, and by Riemann-Hurwitz theorem the genus of $`\mathrm{\Gamma }`$ is $`g=2Nll(l+1)+1.`$ In virtue of the Lax representation the parameters of the spectral curve are integrals of motion. If, in addition, we fix $`c(t)|_{t=0}=const`$ it can be deduced that the poles of vector $`c`$ do not depend on time. In general position the number of such poles is half the number of branch points, and it is equal to $`2Nll(l+1)=g1.`$ The $`\mathrm{\Psi }`$ function has essential singularities at the points over $`z=z_1,z=z_2`$ and additional $`2l`$ poles. Summarizing the analytic properties of the eigenfunction we obtain

###### Theorem 2

The function $`\mathrm{\Psi }`$ is meromorphic on the curve $`\mathrm{\Gamma }`$ of genus $`g`$ outside $`2l`$ points $`P_i`$, and has $`g+2l1`$ poles $`\gamma _1,\mathrm{},\gamma _{g+2l1}`$ which do not depend on $`t,x`$. At the points $`P_i,i=1,\mathrm{},l`$ over $`z=z_1`$ it has the following decomposition

$$\mathrm{\Psi }_\alpha (x,t,P)=(\chi _0^{\alpha i}+\underset{s=1}{\overset{\mathrm{}}{}}\chi _s^{\alpha i}(x,t)z^s)\mathrm{exp}(\lambda _i(x+t))\mathrm{\Psi }_1(0,0,P)$$

and when $`i=l+1,\mathrm{},2l`$ at the points over $`z=z_2`$ we have

$$\mathrm{\Psi }_\alpha (x,t,P)=(\chi _0^{\alpha i}+\underset{s=1}{\overset{\mathrm{}}{}}\chi _s^{\alpha i}(x,t)z^s)\mathrm{exp}(\lambda _i(xt))\mathrm{\Psi }_1(0,0,P)$$

where $`\lambda _i(z)=(1\nu _i)z^1h_i(z)`$, $`h_i`$ are regular functions in the vicinity of $`z=0`$, $`\nu _i`$ are eigenvalues of the matrix $`L`$ at the points over $`z=z_1,z=z_2`$ respectively and $`\chi _0^{\alpha i}`$ are constants.

The similar proposition on the analytic properties of the solution of dual equation is true

###### Theorem 3

The solution of the dual equation $`\mathrm{\Psi }^+`$ is a meromorphic function on the same curve $`\mathrm{\Gamma }`$ of genus $`g`$ outside the points $`P_i,i=1,\mathrm{},2l`$ It has $`g+2l1`$ poles $`\gamma _1^+,\mathrm{},\gamma _{g+2l1}^+`$ which do not depend on $`t,x`$. At the points $`P_i,i=1,\mathrm{},l`$ over $`z=z_1`$ it has the following decomposition

$$\mathrm{\Psi }_\alpha ^+(x,t,P)=(\chi _0^{+,\alpha i}+\underset{s=1}{\overset{\mathrm{}}{}}\chi _s^{+,\alpha i}(x,t)z^s)\mathrm{exp}(\lambda _i(x+t))\mathrm{\Psi }_1^+(0,0,P),$$

and at the points $`P_i`$ over $`z=z_2`$, $`i=l+1,\mathrm{},2l`$ one has

$$\mathrm{\Psi }_\alpha ^+(x,t,P)=(\chi _0^{+,\alpha i}+\underset{s=1}{\overset{\mathrm{}}{}}\chi _s^{+,\alpha i}(x,t)z^s)\mathrm{exp}(\lambda _i(xt))\mathrm{\Psi }_1^+(0,0,P),$$

where $`\chi _0^{+,\alpha i}`$ are constants.

## 3 Reconstruction formulas

In this section we demonstrate that the analytic properties that we found in the previous section are sufficiently restrictive and allow us to reconstruct the eigenfunction in terms of the characteristics of the curve, such as abelian differentials, $`\theta `$-functions and Abel transform. The standard proposition about the uniqueness of the function with the following properties

1. $`\mathrm{\Psi }_\alpha (x,t,P)`$ has $`g+2l1`$ poles $`\gamma _1,\mathrm{},\gamma _{g+2l1}`$ on the curve $`\mathrm{\Gamma }`$ of genus $`g`$
2. it has essential singularities at the points $`P_j,j=1,..,l`$ of the form

$$\mathrm{\Psi }_\alpha (x,t,P)=e^{w_j^1(x+t)}(\delta _{\alpha j}+\underset{s}{}F_{\alpha j}^s(x,t)w_j^s)$$

and at the points $`P_j,j=l+1,\mathrm{},2l`$ of the form

$$\mathrm{\Psi }_\alpha (x,t,P)=e^{w_j^1(xt)}(\delta _{\alpha j}+\underset{s}{}F_{\alpha j}^s(x,t)w_j^s)$$

where $`w_j(P)`$ is the local parameter at $`P_j`$

arises from the Riemann-Roch theorem. And its existence follows from the construction: let $`d\mathrm{\Omega }^1,d\mathrm{\Omega }^2`$ be unique meromorphic differentials, holomorphic outside the points $`P_j`$ such that:

$$d\mathrm{\Omega }^1=d(w_j^1+O(w_j))\text{ }atthepointsP_j$$

$$d\mathrm{\Omega }^2=\{\begin{array}{cc}d(w_j^1+O(w_j))\hfill & \text{at the points }P_j\text{ with }j=1,\mathrm{},l\hfill \\ d(w_j^1+O(w_j))\hfill & \text{at the points }P_j\text{ with }j=l+1,\mathrm{},2l\hfill \end{array}$$

normalized by the condition $`_{a_j}𝑑\mathrm{\Omega }^k=0`$. We set $`U_j^{(k)}=\frac{1}{2\pi i}_{b_j}𝑑\mathrm{\Omega }^k`$; then, up to a normalization constant

$$\mathrm{\Psi }_\alpha (x,t,P)=\frac{\theta (A(P)+U^{(1)}x+U^{(2)}t+Z_\alpha )\theta (Z_0)_{j\alpha }\theta (A(P)+S_j)}{\theta (A(P)+Z_\alpha )\theta (U^{(1)}x+U^{(2)}t+Z_0)_{i=1}^{2l}\theta (A(P)+R_i)}e^{x\mathrm{\Omega }^{(1)}(P)+t\mathrm{\Omega }^{(2)}(P)}$$
(3.1)
where $`A(P)`$ is the Abel transform, $`K`$ is the vector of Riemann constants,

$$R_i=K\underset{s=1}{\overset{g1}{}}A(\gamma _s)A(\gamma _{g1+i}),S_i=K\underset{s=1}{\overset{g1}{}}A(\gamma _s)A(P_i),$$

$$Z_0=K\underset{i=1}{\overset{g+2l1}{}}A(\gamma _i)+\underset{j=1}{\overset{2l}{}}A(P_j),Z_\alpha =Z_0A(P_\alpha ),$$

$`\mathrm{\Omega }^{(k)}=_{P_0}^P𝑑\mathrm{\Omega }^k`$. The dual eigenfunction is constructed analogously

$$\mathrm{\Psi }_\alpha ^+(x,t,P)=\frac{\theta (A(P)U^{(1)}xU^{(2)}t+Z_\alpha ^+)\theta (Z_0^+)_{j\alpha }\theta (A(P)+S_j^+)}{\theta (A(P)+Z_\alpha ^+)\theta (U^{(1)}x+U^{(2)}tZ_0^+)_{i=1}^{2l}\theta (A(P)+R_i^+)}e^{x\mathrm{\Omega }^{(1)}(P)t\mathrm{\Omega }^{(2)}(P)}$$

where

$$S_j^+=KA(P_j)\underset{s=1}{\overset{g1}{}}A(\gamma _s^+),R_i^+=KA(\gamma _{g1+i}^+),$$

$$K_0=\underset{i=1}{\overset{g+2l1}{}}(A(\gamma _i)+A(\gamma _i^+))2\underset{j=1}{\overset{2l}{}}A(P_j),Z_0^+=Z_02KK_0,Z_\alpha ^+=Z_0^+A(P_\alpha ).$$

These functions satisfy the auxiliary linear problems (2.2), (2.5). The finite gap potential of the Dirac operator, i.e. the solution of the matrix Davey-Stewartson equation is given by $`A=\sigma F^1\sigma F^1`$, where $`\sigma `$ is the diagonal matrix with diagonal elements $`\underset{l}{\underset{}{1,\mathrm{},1}},\underset{l}{\underset{}{1,\mathrm{},1}}`$ and $`F_1`$ is the matrix coefficient of the decomposition near the singularities of the FBA. Therefore the solution of the system of poles dynamics (2.7),…,(2.12) is given by the formulas:

1. the poles of the eigenfunction can be found as the solutions of the equation

$$\theta (U^{(1)}x_i(t)+U^{(2)}t+Z_0)=0,$$
2. the spin variables are given by:

$$\stackrel{~}{a}_{i,\alpha }(t)=Q_i^1(t)\frac{\theta (U^{(1)}x_i(t)+U^{(2)}t+Z_\alpha )}{\theta (Z_\alpha )\theta (S_\alpha )},$$

$$\stackrel{~}{b}_i^\alpha (t)=Q_i^1(t)\frac{\theta (U^{(1)}x_i(t)+U^{(2)}tZ_\alpha ^+)}{\theta (Z_\alpha ^+)\theta (S_\alpha ^+)},$$
3. the normalizing factor is:

$$Q_i^2(t)=\frac{1}{2}\underset{\alpha =1}{\overset{2l}{}}\frac{\theta (U^{(1)}x_i(t)+U^{(2)}tZ_\alpha )\theta (U^{(1)}x_i(t)+U^{(2)}tZ_\alpha ^+)}{\theta (Z_\alpha )\theta (Z_\alpha ^+)\theta (S_\alpha )\theta (S_\alpha ^+)}.$$

In the formulas above we used the notation:

$$a_{i,\alpha }^1=2\stackrel{~}{a}_{i,\alpha },a_{i,\alpha }^2=2\stackrel{~}{a}_{i,\alpha +l}$$

$$b_{i}^{\alpha }{}_{}{}^{1}=\stackrel{~}{b}_i^{\alpha +l},b_{i}^{\alpha }{}_{}{}^{2}=\stackrel{~}{b}_i^\alpha .$$

## 4 Universal symplectic structure

The systems obtained as poles dynamics in the finite gap theory were often introduced from physical motivations. Such an explanation leaves the structural questions concerning dynamical systems apart. By structure we mean the algebraic symmetries in the broad sense. The importance of the hamiltonian structure cannot be overestimated in the fundamental problem of quantization, likewise the “action-angle” variables have a crucial meaning in the Seiberg-Witten theory. Heuristically most of the systems involved are hamiltonian and moreover integrable in the sense of Liouville. In  was proposed the universal hamiltonian description of such systems in terms of the flows on the space of quasi-periodic operators generated by the Lax equation. First, this approach was applied to the description of the hamiltonian structure of the nonlinear equations of the soliton theory, but in  the efficiency of this method when applied to the finite-dimensional case was confirmed. The essential point is the algebraic-geometric correspondence which maps the coordinates of the phase space of the system of poles dynamics to the space of spectral data, namely to the set of modular coordinates, the divisor of poles of the Baker-Akhiezer function and the type of singularities, which are given by a pair of Abelian differentials on the curve, in our case by the differentials $`dk,dz`$.

The coordinates on the phase space of the system (2.7),…,(2.12) given by $`x_i,\lambda _i,a_1^j,a_2^j,b_1^j,b_2^j`$ are invariant under the following change of variables

$$a_1^iW_1^1a_1^i,b_{2}^{i}{}_{}{}^{+}b_{2}^{i}{}_{}{}^{+}W_1,$$

$$a_2^iW_2^1a_2^i,b_{1}^{i}{}_{}{}^{+}b_{1}^{i}{}_{}{}^{+}W_2$$

where $`W_k,k=1,2`$ are constant matrices. We also have the symmetry of the system

$$a^i\mu _ia^i;b^i\mu _{i}^{}{}_{}{}^{1}b^i$$

where $`\mu _i`$ are constant scalars. We also recall the condition $`(b_{1}^{j}{}_{}{}^{+}a_2^j)(b_{2}^{j}{}_{}{}^{+}a_1^j)=2.`$ The dimension of the reduced space is $`4Nl+2N2l(l1)NN=4Nl2l(l1).`$ We fix $`W_1,W_2`$ such that the matrices

$$W_1^1\underset{i}{}a_1^ib_{2}^{i}{}_{}{}^{+}W_1,W_2^1\underset{i}{}a_2^ib_{1}^{i}{}_{}{}^{+}W_2$$
(4.1)
are diagonal and the matrices $`W_k`$ leave the vector $`(1,\mathrm{},1)`$ invariant. Further we use the notation $`A_i^k,B_j^m`$ for the new variables. The normalization (4.1) takes the form:

$$\underset{i}{}A_{il}^1B_{im}^2=\delta _{lm}k_l^1,$$

$$\underset{i}{}A_{il}^2B_{im}^1=\delta _{lm}k_l^2.$$

For the eigenvector of the operator $`L`$ the following representation is true:

$$c_i=(zz_1)(c_i^j+O(zz_1))\mathrm{exp}(x_i\zeta (z)),\text{ if }jl,$$

$$c_i=(zz_2)(c_i^j+O(zz_2))\mathrm{exp}(x_i\zeta (z)),\text{ if }2lj>l,$$

$$c_i=(c_i^j+O(z))\mathrm{exp}(x_i\zeta (z)),\text{ if }j>2l.$$

The vectors $`c^j`$ are proportional to the vectors $`B^j`$. We chose for normalizing the following condition

$$\underset{i=1}{\overset{N}{}}A_{im}c_i^j=\delta _{jm},j2l,\underset{i=1}{\overset{N}{}}A_{im}c_i^j=0,j>2l.$$
(4.2)
It implies the customary normalizing of FBA at the points $`P_j`$

$$\mathrm{\Psi }_q=(\delta _{jq}+O(z))\mathrm{exp}(k_jx).$$

Now we define the “action” variables as the poles of the vector $`C`$ which solves the equation $`(k+L)C=0`$ and satisfies the normalizing condition

$$\underset{m=1}{\overset{l}{}}\underset{i=1}{\overset{N}{}}A_{im}c_i\mathrm{\Phi }_1(x_i,z)+\underset{m=l+1}{\overset{2l}{}}\underset{i=1}{\overset{N}{}}A_{im}c_i\mathrm{\Phi }_2(x_i,z)=1$$
(4.3)
which is completely in agreement with the previous normalization due to Riemann-Roch theorem. Indeed, the expression on the left-hand side is a meromorphic function on the spectral curve, it has $`g+2l1`$ poles and it has fixed values $`1`$ at the points $`P_i`$, i.e. additional $`2l`$ conditions. This function has to be equal to $`1`$ identically. Now all the preliminaries are done and we turn to the construction. The Lax operator and its eigenvector normalized by the condition (4.3) are natural coordinates on the phase space and we take as $`\delta L`$ and $`\delta C`$ the corresponding cotangent vectors. We only introduce new variables $`k^{new},L^{new}`$ such that $`k^{new}=\widehat{k}2\zeta (z)+\zeta (zz_1)+\zeta (zz_2)`$ which is a well defined function on the spectral curve, new variables $`\lambda _i^{new}=\lambda _i\zeta (\eta )((B_{2}^{i}{}_{}{}^{+}A_1^i)+(B_{1}^{i}{}_{}{}^{+}A_2^i))`$ and the corresponding expression for the Lax operator:

$$L_{ij}^{new}=\delta _{ij}(\lambda _i^{new}+\zeta (zz_2)((B_{1}^{i}{}_{}{}^{+}A_2^i)1)\zeta (zz_1)((B_{2}^{i}{}_{}{}^{+}A_1^i)+1))+$$

$$(1\delta _{ij})((B_{1}^{i}{}_{}{}^{+}A_2^i)\mathrm{\Phi }_2(x_ix_j,z)(B_{2}^{i}{}_{}{}^{+}A_1^i)\mathrm{\Phi }_1(x_ix_j,z)).$$
(4.4)
The index $`new`$ is omitted in the following. We also introduce the differentials on the total space $`\delta k`$ and $`\delta z`$. To relate it with the canonical differential on the spectral curve we notice that the total space is the fibration on the specific moduli space with fiber which is the symmetric power of the spectral curve. We also notice that the sub-varieties $`z=const`$ are transversal to the fibers. So we could extend differentials on the fiber to the total space defining them zero along the sub-varieties $`z=const`$. It means that we choose the connection which is zero along such sub-varieties. For this choice of a connection we have $`\delta z=dz`$, where we use the notation $`d`$ for the extended differential (for more details see ).

The main object of our concern is the two-form on the phase space

$$\omega =\frac{1}{2}\underset{i=1}{\overset{N}{}}res_{P_i}(<C^{}(\delta L+\delta k)\delta C>)dz,$$
(4.5)
where $`C^{}`$ is the eigenvector of $`L`$ satisfying $`<C^{}C>=1`$.

###### Theorem 4

The two-form $`\omega `$ is

$$\omega =2\underset{s}{}\delta z(\gamma _s)\delta k(\gamma _s).$$
(4.6)
Proof We first clarify the meaning of this formula. $`k`$ is a well defined function on the curve and $`z`$ is a multi-valued function. After the choice of a branch of the function $`z`$, its values at the points $`\gamma _s`$ are functionals on the space of spectral data, i.e. on the phase space of our system. The only thing to justify is that the corresponding cotangent vector $`\delta z(\gamma _s)`$ does not depend on the choice of branch. Indeed, the difference between such functionals on different branches depends only on the modular parameter $`\tau `$ of the base elliptic curve which is a priori fixed. The differential $`\mathrm{\Omega }=<C^{}\delta L\delta C>dz`$ is meromorphic on the spectral curve, the essential singularities reduce. The sum of the residues at the points $`P_i`$ is the opposite of the sum of other residues. $`\mathrm{\Omega }`$ has poles at the points $`\gamma _s`$. Using the condition $`<C^{}C>=1`$ we obtain that the residues at these points are

$$res_{\gamma _s}\mathrm{\Omega }=<C^{}\delta LC>\delta z(\gamma _s)=\delta k(\gamma _s)\delta z(\gamma _s).$$

Other possible poles of $`\mathrm{\Omega }`$ are the branch points of the curve, at these points $`C^{}`$ has poles due to the normalization, but the differential $`dz`$ has simple zeroes at these points. $`\mathrm{}`$

This proposition shows that our two-form is symplectic, i.e. it is closed and non-degenerate. It is non-degenerate because the number of Darboux coordinates is equal to the dimension of the phase space.

###### Theorem 5

The symplectic form $`\omega `$ admits the representation

$$\omega =\underset{i=1}{\overset{N}{}}(\delta \lambda _i\delta x_i+\delta B_{i}^{q}{}_{}{}^{1}\delta A_{iq}^2\delta B_{i}^{q}{}_{}{}^{2}\delta A_{iq}^1).$$
(4.7)
This is proven by a straightforward calculation. We use the notation:

$$A_i=(A_{i1}^2,\mathrm{},A_{il}^2,A_{i1}^1,\mathrm{},A_{il}^1);$$

$$B_i=(B_{i1}^1,\mathrm{},B_{il}^1,B_{i1}^2,\mathrm{},B_{il}^2).$$

Let us decompose $`L,C`$ as $`L=G\stackrel{~}{L}G^1,C=G\stackrel{~}{C}`$ where $`G=Diag(e^{\zeta (z)x_i}).`$ Then

$$\omega =\frac{1}{2}res_{z=0}Tr[\delta \stackrel{~}{L}\delta h+\widehat{C}^1(\delta \stackrel{~}{L}\delta \widehat{C}+[\delta h,\stackrel{~}{L}]\delta \widehat{C}\delta h\widehat{C}\delta \widehat{k})].$$

Here the expression $`\widehat{C}`$ corresponds to the section of the direct image of the line bundle on the spectral curve to the base elliptic curve and $`\delta h=G^1\delta G=diag(\delta x_i\zeta (z)).`$ Using

$$Tr(C^1[\delta h,L]\delta C)=Tr(C^1\delta h(\delta L\widehat{C}\widehat{C}\delta \widehat{k})),$$

we obtain

$$\omega =\frac{1}{2}res_{z=0}Tr(2\delta \stackrel{~}{L}\delta h+\widehat{C}^1\delta L\delta C).$$

The first summand gives $`_i\delta \lambda _i\delta x_i.`$ Calculating the second summand we take into account the normalizing conditions (4.2)

$$\underset{j=1}{\overset{N}{}}(A_{j\alpha }\delta c_j^k+\delta A_{j\alpha }c_j^k)=0,\underset{j=1}{\overset{N}{}}(c_{kj}^{}\delta B_{j\alpha }+\delta c_{kj}^{}B_j^\alpha )=0,$$

which demonstrates the proposition. $`\mathrm{}`$

###### Theorem 6

The system of poles dynamics of the finite gap potentials of the Davey-Stewartson equation (2.7),…,(2.12) is hamiltonian with respect to the symplectic form $`\omega `$ and the Gaudin Hamiltonian

$$H=\frac{1}{2}\underset{i}{}\lambda _i(B_{2}^{i}{}_{}{}^{+}A_1^i+B_{1}^{i}{}_{}{}^{+}A_2^i)\underset{i}{}\zeta (\eta )(B_{2}^{i}{}_{}{}^{+}A_1^i+1)(B_{1}^{i}{}_{}{}^{+}A_2^i1)+$$

$$\underset{ij}{}(B_{1}^{i}{}_{}{}^{+}A_2^j)(B_{2}^{i}{}_{}{}^{+}A_1^i)\frac{\sigma (x_jx_i+\eta )}{\sigma (\eta )\sigma (x_jx_i)}$$
(4.8)
is equal to $`(H_1H_2)/2`$, defined in section 1 for the Gaudin elliptic system.

Proof We substitute the flow defining the dynamics, i.e. the Lax equation, into the symplectic form. Our intention is to show that $`i__t\omega (X)=\omega (X,_t)=dH(X).`$ With the chosen coordinates the vector field is defined by the expressions:

$$_t(L_1L_2)=[L_1,L_2],_tk=0;$$

$$_tC(t,P)=L_2C(t,P)+\mu (t,P)C(t,P),$$

where $`\mu (t,P)`$ is a scalar function.

$$i__t\omega =\frac{1}{2}\underset{i=1}{\overset{2l}{}}res_{P_i}(<C^{}(\delta L\delta k)(L_2+\mu (t,P))C><C^{}[L_2,L_1]\delta C>)dz=$$

$$=\underset{i=1}{\overset{2l}{}}res_{P_i}\delta k\mu (t,P)dz.$$
(4.9)
The function $`\mu `$ has the same expansion at the marked points as the coefficient at $`t`$ in the exponential part of the vector $`\mathrm{\Psi }`$. It means that at $`P_i,i=1\mathrm{}l`$ the principal part of $`\mu `$ is equal to $`\frac{1}{2}k`$ and at the points $`P_i,i=l+1,\mathrm{},2l,\mu \frac{1}{2}k`$.

$$i__t\omega =\frac{1}{4}\underset{i=1}{\overset{l}{}}res_{P_i}\delta k^2dz+\frac{1}{4}\underset{i=l+1}{\overset{2l}{}}res_{P_i}\delta k^2dz=$$

$$=\frac{1}{4}\underset{i=1}{\overset{l}{}}res_{P_i}\delta <C^{}L^2C>dz+\frac{1}{4}\underset{i=l+1}{\overset{2l}{}}res_{P_i}\delta <C^{}L^2C>dz=$$

$$=\frac{1}{4}res_{z=z_2}Tr\delta L^2dz\frac{1}{4}res_{z=z_1}Tr\delta L^2dz.$$

Substituting the expressions for the matrix elements of the Lax operator (4.4) into this formula we find

$$res_{z=z_2}TrL^2dzres_{z=z_1}TrL^2dz=$$

$$\frac{1}{2}\underset{i}{}(\lambda _i((B_{2}^{i}{}_{}{}^{+}A_1^i)+(B_{1}^{i}{}_{}{}^{+}A_2^i))+\zeta (\eta )((B_{2}^{i}{}_{}{}^{+}A_1^i)^2+(B_{1}^{i}{}_{}{}^{+}A_2^i)^2))+$$

$$\underset{i}{}\zeta (\eta )((B_{2}^{i}{}_{}{}^{+}A_1^i)(B_{1}^{i}{}_{}{}^{+}A_2^i)+1)+\underset{ij}{}(B_{1}^{i}{}_{}{}^{+}A_2^j)(B_{2}^{i}{}_{}{}^{+}A_1^i)\frac{\sigma (x_jx_i+\eta )}{\sigma (\eta )\sigma (x_jx_i)}.$$

We have found the Hamiltonian up to a constant. To obtain the Hamiltonian $`H_1H_2`$ of the Gaudin system we have to subtract $`N\zeta (\eta )`$ and use the condition $`(B_{1}^{i}{}_{}{}^{+}A_2^i)(B_{2}^{i}{}_{}{}^{+}A_1^i)=2`$. We only need to note that the system (2.7),…,(2.12) is equivalent to the system of the hamiltonian dynamics with this Hamiltonian, and the equivalence is given by the formulas:

$$a_k^i\mathrm{exp}((\lambda _1^i+\lambda _2^i)/2)a_k^i;$$

$$b_k^i\mathrm{exp}((\lambda _1^i+\lambda _2^i)/2)b_k^i.$$

$`\mathrm{}`$

The “action-angle” variables can be constructed in the traditional way as in . We introduce the differentials $`d\mathrm{\Omega }_j,j=1,\mathrm{},2l1`$ of the third kind with residues $`1`$ and $`1`$ at the singular points $`P_j`$ and $`P_{2l}`$ respectively. On the fundamental domain of our curve we can define the functions $`A_k(Q),k=g+2l1`$ such that

$$A_k(Q)=^Q𝑑\omega _i,i=1,\mathrm{},g,A_{g+i}(Q)=^Q𝑑\mathrm{\Omega }_i,i=1,\mathrm{},2l1.$$

Now we introduce the quantities

$$\varphi _i=\underset{s=1}{\overset{g+2l1}{}}A_i(\gamma _s),I_k=_{a_i}k𝑑z,I_{g+i}=res_{P_i}kdz.$$

The symplectic form is

$$\omega =\underset{i=1}{\overset{g+2l1}{}}\delta \varphi _i\delta I_i.$$

Notice that the quantities $`I_i`$ are integrals of motion. So, the obtained coordinates are the “action-angle” variables for our system. For the $`sl(2)`$-case the separated variables was also realized in .

## 5 Degenerations

The construction of the Hitchin system on degenerate curves is well suited for the analysis of several degenerations of dynamical systems. For example the rational and trigonometric degeneration can be obtained by considering the elliptic curve with cuspidal and nodal singularities respectively instead of the regular elliptic curve. The normalization in both cases is the rational curve and the moduli of holomorphic bundles are parameterized only by the elements gluing the fibers over the singular points. For the Gaudin elliptic system such degenerations exist. The rational one is represented by the hamiltonian

$$H_{rat}=\frac{1}{2}\underset{i}{}\lambda _i(B_{2}^{i}{}_{}{}^{+}A_1^i+B_{1}^{i}{}_{}{}^{+}A_2^i)+\underset{ij}{}(B_{1}^{i}{}_{}{}^{+}A_2^j)(B_{2}^{i}{}_{}{}^{+}A_1^i)\left(\frac{1}{x_ix_j}\frac{1}{\eta }\right)$$
(5.1)
with the canonical symplectic form $`\omega `$ in the representation (4.7). This system involves a rational curve with two simple points and one double point with cuspidal singularity. It also could be solved by the algebraic-geometric direct and inverse problem as it was done in , . The expressions for the $`\mathrm{\Psi }`$ function could be found in  (2.53). The “action-angle” variables for it were also constructed in . The author deals with another Hamiltonian but the separated variables for these systems coincide. There is another type of degenerations of our original system. It is the limit $`z_1z_2,\eta 0.`$ This case is successfully taken care of by the direct and inverse problem as well. The FBA function can be found in  (2.60). The solution provides the double periodic potentials for the Davey-Stewartson equation. The system can be specified by the Hamiltonian

$$H_{ell}=\frac{1}{2}\underset{i}{}\lambda _i(B_{2}^{i}{}_{}{}^{+}A_1^i+B_{1}^{i}{}_{}{}^{+}A_2^i)+\underset{ij}{}(B_{1}^{i}{}_{}{}^{+}A_2^j)(B_{2}^{i}{}_{}{}^{+}A_1^i)\zeta (x_ix_j)$$
(5.2)
and the standard symplectic form (4.7). A priori this system is not integrable. Proceeding with the analogous construction for the “action-angle” variables, we obtain that their number is $`2l(2N2l+1)`$ which is less than the dimension of the phase space. The same degeneracy can be observed by counting the number of Hitchin Hamiltonians. The Hamiltonians of the original system can be represented as $`H_s^{k,l}=res_{z=z_s}((zz_s)^{k1}Tr\mathrm{\Phi }^l),s=1,2.`$ When the two fixed points coincide, there is no more residues of the expressions $`Tr\mathrm{\Phi }^k`$, and the only retrievable quantities are $`H_1^{k,l}+H_2^{k,l}.`$ This case should be investigated with a more delicate geometrical analysis of the jet structure of the holomorphic bundle.

## Conclusion

The main result of this paper is the explicit solution of the Gaudin elliptic system with spin and the uncovered connection between this system and the quasi periodic solutions of the Davey-Stewartson equation. For the methodological aspects we noticed the crucial correspondence between the Hitchin and finite gap descriptions of the dynamical systems. Our subsequent concern is to give a more general illustration of the equivalence of these two settings. The important technical query in this direction is the geometry of the moduli space of semi-stable bundles when we vary the base algebraic curve. It could clarify the Hitchin description of such a system as an elliptic degenerate Gaudin model. Such an attempt was undertaken and the case of the elliptic curve was analyzed in . The construction crucially used the group structure on the elliptic curve and could not be translated directly to the general situation. In order to realize the geometrically established limit procedure a kind of canonical connection on the bundle over the moduli space of algebraic curves with moduli space of holomorphic bundles as the fiber is required. And the role of the Painlevé connection in this context is of immediate concern to us.

Acknowledgements

I would like to thank I. Krichever for introducing me to the beauty of integrable systems and the special spectral method and for his proposal to investigate the special spectral problem for the Davey-Stewartson equation. I would like to thank M. Olshanetsky and A. Levin for fruitful discussions on the nature of the Hitchin system. This work was partially supported by the RFBR 01-01-00546 grant.


text:
# Molecular Gas in Infrared-Excess, Optically-Selected QSOs and the Connection with Infrared Luminous Galaxies

## 1. Introduction

There exists substantial evidence that the bulges of nearby massive galaxies contain quiescent supermassive nuclear black holes ($`M_{}=10^{69}`$ M: Magorrian et al. 1998; van der Marel 1999; Kormendy et al. 1998). The almost ubiquitous presence of such objects, in combination with the possibility that nuclear black hole mass scales with stellar bulge mass (i.e., the Magorrian Relation: Magorrian et al. 1998), indicates that supermassive black holes are an essential part of spheroid formation, and that mass accretion onto nuclear black holes during stellar bulge formation likely resulted in periods of AGN activity.

There also exists growing observational evidence that core111We define the “core” to be the region of a galaxy internal to the radius at which the projected surface brightness is half its central brightness. This differs from the “bulge” of a galaxy, which is the region internal to the radius which encloses half the total luminosity of a galaxy. formation and AGN activity may be linked with ultraluminous infrared galaxy (ULIG)222ULIGs are defined as galaxies with $`L_{\mathrm{IR}}[81000\mu \mathrm{m}]10^{12}`$ L. mergers, lending support to the evolutionary connection between ULIGs and QSOs first proposed by Sanders et al. (1988a). The evidence can be summarized as follows: First, observations of massive concentrations of molecular gas in the nuclear regions of local ULIGs relative to that observed in normal galaxies (e.g. Kormendy & Sanders 1992; Downes & Solomon 1998; Bryant & Scoville 1999) support the highly dissipative nature of ULIGs. Second, spectroscopic evidence that a significant fraction of ULIGs, most of which have “warm” infrared colors (i.e., 25$`\mu `$m to 60$`\mu `$m flux density ratio, $`f_{25}/f_{60}0.20`$), have Seyfert-like emission line spectra (Veilleux et al. 1995; Veilleux, Kim, & Sanders 1999) may indicate that merger events are intimately connected with AGN activity. Third, recent CO($`10`$) observations show the association of molecular gas with the AGN nucleus of several double nucleus ULIGs (Evans et al. 1999; Evans et al. 2000), supporting the necessity of dissipative material for fueling AGN activity. Fourth, the recent evidence of similarities in the space densities of high-redshift submillimeter (SCUBA) galaxies (e.g., Smail, Ivison, & Blain 1997; Barger et al. 1998; Hughes et al. 1998; Barger, Cowie, & Sanders 1999) and local giant elliptical galaxies has been used to infer that ULIGs are connected with massive galaxy formation (Barger et al. 1999; Trentham 2000).333This assumes that high-redshift infrared luminous galaxies are distant analogs of local infrared galaxies.

From the above, it seems reasonable to assume that a significant fraction of QSOs may have molecular gas rich hosts and that molecular gas may play an important role in fueling not only circumnuclear starbursts, but also AGN activity in galaxies experiencing dissipative collapse. As a first step towards establishing the presence of abundant reservoirs of molecular gas in QSOs, a survey of CO($`10`$) emission in a complete sample of low-redshift, optically-selected QSOs has been initiated. In addition to detecting the possible fueling source for both AGN and star formation activity, these observations provide a major missing link in comparing the properties of low-redshifts QSOs and infrared luminous galaxies with both “cool” (i.e., 25$`\mu `$m to 60$`\mu `$m flux density ratio, $`f_{25}/f_{60}<0.20`$) and “warm”, Seyfert-like ($`f_{25}/f_{60}0.2`$) mid-to-far infrared flux density ratios. Specifically, these data make it possible to determine whether the relative molecular gas properties of PG QSOs and ULIGs are consistent with the hypothesis that cool and warm ULIGs represent successive stages in the evolution toward UV-excess QSOs (Sanders et al. 1988a). The CO data presented here are complemented by previous studies of the far-infrared properties (Sanders et al. 1989a) and host galaxy morphologies (Bahcall et al. 1997; Surace 1998; Surace & Sanders 2000: hereafter SS00) of optical QSOs, and thus facilitate a more well-rounded comparison with the wealth of data available of low-redshift ($`z<0.3`$) infrared luminous galaxies.

This paper is divided into 6 sections. Section 2 is a summary of the sample selection. In §3, the general observing procedure is discussed, followed by the results in §4. A comparison of the CO and infrared properties of the PG QSOs with those of local infrared luminous galaxies and high redshift QSOs is provided in §5, along with a summary of molecular gas and infrared properties of the PG QSOs relative to the optical morphologies of their host galaxies. Finally, §6 summarizes the paper. Throughout this paper, $`H_0=75`$ km s-1 Mpc-1 and $`q_0=0.5`$ are assumed.

## 2. Sample Selection

The QSOs observed as part of this CO study were selected from a complete sample of 18 QSOs defined by SS00. The reader is referred to SS00 for a detailed discussion of the sample selection, but it is briefly summarized here for completeness. The selection criteria consisted of all $`0.04<z<0.17`$ QSOs from the Palomar–Green Bright Quasar Survey (PG BGS: Schmidt & Green 1983) having infrared-to-“big blue bump” luminosities, $`L_{\mathrm{IR}}(81000\mathrm{\mu m})/L_{\mathrm{bbb}}(0.11.0\mathrm{\mu m})>0.36`$ and absolute blue magnitudes, $`M_\mathrm{B}22.0`$. The former criteria selects QSOs with infrared excesses that overlap with warm infrared luminous galaxies (see §5.2).

There are two advantages to using the SS00 sample. First, the far-infrared flux densities of these QSOs have been measured using the Infrared Astronomical Satellite (IRAS) data (Sanders et al. 1989a), making direct comparisons of their infrared and molecular gas properties with those of infrared luminous galaxies possible. Second, the molecular gas properties of these QSOs can be compared with pre-existing optical images of the QSO host galaxies to search for correlations of $`L_{\mathrm{CO}}^{}`$ with the morphology of the host galaxies.

While QSOs with infrared excesses may at first seem to be unrepresentative of optical QSOs, it must be remembered that at least two-thirds of all the PG QSOs in the Bright Quasar Survey show clear evidence of infrared emission consistent with thermal dust emission, with the infrared components contributing 20-40% of their total bolometric luminosity (Sanders et al. 1989a). Thus, the sample selected for this survey is simply the optical QSOs with the highest probable dust content.

Of the 18 QSOs in the SS00 sample, two of them (IZw1 and Mrk 1014) have previous confirmed CO($`10`$) detections (e.g. Sanders, Scoville, & Soifer 1988; Barvainis et al. 1989; Solomon et al. 1997; Schinnerer, Eckart, & Tacconi 1998) and two others (PG 0838+770 and PG 1613+658) have unconfirmed CO detections (Alloin et al. 1992). In addition to observing PG 0838+770 and PG 1613+658, eight PG QSOs with $`L_{\mathrm{IR}}/L_{\mathrm{bbb}}`$ ratios in the range of 0.36–1.1 were selected for CO observations. The sample is listed in Table 1.

## 3. Observations

Aperture synthesis maps of CO($`10`$) emission in the sample of PG QSOs were made with the Owens Valley Radio Observatory (OVRO) Millimeter Array during 21 observing periods from 1999 April to 2000 October. The array consists of six 10.4m telescopes, and the longest observed baseline was 103m. Each telescope was configured with 120$`\times `$4 MHz digital correlators. All but one of the observations were done in the low-resolution configuration, providing a $`4.0\mathrm{}`$ (full width at half the maximum intensity, or FWHM) synthesized beam with natural weighting. For each QSO, the receivers were tuned to the expected frequency of the redshifted CO emission line. These redshifts were estimated from published optical emission line redshifts compiled by the NASA Extragalactic Database (NED); typically, the average or the most common redshift was adopted. During the observations, nearby quasars were observed every 25 minutes to monitor phase and gain variations, and several brighter quasars were observed to determine the passband structure. Observations of Uranus were also made for absolute flux calibration; typical uncertainties in the flux calibration are 15-20%. A journal of observations is provided in Table 2.

## 4. Results

The data were reduced and calibrated using the standard Owens Valley array program MMA (Scoville et al. 1992). The data were then exported to the mapping program DIFMAP (Shepherd, Pearson, & Taylor 1995) and to the NRAO software package AIPS for spectral line extraction.

Figure 1 shows the integrated CO($`10`$) intensity maps for the detected PG QSOs (integrated intensity values are listed in Table 3). Six out of the ten QSOs observed were detected in CO; all of the detections have a signal-to-noise ratio, S/N, greater than 4.5$`\sigma `$. A comparison was done between the CO beam centroid positions and the optical QSO positions, the latter of which were extracted from the USNO-A2.0 Catalog. With the exception of PG 1351+640 (1.44″ offset) and PG 1613+658 (1.05″ offset), all the CO centroid and optical positions agreed to within 0.9″. Of the four non-detections, two are of the QSOs with the highest redshifts in the sample (PG 1402+261 and PG 1202+281). For all four non-detections, the 1.0 GHz continuum channel data was examined to search for CO emission over a broader (albeit less sensitive) bandwidth than covered by the spectrometer. Such searches yielded null results.

In the present sample of ten QSOs, PG 0007+106 is the only one having detectable 3mm continuum emission. The unresolved continuum emission is observed to be at RA = 00:10:30.98 dec = +10:58:29.8 (J2000.0), which is within 0.75$`\mathrm{}`$ of the USNO-A2.0 optical position, and the continuum flux density is measured to be 0.1 Jy.

Figure 2 shows the CO($`10`$) spectra of the six detected QSOs. The CO emission line widths are varied, ranging from FWHM of $`\mathrm{\Delta }v_{\mathrm{FWHM}}5090`$ km s-1 for PG 1415+451 and PG 0838+770 to $`500`$ km s-1 for PG 1613+658. All of CO emission line centroids are within 300 km s-1 of the adopted optical redshifts, indicating that the optical emission lines are only moderately red/blueshifted relative to the systemic velocity. Both PG 0838+770 and PG 1613+658 are detected in CO; the emission line widths and redshifts are similar to those reported by Alloin et al. (1992), however, the integrated CO line intensity of PG 0838+770 is only 40% of that reported by Alloin et al. (1992). By comparison, the CO line intensity of PG 1613+658 agrees with the Alloin et al. (1992) measurement to within 12% (i.e. within the error of the OVRO measurements).

While the CO emission line of PG 1415+451 was detected with a signal-to-noise of 5.2$`\sigma `$, the narrow appearance of the emission feature may be an indication that the “line” is actually an artifact; i.e., a spurious interference spike created at the phase center of the interferometer. We believe that the emission feature is authentic based on a combination of the following three facts: First, the emission feature spans five independent (i.e., unbinned) channels. Second, while the phase center of the PG 1415+451 observations is within the beam of the emission feature, and the centroid position of the feature is $`0.9\mathrm{}`$ closer to the optical QSO position (obtained from the USNO–A2.0 catalog) than the phase center is to the optical QSO position. Third, the CO observation of PG 0838+770, which has a CO line width less than a factor of two wider than the emission feature observed for PG 1415+451, was a confirmation of prior observations reported by Alloin et al. (1992), and is thus authentic. It is thus likely that the host galaxies of both QSOs simply have molecular disks that are being observed nearly face-on.

### 4.1. CO Line Luminosities and Molecular Gas Masses

Table 3 summarizes the emission line properties of the ten PG QSOs. The luminosity distance for a source at a given redshift is,

$$D_\mathrm{L}=cH_0^1q_0^2\left\{zq_0+(q_01)\left(\sqrt{2q_0z+1}1\right)\right\}$$

$$[\mathrm{Mpc}].$$
$`(1)`$
Given the measured flux, $`S_{\mathrm{CO}}\mathrm{\Delta }v`$, the CO luminosity of a source at redshift $`z`$ is,

$$L_{\mathrm{CO}}^{}=\left(\frac{c^2}{2k\nu _{\mathrm{obs}}^2}\right)S_{\mathrm{CO}}\mathrm{\Delta }vD_\mathrm{L}^2(1+z)^3,$$
$`(2)`$
where $`c`$ is the speed of light, $`k`$ is the Boltzmann constant, and $`\nu _{\mathrm{obs}}`$ is the observed frequency. In terms of useful units, $`L_{\mathrm{CO}}^{}`$ can be written as,

$$L_{\mathrm{CO}}^{}=1.35\times 10^3\left(\frac{S_{\mathrm{CO}}\mathrm{\Delta }v}{\mathrm{Jy}\mathrm{km}\mathrm{s}^1}\right)\left(\frac{D_\mathrm{L}^2}{\mathrm{Mpc}^2}\right)J^2(1+z)^1$$

$$[\mathrm{K}\mathrm{km}\mathrm{s}^1\mathrm{pc}^2].$$
$`(3)`$
To calculate the mass of molecular gas in the PG QSOs, the reasonable assumption is made that the CO emission is optically thick and thermalized (see Discussion Section), and that it originates in gravitationally bound molecular clouds. Thus, the ratio of the H2 mass and the CO luminosity is given by,

$$\alpha =\frac{M(\mathrm{H}_2)}{L_{\mathrm{CO}}^{}}\frac{\sqrt{n(\mathrm{H}_2)}}{T_\mathrm{b}}\mathrm{M}_{}[\mathrm{K}\mathrm{km}\mathrm{s}^1\mathrm{pc}^2]^1,$$
$`(4)`$
where $`n(`$H$`{}_{2}{}^{})`$ and $`T_\mathrm{b}`$ are the density of H2 and brightness temperature for the CO(1$``$0) transition (Scoville & Sanders 1987; Solomon, Downes, & Radford 1992). Multitransition CO surveys of molecular clouds in the Milky Way (e.g. Sanders et al. 1993), and in nearby starburst galaxies (e.g. Güsten et al. 1993) have shown that hotter clouds tend of be denser such that the density and temperature dependencies tend to cancel each other. The variation in the value of $`\alpha `$ is approximately a factor of 2 for a wide range of kinetic temperatures, gas densities, and CO abundance (e.g. $`\alpha =25M_{}`$ \[K km s-1 pc$`{}_{}{}^{2}]^1`$: Radford, Solomon, & Downes 1991). More recent dynamical mass estimates of a low-redshift infrared galaxy sample observed in CO with the Plateau de Bure Interferometer indicate that $`\alpha `$ can be as low as 1 $`M_{}`$ \[K km s-1 pc$`{}_{}{}^{2}]^1`$ (Downes & Solomon 1998). We adopt a value of 4 $`M_{}`$ (K km s-1 pc$`{}_{}{}^{2})^1`$ for $`\alpha `$, which is similar to the value determined for the bulk of the molecular gas in the disk of the Milky Way (Scoville & Sanders 1987; Strong et al. 1988).

In the four cases where only upper limits to the CO flux were determined, the average FWHM velocity of all the PG QSOs detected in CO to date was adopted (280 km s-1). Conservative upper limits to the CO luminosity and molecular gas mass were then calculated from the 3$`\sigma `$ upper limits to the CO flux density.

## 5. Discussion

The detection of CO emission in the host galaxies of these PG QSOs provides direct evidence that these infrared-excess QSOs live in host galaxies containing significant supplies of molecular gas. It also supports the hypothesis that galaxies with a molecular gas-rich interstellar medium (ISM) may be capable of igniting and sustaining both AGN and star formation activity in QSO host galaxies. The measured CO luminosities for the six detected QSOs in Table 3 translate into molecular gas masses of $`M(\mathrm{H}_2)220\times 10^9`$ M, which, along with the additional two QSOs IZw1 and Mrk 1014, span the full range of masses derived for luminous infrared galaxies (LIGs)444LIGs are defined as galaxies having $`L_{\mathrm{IR}}[81000\mu \mathrm{m}]=10^{11.011.99}`$ L. and ULIGs.

These six PG QSOs are among the faintest CO($`10`$) sources yet detected at redshifts $`z<0.3`$. An interesting issue is whether or not the CO luminosities of these sources, and perhaps the CO luminosities of other PG QSOs as well, could have been predicted apriori. Consider the plot of CO flux versus the 100$`\mu `$m flux density of the PG QSO sample, in addition to a comparative sample of $`z<0.3`$ warm and cool ULIGs observed in CO to date (Figure 3). Solomon et al. (1997) use the tightness of this correlation, in combination with the small sizes of the CO emitting regions in ULIGs, to argue that the dust emission at 100$`\mu `$m and the CO(1$``$0) emission is optically thick (see their Figure 4). The $`4\mathrm{}`$ resolution of the present CO dataset is insufficient to determine how compact the CO emitting regions in these PG QSOs are (e.g., $`4\mathrm{}`$ corresponds to $`6`$ kpc at $`z0.1`$). For our purposes, the fact that all of the PG QSOs detected in CO to date also fall along this correlation indicates that the CO fluxes (and thus, the molecular gas masses) of optical QSOs may potentially be estimated from the 100$`\mu `$m flux density, $`f_{100}`$, to within a factor of a few if $`f_{100}`$ is known to high accuracy. Note, however, that the error in $`f_{100}`$ of any IRAS measurements exceeds 25% below 0.5 Jy; improved far-infrared measurements of QSOs will be possible via observations with the Space Infrared Telescope Facility (SIRTF).

With the detection of molecular gas in the six of the ten optical QSOs listed in Table 3, plus the previous detection of IZw1 and Mrk 1014, physically meaningful comparisons between QSOs and infrared-selected galaxies can be better addressed. The remainder of the discussion section will thus be focussed on comparing the infrared and CO luminosities of the eight QSOs detected in CO to date with infrared luminous galaxies, followed by a brief discussion of recent images of the host galaxies of these PG QSOs (Bahcall et al. 1997; SS00). The four QSOs with CO non-detections will be discussed where applicable.

### 5.1. Comparison with Infrared Luminous Galaxies

Figure 4 is a plot of CO luminosity versus redshift of the 12 QSOs. In order to examine how the CO luminosity of these PG QSOs compares with that of the most CO luminous galaxies presently known at any given redshift, the infrared luminous galaxies observed in CO($`10`$) by the Sanders et al. (1989b), Sanders, Scoville, & Soifer (1991), and Solomon et al. (1997) surveys have also been plotted. All but eight of the infrared luminous galaxies at $`z>0.04`$ are ULIGs. Three features of Figure 4 are readily apparent. First, the CO luminosity range of all twelve QSOs spans more than an order of magnitude. This range is significantly more than that of the cool ULIGs, all of which have $`L_{\mathrm{CO}}^{}10^{10}`$ K km s-1 pc2 over the same redshift range. Second, the three QSOs with the highest $`L_{\mathrm{CO}}^{}`$ (IZw1, PG 1613+658, and Mrk 1014) fall within the scatter of the cool ULIGs distribution, whereas the remaining nine QSOs form a distinct group with lower values of CO luminosity for their given redshift. Third, seven of the eight QSOs detected in CO have a range in $`L_{\mathrm{CO}}^{}`$ similar to that of the warm ULIGs. Thus, while the $`L_{\mathrm{CO}}^{}`$ values of the QSOs are diverse, there is also evidence of overlap in $`L_{\mathrm{CO}}^{}`$ between these optically-selected and infrared-selected samples.

Consider now the infrared luminosity of both samples relative to the CO luminosity. Figure 5 shows two ways in which infrared and CO data of infrared luminous galaxies are commonly plotted. The top panel (Figure 5a) is a plot of $`L_{\mathrm{IR}}`$ versus $`L_{\mathrm{CO}}^{}`$, which shows the relationship between the fuel available for star formation and AGN activity ($`L_{\mathrm{CO}}^{}`$) and the luminosity of the imbedded sources ($`L_{\mathrm{IR}}`$). The bottom panel (Figure 5b) is a plot of $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ versus $`L_{\mathrm{IR}}`$. In starburst galaxies, $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ is commonly referred to as the star formation efficiency, it is a measure of the number of high-mass stars (as traced by $`L_{\mathrm{IR}}`$) created per unit molecular gas mass (as traced by $`L_{\mathrm{CO}}^{}`$). An examination of Figure 5a shows that the quantity $`L_{\mathrm{IR}}`$ increases as $`L_{\mathrm{CO}}^{}`$ increases, and for any given CO luminosity bin, the QSOs are in the upper half of the distribution. In Figure 5b, the QSOs populate the upper half of the LIG-ULIG distribution of $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$. Note also that, in general, the warm infrared galaxies and QSOs populate similar areas of both plots.

The interpretation of Figure 5 depends on whether the infrared emission in QSOs is due to dust heating by AGN or massive stars. First, consider the former possibility: Sanders et al. (1988a) proposed that ULIGs are dust-enshrouded QSOs; in this scenario, the optical QSO is revealed after most of the molecular gas is consumed by star formation and AGN activity, and after most of the dust is destroyed/blown out by winds associated with supernovae and AGN. During this process, ULIGs pass from a cool infrared phase where vigorous star formation and black hole creation/accretion has commenced, through a warm infrared phase (see also Sanders 1988b). This idea is further addressed by Sanders et al. (1989a), who propose that the infrared emission from PG QSOs is thermal emission from dust that has absorbed 20–40% of the emission emanating from the QSO (i.e., $`L_{\mathrm{IR}}=[0.20.4]L_{\mathrm{bol}}`$, where $`L_{\mathrm{bol}}`$ is the bolometric luminosity). Given this, the quantity $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ is not a measure of the star formation efficiency in QSOs, but is an indication that the AGN is contributing significantly to heating the dust in QSOs (see also Evans et al. 1998).

If the dust-enshrouded quasar model of ULIGs is correct and the ratio $`L_{\mathrm{IR}}/L_{\mathrm{bol}}`$ for both ULIGs and QSOs is due to covering factor and not due to orientation effects (e.g. Barthel 1989), optically-selected QSOs must have had a covering factor of unity in their past ULIG phase. Consider now where the QSOs would be on Figure 5a and 5b if the QSO nuclei were completely enshrouded in dust (i.e. $`L_{\mathrm{bol}}=L_{\mathrm{IR}}`$). This increases $`\mathrm{log}L_{\mathrm{IR}}`$ of the QSOs by $``$0.3–0.5 in Figure 5a, and shows to first order a possible range in locations of enshrouded QSOs on this diagram. In Figure 5b, the quantity $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ of the QSOs is increased to extreme values ($`4001600`$). Note that, in both panels, ULIGs with high ($`>160`$) $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ occupy the region covered by the QSOs, similar to what would be expected if these ULIGs were dust-enshrouded QSOs.

As mentioned above, the second possible interpretation of the infrared emission in QSOs is that $`L_{\mathrm{IR}}`$ is primarily due to dust heating by massive stars. Thus, high values of $`L_{\mathrm{IR}}`$ (Figure 5a) and $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ (Figure 5b) would result from an extremely efficient production of high-mass stars per unit molecular gas mass in QSOs relative to most LIGs and ULIGs. The attraction of this scenario is that both stars forming in the central regions of galaxies and AGN contribute significantly to the overall energy budget of their host galaxies, and that AGN activity is somehow strongly dependent on an accompanying starburst. Indeed, observed correlations between the masses of quiescent supermassive nuclear black holes and both the stellar bulge masses (e.g., Magorrian et al. 1998) and velocity dispersions (Ferrarese & Merritt 2000; Gebhardt et al. 2000) of nearby normal galaxies suggests that the accretion of mass onto nuclear black holes may accompany stellar bulge formation. However, whether or not a significant fraction of the stellar bulge population of the present sample of PG QSOs is in the process of being formed cannot be addressed with the present data set.

### 5.2. The Host Galaxies

An obvious issue to consider at this point is whether or not the diversity of $`L_{\mathrm{CO}}^{}`$ and $`L_{\mathrm{IR}}`$ of these 12 QSOs is at all correlated with the galaxy type or evolutionary state of their host galaxies. An examination of optical and near-infrared images of all 12 QSOs (Bahcall et al. 1997; Surace et al. 1998; SS00) reveals that two of the QSOs with the highest $`L_{\mathrm{CO}}^{}`$ (Mrk 1014 and PG 1613+658) have irregular morphologies consistent with an advanced merger of two molecular gas-rich disk galaxies (see simulations by Barnes & Hernquist 1996), and the other (IZw1) has asymmetric spiral arms, suggestive of a strong interaction with a much smaller mass companion. The remaining nine QSOs have host galaxies that appear to be a mixture of barred spiral galaxies, elliptical galaxies, galaxies with indeterminant features, and at least one ongoing merger (i.e., PG 0007+106).

The morphologies are suggestive of a diversity in paths that may ultimately lead to QSO activity (SS00). First, the QSOs having the highest values of $`L_{\mathrm{CO}}^{}`$ appear to be ongoing major mergers (collisions of two or more molecular gas-rich disk galaxies where prominent tidal features are still visible). Thus, in terms of $`L_{\mathrm{CO}}^{}`$, $`L_{\mathrm{IR}}`$, and host morphologies, these QSOs are similar to ULIGs. Second, the less CO-luminous QSOs with barred spiral galaxy hosts have a less clear origin, however, the bar likely acts as the means by which molecular gas in the spiral arms is torqued and drained of enough angular momentum to agglomerate in the nuclear regions (e.g. Sakamoto et al. 1999). These barred spiral galaxies may have lower $`L_{\mathrm{CO}}^{}`$ because only one, massive molecular gas-rich galaxy is involved (though the bar instability may be triggered via an interaction with a satellite galaxy). Third, the QSOs having host galaxies with an indeterminant classification may be spiral galaxies with features too faint to discern with the sensitivity of the SS00 observations, fully evolved elliptical galaxies, or well-evolved major mergers which are evolving into elliptical galaxies (see simulations by Barnes & Hernquist 1996). Thus, these host galaxies may have low $`L_{\mathrm{CO}}^{}`$ because the host galaxy is intrinsically molecular gas poor, or because the molecular gas has already been consumed via both circumnuclear star formation and the QSO. Finally, the QSOs with elliptical host galaxies may be the end stage of the merger one two disk galaxies, or elliptical host galaxies that have cannibalized a satellite galaxy (such as the $`z0.44`$ galaxy IRAS PSC 09104+4109: e.g. Evans et al. 1998).

### 5.3. Comparison with High Redshift QSOs

A significant number of QSOs have been detected at high redshift ($`z>2`$: Barvainis et al. 1994, 1998; Ohta et al. 1996; Omont et al. 1996; Guilloteau et al. 1997; Downes et al. 1999). These QSOs, many of which were initially selected for CO observations based on their strong rest-frame infrared luminosity, are likely among the most molecular gas rich QSOs at their redshift. Figure 6 is a plot of the CO luminosity of the eight PG QSOs detected in CO to date, as well as all of the moderate ($`0.3<z2`$) and high-redshift CO luminous QSOs, versus the redshift. Several of the QSOs at high redshift are known to be gravitationally lensed; in these cases, both observed and the intrinsic CO luminosities are plotted. Based on the present sample, there is no convincing evidence for evolution of molecular gas mass of the most CO luminous QSOs at any given redshift as a function of increasing redshift. The ratio of $`L_{\mathrm{CO}}^{}`$ of BRI 1335-0415 ($`z=4.7`$) to that of the most CO luminous low-redshift QSO (Mrk 1014: $`z=0.16`$) is three. However, given the uncertainty in the $`L_{\mathrm{CO}}^{}`$–to–$`M(`$H$`{}_{2}{}^{})`$ conversion factor, $`\alpha `$, the factor of three is also consistent with no evolution.

The derived molecular gas masses of the most CO luminous QSOs at any given redshift are significantly below the average stellar mass of present-day elliptical galaxies ($``$ several$`\times 10^{11}`$ M). Given that massive elliptical galaxies are a likely end stage for many QSO host galaxies (e.g., as inferred from Bahcall et al. 1997 or the Magorrian Relation), the mass discrepancy can be rectified by concluding that present-day massive galaxies are created from the merger of smaller galaxies, and/or that galaxies with molecular gas masses in excess of several$`\times 10^{10}`$ M quickly consume the gas by forming high mass stars (e.g. Evans et al. 1996, 1998).

## 6. Summary

Millimeter (CO) observations of 10 infrared-excess, optically-selected Palomar-Green QSOs with redshifts in the range $`0.04<z<0.17`$ have been presented. These observations represent the first step towards detecting molecular gas in a complete sample of optical QSOs. The following conclusions have been reached:

1) Six out of 10 QSOs in the sample were detected in CO$`(10)`$. The CO detections in PG 1613+658 and PG 0838+770 confirm previously reported detections by Alloin et al. (1992). These detections, combined with previous detections of CO in the PG QSOs IZw1 and Mrk 1014, provide strong support for the possibility that molecular gas is a likely source of fuel for both AGN and star formation in QSO host galaxies.

2) All of the QSOs have high values of $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ relative to most $`z<0.3`$ infrared luminous galaxies of the same $`L_{\mathrm{IR}}`$. This is because QSOs contribute significantly to the heating of dust in their host galaxies, and/or because the dust is heated by massive stars formed with high efficiency (i.e., per unit molecular gas mass) relative to most infrared luminous galaxies. If large values of $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ are achieved primarily through dust heating by QSOs, a significant fraction of ULIGs with similar values of $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ may also contain buried AGN. Alternatively, if high $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ is due primarily to star formation, then high-mass star formation may be intimately connected to the QSO phenomenon.

3) The warm infrared luminous galaxies, which by definition have Seyfert-like mid-to-far infrared flux density ratios (i.e., $`f_{25}/f_{60}0.2`$), have similar ranges of $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ and $`L_{\mathrm{CO}}^{}`$ as the infrared-excess PG QSOs.

4) A comparison of the CO luminosities and optical morphologies of the QSO host galaxies indicates that QSOs with the highest infrared and CO luminosities reside in morphologically disturbed galaxies, whereas the rest have a mixture of barred spiral host galaxies and “featureless” host galaxies.

5) There is no convincing evidence of evolution in the molecular gas content of the most CO luminous QSOs in the redshift range $`z=05`$. The most CO luminous high-redshift QSO (BRI 1335-0415) has an $`L_{\mathrm{CO}}^{}`$ a factor of three greater than the most CO luminous low-redshift QSO (Mrk 1014). However, given the uncertainty in the $`L_{\mathrm{CO}}^{}`$–to–$`M(`$H$`{}_{2}{}^{})`$ conversion factor, $`\alpha `$, the factor of three is also consistent with no evolution.

We thank the staff and postdoctoral scholars of the Owens Valley Millimeter array for their support both during and after the observations were obtained. ASE also thanks J. Carpenter and A. Baker for useful discussions and assistance during the preparation of this manuscript, and the referee for comments that led to the clarification of several key points. ASE was supported by RF9736D and AST 0080881. JAS was supported by the Jet Propulsion Laboratory, California Institute of Technology, under contract with NASA. The Owens Valley Millimeter Array is a radio telescope facility operated by the California Institute of Technology and is supported by NSF grant AST 9981546. This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory.
Figure Captions

Figure 1. Integrated intensity maps of the six PG QSOs detected in CO($`10`$). a) PG 0838+770 - contours are plotted as 1$`\sigma \times (1.6,1.6,2.6,3.6,4.6,5.6,6.6)`$. The peak intensity is 0.017 Jy beam-1 and corresponds to the position RA=08:44:45.22 dec=+76:53:08.96 (J2000.0). b) PG 1119+120 - contours are plotted as 1$`\sigma \times (2.3,2.3,3.3,4.3,5.3)`$. The peak intensity is 0.014 Jy beam-1, and corresponds to the position RA=11:21:47.12 dec=+11:44:18.30 (J2000.0). c) PG 1351+640 - contours are plotted as 1$`\sigma \times (1.6,1.6,2.6,3.6,4.6)`$; the peak intensity is 0.0055 Jy beam-1, and corresponds to the position RA=13:53:15.62 dec=+63:45:45.72 (J2000.0). d) PG 1415+451 - contours are plotted as 1$`\sigma \times (2.2,2.2,3.2,4.2,5.2)`$; the peak intensity is 0.051 Jy beam-1, and corresponds to the position RA=14:17:00.76 dec=+44:56:06.50 (J2000.0). e) PG 1440+356 - contours are plotted as 1$`\sigma \times (2.3,2.3,3.3,4.3,5.3,6.3)`$; the peak intensity is 0.016 Jy beam-1, and corresponds to the position RA=14:42:07.48 dec=+35:26:22.33 (J2000.0). f) PG 1613+658 - contours are plotted as 1$`\sigma \times (1.6,1.6,2.6,3.6,4.6,5.6,6.6,7.6,8.6)`$; the peak intensity is 0.012 Jy beam-1, and corresponds to the position RA=16:13:57.15 dec=+65:43:09.62 (J2000.0).

Figure 2. CO($`10`$) spectra of the six detected PG QSOs. For PG 0838+770 and PG 1415+451, the data have been smoothed to a resolution of 16 MHz with a 8 MHz sampling. For PG 1119+120, the data have been smoothed to a resolution of 32 MHz with a 16 MHz sampling. For PG 1351+640, the data have been smoothed to a resolution of 40 MHz with a 20 MHz sampling. For PG 1440+356 and PG 1613+658, the data have been smoothed to a resolution of 24 MHz with a 12 MHz sampling.

Figure 3. A plot of the CO flux vs. the 100 $`\mu `$m flux density for PG QSOs and ULIGs detected to date in CO. The CO data for IZw1 and Mrk 1014 were obtained from Barvainis et al. (1989) and Solomon et al. (1997), respectively. The ULIGs CO data were obtained from Sanders et al. (1989b), Sanders, Scoville, & Soifer (1988, 1991), Solomon et al. (1997), and Evans et al. (1999). The 100$`\mu `$m flux densities were obtained from Sanders et al. (1989a, 1991) and Solomon et al. (1997). Arrows denote 3$`\sigma `$ upper limits on the CO luminosity of PG 1126-041, PG 1202+281, and PG 1402+261. Only upper limits on $`L_{\mathrm{CO}}^{}`$ and $`f_{100}`$ currently exist for PG 0007+106, and thus PG 0007+106 is not included on the plot. Adapted from Figure 4 of Solomon et al. (1997).

Figure 4. A plot of $`\mathrm{log}L_{\mathrm{CO}}^{}`$ versus redshift ($`z`$) for the low-$`z`$ PG QSO sample, a flux-limited sample ($`f_{60}>5.24`$ Jy) of infrared luminous galaxies and a sample of ultraluminous infrared galaxies. The vertical dashed lines represent the upper and lower $`z`$ boundaries of the PG QSO sample. The LIG and ULIGs data have been obtained from the same sources as in Figure 3. Arrows denote 3$`\sigma `$ upper limits on the CO luminosity of PG 1126-041, PG 1202+281, and PG 1402+261.

Figure 5. a) A plot of $`L_{\mathrm{IR}}`$ vs. $`L_{\mathrm{CO}}^{}`$ for the low-$`z`$ QSO sample, a flux-limited sample ($`f_{60\mu \mathrm{m}}>5.24`$ Jy) of infrared luminous galaxies and a sample of ultraluminous infrared galaxies. Arrows denote 3$`\sigma `$ upper limits on the CO luminosity of PG 1126-041, PG 1202+281, and PG 1402+261. b) A plot of $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ vs. $`L_{\mathrm{IR}}`$ for the low-$`z`$ QSO sample, a flux-limited sample ($`f_{60\mu \mathrm{m}}>5.24`$ Jy) of infrared luminous galaxies and a sample of ultraluminous infrared galaxies. Arrows denote 3$`\sigma `$ lower limits on $`L_{\mathrm{IR}}/L_{\mathrm{CO}}^{}`$ of PG 1126-041, PG 1202+281, and PG 1402+261. For simplicity, all of the cool infrared luminous galaxies are plotted as plus signs.

Figure 6. A plot of $`\mathrm{log}L_{\mathrm{CO}}^{}`$ versus redshift ($`z`$) for the low-$`z`$ PG QSO sample and the high-redshift QSOs detected in CO to date. The moderate and high-redshift QSO data have been obtained from the following references: 3C 48 ($`z=0.37`$: Scoville et al. 1993), H1413+117 ($`z=2.56`$: Barvainis et al. 1994), MG 0414+0534 ($`z=2.64`$: Barvainis et al. 1998), APM 08279+5255 ($`z=3.91`$: Downes et al. 1999), BRI 1335-0415 ($`z=4.41`$: Guilloteau et al. 1997), and BR 1202-0725 ($`z=4.69`$: Omont et al. 1996). For gravitationally lensed QSOs, $`L_{\mathrm{CO}}^{}`$ is plotted in terms of both the observed value and the intrinsic value. For MG 0414+0534, the amplification is unknown, thus MG 0414+0534 is plotted with a solid line extending downward from the observed $`L_{\mathrm{CO}}^{}`$. Adapted from Figure 3 of Frayer et al. (1999).


text:
# Untitled Document

ON THE GREEN FUNCTIONS OF

GRAVITATIONAL RADIATION THEORY

Giampiero Esposito

Istituto Nazionale di Fisica Nucleare, Sezione di Napoli, Complesso Universitario di Monte S. Angelo, Via Cintia, Edificio N’, 80126 Napoli, Italy Università di Napoli Federico II, Dipartimento di Scienze Fisiche, Complesso Universitario di Monte S. Angelo, Via Cintia, Edificio N’, 80126 Napoli, Italy

Abstract. Previous work in the literature has studied gravitational radiation in black-hole collisions at the speed of light. In particular, it had been proved that the perturbative field equations may all be reduced to equations in only two independent variables, by virtue of a conformal symmetry at each order in perturbation theory. The Green function for the perturbative field equations is here analyzed by studying the corresponding second-order hyperbolic operator with variable coefficients, instead of using the reduction method from the retarded flat-space Green function in four dimensions. After reduction to canonical form of this hyperbolic operator, the integral representation of the solution in terms of the Riemann function is obtained. The Riemann function solves a characteristic initial-value problem for which analytic formulae leading to the numerical solution are derived.

1. Introduction

The construction of suitable inverses of differential operators lies still at the very heart of many profound properties in classical and quantum field theory. For example, the theory of small disturbances in local field theory can only be built if suitable invertible operators are considered . In a path-integral formulation, these correspond to the gauge-field and ghost operators, respectively . Moreover, the Peierls bracket on the space of physical observables, which is a Poisson bracket preserving the invariance under the full infinite-dimensional symmetry group of the theory, is obtained from the advanced and retarded Green functions of the theory via the supercommutator function \[1–4\], and leads possibly to a deeper approach to quantization. Last, but not least, a perturbation approach to classical general relativity relies heavily on a careful construction of Green functions of operators of hyperbolic \[5–8\] and elliptic  type. In particular, following \[5–8\], we shall be concerned with the axisymmetric collision of two black holes travelling at the speed of light, each described in the centre-of-mass frame before the collision by an impulsive plane-fronted shock wave with energy $`\mu `$. One then passes to a new frame to which a large Lorentz boost is applied. There the energy $`\nu =\mu e^\alpha `$ of the incoming shock $`1`$ obeys $`\nu >>\lambda `$, where $`\lambda =\mu e^\alpha `$ is the energy of the incoming shock $`2`$ and $`\mathrm{e}^\alpha \sqrt{\frac{1+\beta }{1\beta }}`$ ($`\beta `$ being the usual relativistic parameter). In the boosted frame, to the future of the strong shock $`1`$, the metric can be expanded in the form

$$g_{ab}\nu ^2\left[\eta _{ab}+\underset{i=1}{\overset{\mathrm{}}{}}\left(\frac{\lambda }{\nu }\right)^ih_{ab}^{(i)}\right]$$
$`(1.1)`$
where $`\eta _{ab}`$ is the standard notation for the Minkowski metric. The task of solving the Einstein field equations becomes then a problem in singular perturbation theory, having to find $`h_{ab}^{(1)},h_{ab}^{(2)},\mathrm{}`$ by solving the linearized field equations at first, second, … order respectively in $`\frac{\lambda }{\nu }`$, once that characteristic initial data are given just to the future of the strong shock $`1`$. The perturbation series (1.1) is physically relevant because, on boosting back to the centre-of-mass frame, it is found to give an accurate description of space-time geometry where gravitational radiation propagates at small angles away from the forward symmetry axis $`\widehat{\theta }=0`$. The news function $`c_0`$, which describes gravitational radiation arriving at future null infinity in the centre-of-mass frame, is expected to have the convergent series expansion

$$c_0(\widehat{\tau },\widehat{\theta })=\underset{n=0}{\overset{\mathrm{}}{}}a_{2n}(\widehat{\tau }/\mu )(\mathrm{sin}\widehat{\theta })^{2n}$$
$`(1.2)`$
with $`\widehat{\tau }`$ a suitable retarded time coordinate, and $`\mu `$ the energy of each incoming black hole in the centre-of-mass frame. In  a very useful analytic expression of $`a_2(\widehat{\tau }/\mu )`$ was derived, exploiting the property that perturbative field equations may all be reduced to equations in only two independent variables, by virtue of a remarkable conformal symmetry at each order in perturbation theory. The Green function for perturbative field equations was then found by reduction from the retarded flat-space Green function in four dimensions.

However, a direct approach to the evaluation of Green functions appears both desirable and helpful in general, and it has been our aim to pursue such a line of investigation. For this purpose, reduction to two dimensions with the associated hyperbolic operator is studied again in section 2. Section 3 performs reduction to canonical form with the associated Riemann function. Equations for the Goursat problem obeyed by the Riemann function are derived in section 4, while the corresponding numerical algorithm is discussed in section 5.

2. Reduction to two dimensions and the associated operator

As is well known from the work in  and , the field equations for the first-order correction $`h_{ab}^{(1)}`$ in the expansion (1.1) are particular cases of the general system given by the flat-space wave equation (here $`u\frac{1}{\sqrt{2}}(z+t),v\frac{1}{\sqrt{2}}(zt)`$)

$$\text{ }\text{ / }\text{ }\psi =2\frac{^2\psi }{uv}+\frac{1}{\rho }\frac{}{\rho }\left(\rho \frac{\psi }{\rho }\right)+\frac{1}{\rho ^2}\frac{^2\psi }{\varphi ^2}=0$$
$`(2.1)`$
supplemented by the boundary condition

$$\psi (u=0)=\mathrm{e}^{im\varphi }\rho ^nf[8\mathrm{log}(v\rho )\sqrt{2}v]$$
$`(2.2a)`$
$$f(x)=0x<0.$$
$`(2.2b)`$
Moreover, $`\psi `$ should be of the form $`\mathrm{e}^{im\varphi }\rho ^n\chi (q,r)`$ for $`u0`$, where

$$qu\rho ^2$$
$`(2.3)`$
$$r8\mathrm{log}(\nu \rho )\sqrt{2}v.$$
$`(2.4)`$
For the homogeneous wave equation (2.1) there is no advantage in eliminating $`\rho `$ and $`\varphi `$ from the differential equation. However, the higher-order metric perturbations turn out to obey inhomogeneous flat-space wave equations of the form

$$\text{ }\text{ / }\text{ }\psi =S$$
$`(2.5)`$
where $`S`$ is a source term equal to $`\mathrm{e}^{im\varphi }\rho ^{(n+2)}H(q,r)`$. This leads to the following equation for $`\chi \mathrm{e}^{im\varphi }\rho ^n\psi `$:

$$_{m,n}\chi (q,r)=H(q,r)$$
$`(2.6)`$
where $`_{m,n}`$ is an hyperbolic operator in the independent variables $`q`$ and $`r`$, and takes the form

$$\begin{array}{ccc}\hfill _{m,n}& =(2\sqrt{2}+32q)\frac{^2}{qr}+4q^2\frac{^2}{q^2}+64\frac{^2}{r^2}\hfill & \\ & +4(n+1)q\frac{}{q}16n\frac{}{r}+n^2m^2.\hfill & (2.7)\hfill \end{array}$$

The proof of hyperbolicity of $`_{m,n}`$, with the associated normal hyperbolic form, can be found in section 3 of , and in . The advantage of studying Eq. (2.6) is twofold: to evaluate the solution at some space-time point one has simply to integrate the product of $`H`$ and the Green function $`G_{m,n}`$ of $`_{m,n}`$:

$$\chi (q,r)=G_{m,n}(q,r;q_0,r_0)H(q_0,r_0)𝑑q_0𝑑r_0$$
$`(2.8)`$
and the resulting numerical calculation of the solution is now feasible .

Since we are interested in a direct approach to the evaluation of the Green function $`G_{m,n}`$ in the $`(q,r)`$ coordinates, we begin by noticing that, in all derivatives which are not mixed, the operator $`_{m,n}`$ can be made a constant coefficient operator upon setting

$$\alpha \mathrm{log}(q)$$
$`(2.9)`$
which implies

$$q\frac{}{q}=\frac{}{\alpha }$$
$`(2.10)`$
$$q^2\frac{^2}{q^2}=\frac{^2}{\alpha ^2}\frac{}{\alpha }$$
$`(2.11)`$
and hence

$$\begin{array}{ccc}\hfill _{m,n}& =(2\sqrt{2}\mathrm{e}^\alpha +32)\frac{^2}{\alpha r}+4\frac{^2}{\alpha ^2}+64\frac{^2}{r^2}\hfill & \\ & +4n\frac{}{\alpha }16n\frac{}{r}+n^2m^2.\hfill & (2.12)\hfill \end{array}$$

This operator is further simplified upon defining the variable

$$R\frac{r}{4}$$
$`(2.13)`$
which implies that

$$\begin{array}{ccc}\hfill _{m,n}& =\left(\frac{1}{\sqrt{2}}\mathrm{e}^\alpha +8\right)\frac{^2}{\alpha R}+4\left(\frac{^2}{\alpha ^2}+\frac{^2}{R^2}\right)\hfill & \\ & +4n\left(\frac{}{\alpha }\frac{}{R}\right)+n^2m^2.\hfill & (2.14)\hfill \end{array}$$

This suggests defining yet new variables

$$X\alpha +R$$
$`(2.15)`$
$$Y\alpha R$$
$`(2.16)`$
so that

$$\frac{^2}{\alpha R}=\frac{^2}{X^2}\frac{^2}{Y^2}$$
$`(2.17)`$
$$\frac{^2}{\alpha ^2}+\frac{^2}{R^2}=2\left(\frac{^2}{X^2}+\frac{^2}{Y^2}\right)$$
$`(2.18)`$
$$\frac{}{\alpha }\frac{}{R}=2\frac{}{Y}$$
$`(2.19)`$
and hence $`_{m,n}`$ reads eventually

$$\begin{array}{ccc}\hfill T_{m,n}& =16\frac{^2}{Y^2}+8n\frac{}{Y}+n^2m^2\hfill & \\ & \frac{1}{\sqrt{2}}\mathrm{e}^{(X+Y)/2}\left(\frac{^2}{X^2}\frac{^2}{Y^2}\right)\hfill & (2.20)\hfill \end{array}$$

with Green function satisfying the equation

$$\begin{array}{ccc}& T_{m,n}G_{m,n}(\mathrm{e}^{\frac{X+Y}{2}},2(XY);\mathrm{e}^{\frac{X_0+Y_0}{2}},2(X_0Y_0))\hfill & \\ & =\frac{1}{2}\delta \left(\mathrm{e}^{\frac{X+Y}{2}}\mathrm{e}^{\frac{X_0+Y_0}{2}}\right)\delta ((XY)(X_0Y_0)).\hfill & (2.21)\hfill \end{array}$$

The operator $`T_{m,n}`$ is the sum of an elliptic operator in the $`Y`$ variable and a two-dimensional wave operator ‘weighted’ with the exponential $`\mathrm{e}^{(X+Y)/2}`$, which is the main source of technical complications in these variables.

3. Reduction to canonical form and the Riemann function

It is therefore more convenient, in our general analysis, to reduce first Eq. (2.6) to canonical form, and then find an integral representation of the solution. Reduction to canonical form means that new coordinates $`x=x(q,r)`$ and $`y=y(q,r)`$ are introduced such that the coefficients of $`\frac{^2}{x^2}`$ and $`\frac{^2}{y^2}`$ vanish. As is shown in , this is achieved if

$$\frac{x}{r}=\frac{y}{r}=1$$
$`(3.1)`$
$$\frac{x}{q}=\frac{1+8q\sqrt{2}+\sqrt{1+16q\sqrt{2}}}{2\sqrt{2}q^2}$$
$`(3.2a)`$
$$\frac{y}{q}=\frac{1+8q\sqrt{2}\sqrt{1+16q\sqrt{2}}}{2\sqrt{2}q^2}.$$
$`(3.3a)`$
The resulting formulae are considerably simplified if one defines

$$t\sqrt{1+16q\sqrt{2}}=t(x,y).$$
$`(3.4)`$
The dependence of $`t`$ on $`x`$ and $`y`$ is obtained implicitly by solving the system

$$x=r+8\mathrm{log}\left(\frac{t1}{2}\right)\frac{8}{(t1)}4$$
$`(3.5)`$
$$y=r+8\mathrm{log}\left(\frac{t+1}{2}\right)+\frac{8}{(t+1)}4.$$
$`(3.6)`$
This leads to the equation

$$\mathrm{log}\frac{(t1)}{(t+1)}\frac{2t}{(t^21)}=\frac{(xy)}{8}$$
$`(3.7a)`$
which can be cast in the form

$$\frac{(t1)}{(t+1)}\mathrm{e}^{\frac{2t}{(1t^2)}}=\mathrm{e}^{\frac{(xy)}{8}}.$$
$`(3.7b)`$
This suggests defining

$$w\frac{(t1)}{(t+1)}$$
$`(3.8)`$
so that one first has to solve the transcendental equation

$$w\mathrm{e}^{\frac{(w^21)}{2w}}=\mathrm{e}^{\frac{(xy)}{8}}$$
$`(3.9)`$
to obtain $`w=w(xy)`$, from which one gets

$$t=\frac{(1+w)}{(1w)}=t(xy).$$
$`(3.10)`$
On denoting by $`g(w)`$ the left-hand side of Eq. (3.9), one finds that, in the plane $`(w,g(w))`$, the right-hand side of Eq. (3.9) is a line parallel to the $`w`$-axis, which intersects $`g(w)`$ at no more than one point for each value of $`xy`$. For example, when $`w=1`$, $`g(w)`$ intersects the line taking the constant value $`1`$, for which $`xy=0`$. The function

$$g:wg(w)=w\mathrm{e}^{\frac{(w^21)}{2w}}$$

is asymmetric and has the limiting behaviour described by

$$\underset{w0^{}}{lim}g(w)=\mathrm{}\underset{w0^+}{lim}g(w)=0$$
$`(3.11)`$
$$\underset{w\mathrm{}}{lim}g(w)=0\underset{w+\mathrm{}}{lim}g(w)=\mathrm{}.$$
$`(3.12)`$
Thus, in the lower half-plane, $`g`$ has an horizontal asymptote given by the $`w`$-axis, and a vertical asymptote given by the line $`w=0`$, while it has no asymptotes in the upper half-plane, since

$$\underset{w\mathrm{}}{lim}\frac{g(w)}{w}=\mathrm{}$$

in addition to (3.12). The first derivative of $`g`$ reads

$$g^{}(w)=\frac{(w+1)^2}{2w}\mathrm{e}^{\frac{(w^21)}{2w}}.$$
$`(3.13)`$
One therefore has $`g^{}(w)>0`$ for all $`w>0`$, and $`g^{}(w)<0`$ for all $`w(\mathrm{},0)\left\{1\right\}`$, and $`g`$ is monotonically decreasing for negative $`w`$ and monotonically increasing for positive $`w`$. The point $`w=1`$, at which $`g^{}(w)`$ vanishes, is neither a maximum nor a minimum point, because

$$g^{\prime \prime }(w)=\left(\frac{1}{4w^3}+\frac{1}{2w}+1+\frac{w}{4}\right)\mathrm{e}^{\frac{(w^21)}{2w}}$$
$`(3.14)`$
$$g^{\prime \prime \prime }(w)=\left(\frac{1}{8w^5}\frac{3}{4w^4}+\frac{3}{8w^3}+\frac{3}{8w}+\frac{3}{4}+\frac{w}{8}\right)\mathrm{e}^{\frac{(w^21)}{2w}}.$$
$`(3.15)`$
These formulae imply that $`g^{\prime \prime }(1)=0`$ but $`g^{\prime \prime \prime }(1)=10`$, and hence $`w=1`$ yields a flex of $`g(w)`$.

In the $`(x,y)`$ variables, the operator $`_{m,n}`$ reads therefore

$$_{m,n}=f(x,y)\frac{^2}{xy}+g(x,y)\frac{}{x}+h(x,y)\frac{}{y}+n^2m^2$$
$`(3.16)`$
where, exploiting the formulae

$$\frac{x}{q}=\frac{64\sqrt{2}}{(t1)^2}$$
$`(3.2b)`$
$$\frac{y}{q}=\frac{64\sqrt{2}}{(t+1)^2}$$
$`(3.3b)`$
one finds

$$\begin{array}{ccc}\hfill f(x,y)& =(2\sqrt{2}+32q)\left(\frac{x}{q}+\frac{y}{q}\right)+8q^2\frac{x}{q}\frac{y}{q}+128\hfill & \\ & =256\left[1\frac{2t^2(t^2+1)}{(t1)^2(t+1)^2}\right]\hfill & (3.17)\hfill \end{array}$$

$$g(x,y)=4(n+1)q\frac{x}{q}16n=16\left[1+\frac{2(n+1)}{(t1)}\right]$$
$`(3.18)`$
$$h(x,y)=4(n+1)q\frac{y}{q}16n=16\left[1\frac{2(n+1)}{(t+1)}\right].$$
$`(3.19)`$
The resulting canonical form of Eq. (2.6) is

$$\begin{array}{ccc}\hfill L[\chi ]& =\left(\frac{^2}{xy}+a(x,y)\frac{}{x}+b(x,y)\frac{}{y}+c(x,y)\right)\chi (x,y)\hfill & \\ & =\stackrel{~}{H}(x,y)\hfill & (3.20)\hfill \end{array}$$

where

$$a(x,y)\frac{g(x,y)}{f(x,y)}=\frac{1}{16}\frac{(1t)(t+1)^2(2n+1+t)}{(t^4+4t^21)}$$
$`(3.21)`$
$$b(x,y)\frac{h(x,y)}{f(x,y)}=\frac{1}{16}\frac{(t+1)(t1)^2(2n+1t)}{(t^4+4t^21)}$$
$`(3.22)`$
$$c(x,y)\frac{n^2m^2}{f(x,y)}=\frac{(m^2n^2)}{256}\frac{(t1)^2(t+1)^2}{(t^4+4t^21)}$$
$`(3.23)`$
$$\stackrel{~}{H}(x,y)\frac{H(x,y)}{f(x,y)}=\frac{H(x,y)}{256}\frac{(t1)^2(t+1)^2}{(t^4+4t^21)}.$$
$`(3.24)`$
Note that $`a(t)=b(t),b(t)=a(t),c(t)=c(t),\stackrel{~}{H}(t)=\stackrel{~}{H}(t)`$.

For an hyperbolic equation in the form (3.20), we can use the Riemann integral representation of the solution. For this purpose, recall from  that, on denoting by $`L^{}`$ the adjoint of the operator $`L`$ in (3.20), which acts according to

$$L^{}[\chi ]=\chi _{xy}(a\chi )_x(b\chi )_y+c\chi $$
$`(3.25)`$
one has to find a ‘function’ $`R(x,y;\xi ,\eta )`$ (actually a kernel) subject to the following conditions ($`(\xi ,\eta )`$ being the coordinates of a point $`P`$ such that characteristics through it intersect a curve $`C`$ at points $`A`$ and $`B`$, $`AP`$ being a segment with constant $`y`$, and $`BP`$ being a segment with constant $`x`$): (i) As a function of $`x`$ and $`y`$, $`R`$ satisfies the adjoint equation

$$L_{(x,y)}^{}[R]=0$$
$`(3.26)`$
(ii) $`R_x=bR`$ on $`AP`$, i.e.

$$R_x(x,y;\xi ,\eta )=b(x,\eta )R(x,y;\xi ,\eta )\mathrm{on}y=\eta $$
$`(3.27)`$
and $`R_y=aR`$ on $`BP`$, i.e.

$$R_y(x,y;\xi ,\eta )=a(\xi ,y)R(x,y;\xi ,\eta )\mathrm{on}x=\xi $$
$`(3.28)`$
(iii) $`R`$ equals $`1`$ at $`P`$, i.e.

$$R(\xi ,\eta ;\xi ,\eta )=1.$$
$`(3.29)`$
It is then possible to express the solution of Eq. (3.20) in the form

$$\begin{array}{ccc}\hfill \chi (P)& =\frac{1}{2}[\chi (A)R(A)+\chi (B)R(B)]\hfill & \\ & +_{AB}([\frac{R}{2}\chi _x+(bR\frac{1}{2}R_x)\chi ]dx\hfill & \\ & [\frac{R}{2}\chi _y+(aR\frac{1}{2}R_y)\chi ]dy)\hfill & \\ & +_\mathrm{\Omega }R(x,y;\xi ,\eta )\stackrel{~}{H}(x,y)𝑑x𝑑y\hfill & (3.30)\hfill \end{array}$$

where $`\mathrm{\Omega }`$ is a domain with boundary.

Note that Eqs. (3.27) and (3.28) are ordinary differential equations for the Riemann function $`R(x,y;\xi ,\eta )`$ along the characteristics parallel to the coordinate axes. By virtue of (3.29), their integration yields

$$R(x,\eta ;\xi ,\eta )=\mathrm{exp}_\xi ^xb(\lambda ,\eta )𝑑\lambda $$
$`(3.31)`$
$$R(\xi ,y;\xi ,\eta )=\mathrm{exp}_\eta ^ya(\lambda ,\xi )𝑑\lambda $$
$`(3.32)`$
which are the values of $`R`$ along the characteristics through $`P`$. Equation (3.30) yields instead the solution of Eq. (3.20) for arbitrary initial values given along an arbitrary non-characteristic curve $`C`$, by means of a solution $`R`$ of the adjoint equation (3.26) which depends on $`x,y`$ and two parameters $`\xi ,\eta `$. Unlike $`\chi `$, the Riemann function $`R`$ solves a characteristic initial-value problem.

4. Goursat problem for the Riemann function

By fully exploiting the reduction to canonical form of Eq. (2.6) we have considered novel features with respect to the analysis in , because the Riemann formula (3.30) contains also the integral along the piece of curve $`C`$ from $`A`$ to $`B`$, and the term $`\frac{1}{2}[\chi (A)R(A)+\chi (B)R(B)]`$. This representation of the solution might be more appropriate for the numerical purposes considered in , but the task of finding the Riemann function $`R`$ remains extremely difficult. One can however use approximate methods for solving Eq. (3.26). For this purpose, we first point out that, by virtue of Eq. (3.25), Eq. (3.26) is a canonical hyperbolic equation of the form

$$\left(\frac{^2}{xy}+A\frac{}{x}+B\frac{}{y}+C\right)R(x,y;\xi ,\eta )=0$$
$`(4.1)`$
where

$$Aa$$
$`(4.2)`$
$$Bb$$
$`(4.3)`$
$$Cca_xb_y.$$
$`(4.4)`$
Thus, on defining

$$UR$$
$`(4.5)`$
$$VR_x+BR$$
$`(4.6)`$
the equation (4.1) for the Riemann function is equivalent to the hyperbolic canonical system

$$U_x=f_1(x,y)U+f_2(x,y)V$$
$`(4.7)`$
$$V_y=g_1(x,y)U+g_2(x,y)V$$
$`(4.8)`$
where

$$f_1B=b$$
$`(4.9)`$
$$f_2=1$$
$`(4.10)`$
$$g_1ABC+B_y=abc+a_x$$
$`(4.11)`$
$$g_2A=a.$$
$`(4.12)`$
For the system described by Eqs. (4.7) and (4.8) with boundary data (3.31) and (3.32) an existence and uniqueness theorem holds (see  for the Lipschitz conditions on boundary data), and we can therefore exploit the finite differences method to find approximate solutions for the Riemann function $`R(x,y;\xi ,\eta )`$, and eventually $`\chi (P)`$ with the help of the integral representation (3.30).

5. Concluding remarks

The inverses of hyperbolic operators  and the Cauchy problem for hyperbolic equations with polynomial coefficients  have always been the object of intensive investigation in the mathematical literature. We have here considered the application of such issues to axisymmetric black hole collisions at the speed of light, relying on the work in \[5–8\]. We have pointed out that, for the inhomogeneous equations (2.6) occurring in the perturbative analysis, the task of inverting the operator (2.7) can be accomplished with the help of the Riemann integral representation (3.30), after solving Eq. (4.1) for the Riemann function. One has then to solve a characteristic initial-value problem for a homogeneous hyperbolic equation in canonical form in two independent variables, for which we have developed formulae to be used for the numerical solution with the help of a finite differences scheme. For this purpose one studies the canonical system (cf (4.7) and (4.8))

$$U_x=F(x,y,U,V)$$
$`(5.1)`$
$$V_y=G(x,y,U,V)$$
$`(5.2)`$
in the rectangle $`\{x,y:x[x_0,x_0+a],y[y_0,y_0+b]\}`$ with known values of $`U`$ on the vertical side $`AD`$ where $`x=x_0`$, and known values of $`V`$ on the horizontal side $`AB`$ where $`y=y_0`$. The segments $`AB`$ and $`AD`$ are then divided into $`m`$ and $`n`$ equal parts, respectively. On setting $`\frac{a}{m}h`$ and $`\frac{b}{n}k`$, the original differential equations become equations relating values of $`U`$ and $`V`$ at three intersection points of the resulting lattice, i.e.

$$\frac{U(P_{r,s+1})U(P_{rs})}{h}=F$$
$`(5.3a)`$
$$\frac{V(P_{r+1,s})V(P_{rs})}{k}=G.$$
$`(5.4a)`$
It is now convenient to set $`U_{rs}U(P_{rs}),V_{rs}V(P_{rs})`$, so that these equations read

$$U_{r,s+1}=U_{rs}+hF(P_{rs},U_{rs},V_{rs})$$
$`(5.3b)`$
$$V_{r+1,s}=V_{rs}+kG(P_{rs},U_{rs},V_{rs}).$$
$`(5.4b)`$
Thus, if both $`U`$ and $`V`$ are known at $`P_{rs}`$, one can evaluate $`U`$ at $`P_{r,s+1}`$ and $`V`$ at $`P_{r+1,s}`$. The evaluation at subsequent intersection points of the lattice goes on along horizontal or vertical segments. In the former case, the resulting algorithm is

$$U_{rs}=U_{r0}+h\underset{i=1}{\overset{s1}{}}F(P_{ri},U_{ri},V_{ri})$$
$`(5.5)`$
$$V_{rs}=V_{r1,s}+kG(P_{r1,s},U_{r1,s},V_{r1,s})$$
$`(5.6)`$
while in the latter case one obtains the algorithm expressed by the equations

$$V_{rs}=V_{0s}+\underset{i=1}{\overset{r1}{}}G(P_{is},U_{is},V_{is})$$
$`(5.7)`$
$$U_{rs}=U_{r,s1}+hF(P_{r,s1},U_{r,s1},V_{r,s1}).$$
$`(5.8)`$
Stability of such solutions is closely linked with the geometry of the associated characteristics, and the criteria to be fulfilled are studied in section 13.2 of  (stability depends crucially on whether or not $`\frac{h}{k}1`$).

To sum up, one solves numerically Eq. (3.9) for $`w=w(x,y)=w(xy)`$, from which one gets $`t(xy)`$ with the help of (3.10), which is a fractional linear transformation. This yields $`a,b,c`$ and $`\stackrel{~}{H}`$ as functions of $`(x,y)`$ according to (3.21)–(3.24), and hence $`A,B`$ and $`C`$ in the equation for the Riemann function are obtained according to (4.2)–(4.4), where derivatives with respect to $`x`$ and $`y`$ are evaluated numerically. Eventually, the system given by (4.7) and (4.8) is solved according to the finite-differences scheme of the present section, with

$$F=f_1U+f_2V=f_1R+f_2(R_x+BR)$$
$`(5.9)`$
$$G=g_1U+g_2V=g_1R+g_2(R_x+BR).$$
$`(5.10)`$
Once the Riemann function $`R=U`$ is obtained with the desired accuracy, numerical evaluation of the integral (3.30) yields $`\chi (P)`$, and $`\chi (q,r)`$ is obtained upon using Eqs. (3.5) and (3.6) for the characteristic coordinates. Our steps are conceptually desirable since they rely on well established techniques for the solution of hyperbolic equations in two independent variables , and provide a viable alternative to the numerical analysis performed in , because all functions should be evaluated numerically. Our method is not obviously more powerful than the one used in \[5–8\], but is well suited for a systematic and lengthy numerical analysis, while its analytic side provides an interesting alternative for the evaluation of Green functions both in black hole physics and in other problems where hyperbolic operators with variable coefficients might occur. This task remains very important because a strong production of gravitational radiation is mainly expected in the extreme events studied in \[5–8\] and which motivated our paper. Any viable way of looking at mathematical and numerical aspects of the problem is therefore of physical interest for research planned in the years to come .

Acknowledgments

The work of G Esposito has been partially supported by the National project SINTESI 2000.

References

 DeWitt B S 1965 Dynamical Theory of Groups and Fields (New York: Gordon and Breach)

 DeWitt B S 1984 in Relativity, Groups and Topology II eds B S DeWitt and R Stora (Amsterdam: North–Holland)

 Esposito G 2001 Quantum Gravity in Four Dimensions (New York: Nova Science)

 Cartier P and DeWitt–Morette C 2000 J. Math. Phys. 41 4154

 D’Eath P D and Payne P N 1992 Phys. Rev. D 46 658

 D’Eath P D and Payne P N 1992 Phys. Rev. D 46 675

 D’Eath P D and Payne P N 1992 Phys. Rev. D 46 694

 D’Eath P D 1996 Black Holes: Gravitational Interactions (Oxford: Clarendon Press)

 Esposito G and Stornaiolo C 2000 Class. Quantum Grav. 17 1989

 Esposito G and Stornaiolo C 2000 Found. Phys. Lett. 13 279

 Courant R and Hilbert D 1961 Methods of Mathematical Physics. II. Partial Differential Equations (New York: Interscience)

 Leray J 1953 Hyperbolic Differential Equations (Princeton University Press)

 Leray J 1956 C. R. Acad. Sci. Paris 242 953

 Garabedian P R 1964 Partial Differential Equations (New York: Chelsea)

 Allen B and Ottewill A (2000) Gen. Rel. Grav. 32 385