format: text

text:
# Cluster physics from joint weak gravitational lensing and Sunyaev-Zel’dovich data

## 1 Introduction

Whereas clusters of galaxies, as the largest gravitationnaly bound structures of the universe, form natural probe of cosmology, observations, numerical simulations as well as timing arguments provide compelling evidences that most of them are young and complex systems. Interaction with large-scale structures, merging processes and coupling of dark matter with the intra-cluster medium complicate the interpretation of observations and the modeling of each of its components. Since they are composed of dark mater (DM), galaxies and a hot dilute X-ray emitting gas (Intra cluster medium, ICM) accounting respectively for $`85\%`$, $`15\%`$ and $`5\%`$ of their mass, the physics of the ICM bounded in a dark matter gravitational potential plays a major role in cluster formation and evolution. This variety of components can be observed in many various ways. In particular, gravitational lensing effects (the weak-lensing regime here, WL) Mellier (2000); Bartelmann and Schneider (2001), Sunyaev-Zel’dovich (SZ) effect Sunyaev and Zel’dovich (1972); Birkinshaw (1999) and X-ray emission (X) Sarrazin (1988). Whereas the former probes mostly the dark matter component, both the latter probe the baryons of the gravitationally bound ICM.

Due to observational progress, increasingly high quality data are delivered which enables multi-wavelength investigation of clusters on arcminute scale (the most recent is the spectacular progress in SZ measurements, e.g. Reese et al.  (2000); Désert et al.  (1998)) and we therefore think it is timely to explore how we should perform some joint analysis of these high quality data sets and exploit them at best their complementarity. This challenge has already been tackled by several groups Zaroubi et al.  (1998); Grego et al.  (1999); Reblinsky (2000); Zaroubi et al.  (2000); Castander et al.  (2000); Holder et al.  (2000). Zaroubi et al. and Reblinsky et al. attempted a full deprojection by assuming isothermality and axial symmetry, using respectively a least square minimization or a Lucy-Richardson algorithm , Grego et al. compare SZ derived gas mass to WL derived total mass by fitting a spheroidal $`\beta `$ model. But whereas these methods give reasonable results it has been illustrated, e.g. by Inagaki et al. 1995 in the context of $`H_0`$ measurement from SZ and X-ray observations, that both non isothermality and asphericity analysis can trigger systematic errors as high as $`20\%`$. Therefore, we aim at exploring an original approach which allows to get rid of both isothermality and departure from sphericity. Based on a self-consistent use of both observables, and based on a perturbative development of general physical hypothesis, this method allow us to test some very general physical hypothesis of the gas (hydrostatic equilibrium, global thermodynamic equilibrium) and also provide naturally some X observation predictions.

Observations only provide us with $`2D`$ projected quantities (e.g. mass, gas pression,…). This quantities are related by some physical hypothesis which are explicited in $`3D`$ equalities (e.g. hydrostatic equilibrium, equation of state). The point is that these $`3D`$ equalities do not have any tractable equivalent relating projected $`2D`$ quantities: in particular, projection along the line of sight does not provide an equation of state or a projected hydrostatic equilibrium equation. Therefore as soon as we want to compare this data (WL, SZ, X) we have to deproject the relevant physical quantities ($`P_g,T_g,\rho _g`$…). This can be done only using strong assumptions, either by using parametric models (e.g. a $`\beta `$ model Cavaliere & Fusco-Femaino (1976)) or by assuming mere geometrical hypothesis (the former necessarily encompassing the latter) Fabian et al.  (1981); Yoshikawa and Suto (1999). We choose the geometric approach in order to use as general physical grounds as possible and to avoid as many theoretical biases as possible.

This simplest choice might be naturally motivated first by looking at some images of observed clusters Désert et al.  (1998); Grego et al.  (1999). Their regularity is striking : some have almost circular or ellipsoidal appearance as we expect for fully relaxed system. Then since relaxed clusters are expected to be spheroidal in favored hierarchical structure formation scenario, it is natural to try to relate the observed quasi-circularity (quasi-sphericity) to the $`3D`$ quasi-sphericity (quasi-spheroidality). We perform this using some linearly perturbed spherical (spheroidal) symmetries in a self-consistent approach.

We proceed as follows: in section 2 we defined our physical hypothesis and our notations. The method is precisely described in section 3. We consider both the spherical as well as spheroidal cases and obtain a predicted X surface brightness map from a SZ decrement map and a WL gravitational distortion map. In section 4 a demonstration with simulated clusters is presented before discussing its application to genuine data as well as further developments in section 5.

## 2 Hypothesis, Sunyaev-Zel’dovich effect and the Weak lensing

We now briefly describe our notations as well as our physical hypothesis.

### 2.1 General hypothesis

Following considerations fully detailed in Sarrazin (1988) the ICM can be regarded as a hot and dilute plasma constituted from ions and electrons, whose respective kinetic temperatures $`T_p`$ and $`T_e`$ will be considered as equal $`T_p=T_eT_g`$. This is the *global thermodynamic equilibrium hypothesis* which is expected to hold up to $`r_{virial}`$ ( see Teyssier et al.  (1997); Chièze et al.  (1998) for a precise discussion). Given the low density (from $`n_e10^1\mathrm{cm}^3`$ in the core to $`10^5\mathrm{cm}^3`$ in the outer part) and high temperature of this plasma ($`10\mathrm{k}\mathrm{e}\mathrm{V}`$), it can be treated as a perfect gas satisfying the equation of state :

$$P_g=\frac{\rho _gk_BT_g}{\mu _em_p}=\beta \rho _gT_g$$
(1)
with $`\beta \frac{k_B}{\mu _em_p}`$. Let us neglect then the gas mass with regards to the dark matter mass, and assume *stationarity* (no gravitational potential variation on time scale smaller than the hydrodynamic time scale, e.g. no recent mergers). Then the gas assumed to be in hydrostatic equilibrium in the dark matter gravitational potential satisfies:

$`(\rho _g𝐯_𝐠)`$ $`=`$ $`0`$ (2)
$`P_g`$ $`=`$ $`\rho _g\mathrm{\Phi }_{DM}.`$ (3)
At this point there is no need to assume isothermality.

### 2.2 Sunyaev-Zel’dovich effect and weak lensing

Inverse Compton scattering of cosmic background (CMB) photons by the electrons in the ICM modifies the CMB spectrum Zel’dovich and Sunyaev (1969); Sunyaev and Zel’dovich (1972, 1980). The amplitude of the SZ temperature decrement $`\frac{\mathrm{\Delta }T_{SZ}}{T_{CMB}}`$ is directly proportional to the Comptonisation parameter $`y`$ which is given by :

$`y={\displaystyle \frac{\sigma _T}{m_ec^2}}{\displaystyle 𝑑ln_ek_BT_e}`$ $`=`$ $`{\displaystyle \frac{\sigma _T}{m_ec^2}}{\displaystyle 𝑑lp_e}`$ (4)
$`={\displaystyle \frac{\sigma _T}{m_ec^2}}{\displaystyle 𝑑l\frac{\rho _gk_BT_g}{\mu _em_p}}`$ $`=`$ $`\alpha {\displaystyle 𝑑lP_g}.`$ (5)
where $`\alpha \frac{\sigma _T}{m_ec^2}`$, $`k_B`$ is the Boltzmann’s constant, $`\sigma _T`$ is the Thomson scattering cross section and $`dl`$ is the physical line-of-sight distance. $`m_e`$, $`n_e`$, $`T_e`$ and $`p_e`$ are the mass, the number density, the temperature and the thermal pressure of electrons. $`\rho _g`$ and $`T_g`$ respectively denote the gas density and temperature, and $`\mu _e`$ is the number of electrons per proton mass. Some further corrections to this expression can be found in Rephaeli (1995); Birkinshaw (1999).

In parallel to this spectral distortions, the statistical determination of the shear field $`\kappa `$ affecting the images of background galaxies enable, in the weak lensing regime, to derive the dominant projected gravitational potential of the lens (the clustered dark matter) : $`\varphi _{DM}`$ in our general hypothesis (see Mellier (2000) for details).

## 3 Method

### 3.1 Principle

We now answer the question : how should we co-analyze these various data set ? Our first aim is to develop a method which allows us to get maps of projected thermodynamical quantities with as few physical hypothesis as possible.

Our method is the following. Let us suppose we have for a given cluster a set of data a SZ and WL data which enables us to construct a $`2D`$ map of projected gas pressure as well as a $`2D`$ projected gravitational potential map. Let us suppose as well that these maps exhibit an approximate spherical symmetry as it is the case for a vast class of experimental observations as e.g. in figure 1. More precisely, let us suppose that the projected gas pressure $`y`$ as well as the observed projected gravitational potential $`\varphi _{DM}`$ can be well fitted by the following type of functions :

$`y(R,\phi )`$ $`=`$ $`y_0(R)+\epsilon y_1(R)m(\phi )`$ (6)
$`\varphi _{DM}(R,\phi )`$ $`=`$ $`\varphi _{DM,0}(R)+\epsilon \varphi _{DM,1}(R)n(\phi )`$ (7)
where $`\epsilon 1`$, $`(R,\phi )`$ denotes polar coordinates in the image plane and $`m`$ and $`n`$ are some particular functions. This description means first of all that the images we see are linear perturbations from some perfect circularly symmetric images, and second that the perturbation might be described conveniently by the product of a radial function and an angular function. Equivalently we can assert that to first order in $`\epsilon `$ our images are circularly symmetric but they admit some corrections to second order in $`\epsilon `$.
We then assume that these observed perturbed symmetries are a consequence of an intrinsic $`3D`$ spherical symmetry linearly perturbed too. This point constitutes our key hypothesis. It means that to first order in a certain parameter (e.g. $`\epsilon `$) our clusters are regular objects with a strong circular symmetry but they admit some second order linear perturbations away from this symmetry. As a consequence of these assumptions we will make use of this linearly perturbed symmetry to get a map of some complementary projected thermodynamical quantities, the gas density $`D_g`$ and the gas temperature $`\zeta _g`$, successively to first and second order in $`\epsilon `$.

Formulated this way, the problem yields a natural protocol :

* Looking at some maps with this kind of symmetry, we compute a zero-order map ($`y_0(R)`$, $`\varphi _0(R)`$) with a perfect circular symmetry by averaging over some concentric annulus. A correction for the bias introduced by perturbations is included. These first order quantities allow us to derive some first order maps of $`D_{g,0}(R)`$ and $`\zeta _{g,0}(R)`$ with a perfect circular symmetry.
* We then take into account the first order corrections to this perfect symmetry ($`y_1(R)m(\phi )`$, $`\varphi _1(R)m(\phi )`$) and infer from them first order correction terms to the zeroth order maps: $`D_{g,1}(R,\phi )`$ and $`\zeta _{g,1}(R,\phi )`$.

Even if for clarity’s sake we formulate our method assuming a perturbed circular symmetry, it applies equivalently to a perturbed elliptical symmetry as it will be shown below. In this more general case, we assume that the cluster exhibit a linearly perturbed spheroidal symmetry.

### 3.2 The spherically symmetric case : from observations to predictions

Let us now apply the method to the case where the projected gas density (SZ data) and the projected gravitational potential (WL data) exhibit some approximate circular symmetry. These observations lead us to suppose that the $`3D`$ gas pressure, the gravitational potential, the gas density and the gas temperature can be well described by the following equations:

$$\{\begin{array}{cccccc}P_g(r,\theta ,\phi )& =& P_{g,0}(r)& +& \epsilon & P_{g,1}(r)f(\theta ,\phi )\\ \mathrm{\Phi }_{DM}(r,\theta ,\phi )& =& \mathrm{\Phi }_{DM,0}(r)& +& \epsilon & \mathrm{\Phi }_{DM,1}(r)g(\theta ,\phi )\\ \rho _g(r,\theta ,\phi )& =& \rho _{g,0}(r)& +& \epsilon & \rho _{g,1}(r)h(\theta ,\phi )\\ T_g(r,\theta ,\phi )& =& T_{g,0}(r)& +& \epsilon & T_{g,1}(r)k(\theta ,\phi )\end{array}$$
(8)
where $`(r,\theta ,\phi )`$ are spherical coordinates centered on the cluster.

#### 3.2.1 The hydrostatic equilibrium

If we first apply the hydrostatic equilibrium equation $`P_g=\rho _g\mathrm{\Phi }_{DM}`$ we get the following equations. To first order in $`\epsilon `$ we have

$$P_{g,0}^{}(r)=\rho _{g,0}(r)\mathrm{\Phi }_{DM,0}^{}(r),$$
(9)
and to second order in $`\epsilon `$ :

$$\{\begin{array}{cccc}P_{g,1}^{}(r)f(\theta ,\phi )& =& \rho _{g,0}(r)\mathrm{\Phi }_{DM,1}^{}(r)h(\theta ,\phi )& \\ & & \rho _{g,1}(r)\mathrm{\Phi }_{DM,0}^{}(r)g(\theta ,\phi )& (a)\\ P_{g,1}(r)_\theta f(\theta ,\phi )& =& \rho _{g,0}(r)\mathrm{\Phi }_{DM,1}(r)_\theta h(\theta ,\phi )& (b)\\ P_{g,1}(r)_\phi f(\theta ,\phi )& =& \rho _{g,0}(r)\mathrm{\Phi }_{DM,1}(r)_\phi h(\theta ,\phi )& (c)\end{array}$$
(10)
where “ ’ “ denotes the derivative with regards to $`r`$.
Combining equations (10.b) and (10.c) we get

$`f(\theta ,\phi )=\lambda _1h(\theta ,\phi )+\lambda _2`$ (11)
where $`\lambda _{1,2}`$ are some constants. Then, using equation (10.a) we can write

$`f(\theta ,\phi )=\gamma _1g(\theta ,\phi )+\gamma _2`$ (12)
where $`\gamma _{1,2}`$ are some constants as well. At this point, we can get rid of $`\lambda _2`$ and $`\gamma _2`$ by absorbing them in the order $`1`$ mere radial term (i.e. $`\rho _{g,0}(r)`$ and $`\mathrm{\Phi }_{DM,0}(r)`$). This means we can consider $`\lambda _2=0`$ and $`\gamma _2=0`$ . Similarly we choose to rescale $`\rho _{g,1}(r)`$ and $`\mathrm{\Phi }_{DM,1}(r)`$ so that we can take $`\gamma _1=\lambda _1=1`$. These simple equalities lead us to assume from now on :

$$f(\theta ,\phi )=h(\theta ,\phi )=g(\theta ,\phi ).$$
(13)
This is in no way a restriction since it simply means that we absorb integration constants by redefining some terms. This is possible since the relevant part of $`f`$ (and thus $`h`$) will be fitted on observations as will be shown below. Taking equation (13) into account, equation (10) simplifies to :

$`P_{g,0}^{}(r)`$ $`=`$ $`\rho _{g,0}(r)\mathrm{\Phi }_{DM,0}^{}(r)`$ (14)
$`P_{g,1}^{}(r)`$ $`=`$ $`\rho _{g,0}(r)\mathrm{\Phi }_{DM,1}^{}(r)\rho _{g,1}(r)\mathrm{\Phi }_{DM,0}^{}(r)`$ (15)
$`P_{g,1}(r)`$ $`=`$ $`\rho _{g,0}(r)\mathrm{\Phi }_{DM,1}(r).`$ (16)
#### 3.2.2 The equation of state

We have now identified the angular part to the first order correction of $`P_g`$, $`\mathrm{\Phi }_{DM}`$ and $`\rho _g`$. We still have to link those quantities to the angular dependent part of the temperature $`T_g`$, namely $`k(\theta ,\phi )`$. This is done naturally using the equation of state (1), which directly provide to first and second order in $`\epsilon `$ :

$`P_{g,0}(r)`$ $`=`$ $`\beta \rho _{g,0}(r)T_{g,0}(r)`$ (17)
$`P_{g,1}(r)f(\theta ,\phi )`$ $`=`$ $`\beta \rho _{g,1}(r)T_{g,0}(r)f(\theta ,\phi )`$
$`+`$ $`\beta \rho _{g,0}(r)T_{g,1}(r)k(\theta ,\phi )`$ (18)
This last equation leads naturally to $`f(\theta ,\phi )=k(\theta ,\phi )`$ if we decide once again to absorb any multiplicative factor in the radial part. This way we see that our choice of separating the radial and angular part is in no way a restriction. We eventually get

$`P_{g,0}(r)`$ $`=`$ $`\beta \rho _{g,0}(r)T_{g,0}(r)`$ (19)
$`P_{g,1}(r)`$ $`=`$ $`\beta \rho _{g,1}(r)T_{g,0}(r)+\beta \rho _{g,0}(r)T_{g,1}(r).`$ (20)
#### 3.2.3 The observations

Given this description of the cluster hot gas, the experimental SZ and WL data which respectively provide us with the projected quantities $`y(R,\phi )`$ and $`\varphi _{DM}(R,\phi )`$ write

$`y(R,\phi )`$ $`=`$ $`\alpha {\displaystyle P_{g,0}(r)𝑑l}+\epsilon \alpha {\displaystyle P_{g,1}(r)f(\theta ,\phi )𝑑l}`$ (21)
$``$ $`y_0(R)+\epsilon y_1(R)m(\phi )`$
$`\varphi _{DM}(R,\phi )`$ $`=`$ $`{\displaystyle \mathrm{\Phi }_{DM,0}(r)𝑑l}+\epsilon {\displaystyle \mathrm{\Phi }_{DM,1}(r)f(\theta ,\phi )𝑑l}`$ (22)
$``$ $`\varphi _{DM,0}(R)+\epsilon \varphi _{DM,1}(R)m(\phi ).`$
Note that in order to get this set of definitions we choose the polar axis of the cluster along the line of sight so that the same azimuthal angle $`\phi `$ is used for $`2D`$ and $`3D`$ quantities.

Our aim is now to derive both a projected gas density map and projected temperature map that we define this way :

$`D_g(R,\phi )`$ $`=`$ $`{\displaystyle \rho _g(r,\phi )𝑑l}`$ (23)
$`=`$ $`{\displaystyle \rho _{g,0}(r)𝑑l}+\epsilon {\displaystyle \rho _{g,1}(r)f(\theta ,\phi )𝑑l}`$ (24)
$``$ $`D_{g,0}(R)+D_{g,1}(R,\phi )`$ (25)
$`\zeta _g(R,\phi )`$ $`=`$ $`{\displaystyle T_g(r,\phi )𝑑l}`$ (26)
$`=`$ $`{\displaystyle T_{g,0}(r)𝑑l}+\epsilon {\displaystyle T_{g,1}(r)f(\theta ,\phi )𝑑l}`$ (27)
$``$ $`\zeta _{g,0}(R)+\zeta _{g,1}(R,\phi ).`$ (28)
#### 3.2.4 A projected gas density map to first order…

Now that we have expressed our observables in terms of $`3D`$ physical quantities, it is easy to infer a gas density map successively to first and second order in $`\epsilon `$. To first order the hydrostatic equilibrium condition (9) states that

$$P_{g,0}^{}(r)=\rho _{g,0}(r)\mathrm{\Phi }_{DM,0}^{}(r).$$
(29)
In order to use it we need to deproject the relevant quantities. From the well known spherical deprojection formula Binney and Tremaine (1987) based on Abel’s transform we have :

$`\alpha P_{g,0}(r)`$ $`=`$ $`{\displaystyle \frac{1}{\pi }}{\displaystyle _r^{\mathrm{}}}y_0^{}(R){\displaystyle \frac{dR}{(R^2r^2)^{\frac{1}{2}}}}`$ (30)
$`=`$ $`{\displaystyle \frac{1}{\pi }}{\displaystyle _0^{\mathrm{}}}y_0^{}(r\mathrm{cosh}u)𝑑u`$ (31)
where $`R=r\mathrm{cosh}u`$. Thus, we can write

$`\alpha P_{g,0}^{}(r)`$ $`=`$ $`{\displaystyle \frac{1}{\pi }}{\displaystyle _0^{\mathrm{}}}\mathrm{cosh}uy^{\prime \prime }(r\mathrm{cosh}u)𝑑u`$ (32)
$`=`$ $`{\displaystyle \frac{1}{\pi }}{\displaystyle _r^{\mathrm{}}}{\displaystyle \frac{1}{r}}{\displaystyle \frac{R}{(R^2r^2)^{\frac{1}{2}}}}y_0^{\prime \prime }(R)𝑑R.`$ (33)
Similarly,

$$\mathrm{\Phi }_{DM,0}^{}(r)=\frac{1}{\pi }_r^{\mathrm{}}\frac{1}{r}\frac{R}{(R^2r^2)^{\frac{1}{2}}}\varphi _0^{\prime \prime }(R)𝑑R.$$
(34)
We then get for the projected gas density

$`D_{g,0}(R)=2{\displaystyle _R^{\mathrm{}}}{\displaystyle \frac{rdr}{(r^2R^2)^{\frac{1}{2}}}}{\displaystyle \frac{P_{g,0}^{}(r)}{\mathrm{\Phi }_{DM,0}^{}(r)}}`$ (36)
$`={\displaystyle \frac{2}{\alpha }}{\displaystyle _R^{\mathrm{}}}{\displaystyle \frac{rdr}{(r^2R^2)^{\frac{1}{2}}}}\left({\displaystyle \frac{_r^{\mathrm{}}\frac{sds}{r(s^2r^2)^{\frac{1}{2}}}y_0^{\prime \prime }(s)}{_r^{\mathrm{}}\frac{sds}{r(s^2r^2)^{\frac{1}{2}}}\varphi _0^{\prime \prime }(s)}}\right).`$
#### 3.2.5 …and a projected gas temperature map to first order

Once we built this projected gas density map, we can recover the projected gas temperature map. If we apply the equation of state (17) we get :

$`\zeta _{g,0}(R)`$ $`=`$ $`{\displaystyle \frac{1}{\beta }}{\displaystyle \frac{P_{g,0}(r)}{\rho _{g,0}(r)}𝑑l}`$ (37)
$`=`$ $`{\displaystyle \frac{1}{\beta }}{\displaystyle \frac{P_{g,0}(r)}{P_{g,0}^{}(r)}\mathrm{\Phi }_{DM,0}^{}(r)𝑑l}`$ (38)
$`=`$ $`{\displaystyle \frac{1}{\pi \beta }}{\displaystyle _R^{\mathrm{}}}{\displaystyle \frac{P_{g,0}(r)}{P_{g,0}^{}(r)}}\mathrm{\Phi }_{DM,0}^{}(r){\displaystyle \frac{rdr}{(r^2R^2)^{\frac{1}{2}}}}.`$ (39)
Since all the required functions ($`P_{g,0}`$, $`P_{g,0}^{}`$, $`\mathrm{\Phi }_{DM,0}^{}`$) have been derived in the previous section (equation (31) and (33)) we can get this way a projected gas temperature map.

#### 3.2.6 Corrections from departure to spherical symmetry : a projected gas density map to second order…

We now reach the core of our method, namely we aim at deriving the quantity $`D_{g,1}`$ defined by (25), i.e. the second order correction to the perfectly circular term :

$`D_g(R,\phi )`$ $`=`$ $`D_{g,0}(R)+\epsilon D_{g,1}(R,\phi )`$ (40)
$`=`$ $`{\displaystyle \rho _{g,0}(r)𝑑l}+\epsilon {\displaystyle \rho _{g,1}(r)f(\theta ,\phi )𝑑l}.`$ (41)
If we derive equation (16) and combine it with equation (15) we note that

$`\rho _{g,0}^{}(r)\mathrm{\Phi }_{DM,1}(r)`$ $`=`$ $`\rho _{g,1}(r)\mathrm{\Phi }_{DM,0}^{}(r).`$ (42)
Therefore we can write

$`{\displaystyle \rho _{g,1}(r)f(\theta ,\phi )𝑑l}={\displaystyle \frac{\rho _{g,0}^{}(r)}{\mathrm{\Phi }_{DM,0}^{}(r)}\mathrm{\Phi }_{DM,1}(r)f(\theta ,\phi )𝑑l}.`$ (43)
At this point we want to express this quantity either in terms of WL data or in terms of SZ data depending on the quality of them, or even better in terms of an optimal combination of them.

On one hand, WL data provide us with a straightforward access to the function $`\varphi _1(R)m(\phi )=\mathrm{\Phi }_{DM,1}(r)f(\theta ,\phi )𝑑l`$ thus we choose to approximate (43) by

$`{\displaystyle \rho _{g,1}(r)f(\theta ,\phi )𝑑l}{\displaystyle \frac{\rho _{g,0}^{}(R)}{\mathrm{\Phi }_{DM,0}^{}(R)}}{\displaystyle \mathrm{\Phi }_{DM,1}(r)f(\theta ,\phi )𝑑l}`$ (44)
$`{\displaystyle \frac{\rho _{g,0}^{}(R)}{\mathrm{\Phi }_{DM,0}^{}(R)}}\varphi _1(R)m(\phi )`$
$`{\displaystyle \frac{\rho _{g,0}^{}(R)}{\mathrm{\Phi }_{DM,0}^{}(R)}}\left(\varphi _{DM}(R,\phi )\varphi _0(R)\right)`$
where we used the definitions of section (3.2.3) and where $`R`$ corresponds to the radius observed in the image plane, i.e. the radius $`r`$ equal to the distance between the line of sight and the center of the cluster. We will discuss this approximation in more details in section (3.2.8) and validate it through a practical implementation on simulations in section (4). But we already can make the following statements: would the line of sight follows a line of constant $`r`$ throughout the domain of the perturbation, this expression would be rigorously exact. Moreover it turns out to be a good approximation because of the finite extent of the perturbation.

On the other hand SZ data provide us with a measurement of the function $`y_1(R)m(\phi )=P_{g,1}(r)f(\theta ,\phi )𝑑l`$ therefore we can use equation (16) and (14) to write

$`{\displaystyle \rho _{g,1}(r)f(\theta ,\phi )𝑑l}={\displaystyle \frac{\rho _{g,0}^{}(r)}{P_{g,0}^{}(r)}P_{g,1}(r)f(\theta ,\phi )𝑑l}`$ (45)
$``$ $`{\displaystyle \frac{\rho _{g,0}^{}(R)}{P_{g,0}^{}(R)}}{\displaystyle P_{g,1}(r)f(\theta ,\phi )𝑑l}`$ (46)
$``$ $`{\displaystyle \frac{\rho _{g,0}^{}(R)}{P_{g,0}^{}(R)}}y_1(R)m(\phi )`$ (47)
$``$ $`{\displaystyle \frac{\rho _{g,0}^{}(R)}{P_{g,0}^{}(R)}}\left(y(R,\phi )y_0(R)\right).`$ (48)
Here again we used the same notation and approximation as in equation (44). Note however that as soon as we assumed isothermality, the ratio $`\rho _{g,0}^{}/P_{g,0}^{}`$ is constant therefore this last step is exact. Were we not assuming isothermality, the departure from isothermality is expected to be weak thus this last approximation should be reasonable.

This last two alternative steps are crucial to our method since these approximations link the non spherically symmetric components of various quantities. They are reasonable as will be discussed in section (3.2.8) and will be numerically tested in section (4).
Of course, only well-known quantities appear in equation (44) and (48): $`y`$, $`y_0`$, $`\varphi _{DM}`$ and $`\varphi _0`$ are direct observational data whereas $`P_{g,0}(r)`$ and $`\rho _{g,0}(r)`$ are zeroth order quantities previously derived.

#### 3.2.7 …and a projected gas temperature map to second order

The projected temperature map can be obtained the same way as before. Using first the equation of state we can write :

$`T_{g,0}(r)+\epsilon T_{g,1}(r)f(\theta ,\phi )={\displaystyle \frac{1}{\beta }}{\displaystyle \frac{(P_{g,0}(r)+\epsilon P_{g,1}(r)f(\theta ,\phi ))}{(\rho _{g,0}(r)+\epsilon \rho _{g,1}(r)f(\theta ,\phi ))}}`$ (49)
$``$ $`{\displaystyle \frac{1}{\beta }}\left({\displaystyle \frac{P_{g,0}(r)}{\rho _{g,0}(r)}}+\epsilon P_{g,1}(r){\displaystyle \frac{\rho _{g,0}(r)\rho _{g,1}(r)}{\rho _{g,0}^2(r)}}f(\theta ,\phi )\right).`$
Hence, since

$`\zeta (R,\phi )`$ $`=`$ $`\zeta _0(R,\phi )+\epsilon \zeta _1(R,\phi )`$ (50)
$`=`$ $`{\displaystyle T_{g,0}(r)𝑑l}+\epsilon {\displaystyle T_{g,1}(r)f(\theta ,\phi )𝑑l}`$ (51)
we have

$`\zeta _1(R,\phi )={\displaystyle \frac{\rho _{g,0}(r)\rho _{g,1}(r)}{\rho _{g,0}^2(r)}P_{g,1}(r)f(\theta ,\phi )𝑑l}.`$ (52)
Here we choose to approximate the last integral as previously discussed in order to make use of observational SZ data. Therefore we rewrite this last equation as :

$`\zeta _1(R,\phi ){\displaystyle \frac{\rho _{g,0}(R)\rho _{g,1}(R)}{\rho _{g,0}^2(R)}}{\displaystyle P_{g,1}(r)f(\theta ,\phi )𝑑l}`$ (53)
$``$ $`{\displaystyle \frac{\rho _{g,0}(R)\rho _{g,1}(R)}{\rho _{g,0}^2(R)}}y_1(R)m(\phi )`$
$``$ $`{\displaystyle \frac{\rho _{g,0}(R)\rho _{g,1}(R)}{\rho _{g,0}^2(R)}}\left(y(R,\phi )y_0(R)\right).`$ (54)
We obtain this way an expression to second order for the projected temperature in terms of either observed quantities or previously derived functions.

#### 3.2.8 Why the previous approximation is reasonable on intuitive grounds?

Our previous approximations can be justified on intuitive grounds even if we will take care of validating it numerically in section (4) below. It relies on the fact that perturbations have by definition a finite extent, i.e. the first order correction to the perfectly circular (spherical) term is non zero only within a finite range. The typical size and the amplitude of the perturbation can be easily scaled from the SZ and WL data set. This guarantees the validity of our assumptions on observational grounds. The key point is that the perturbation itself has a kind of axial symmetry, whose axis goes through the center of the cluster and the peak of the perturbation. This is reasonable if the perturbation originates in e.g. an incoming filament but not for a substructure. The latter would therefore have to be treated separately by superposition (see section (5)). This leads naturally to the statement that the typical angle we observe in the image plane is equal to the one we would observe if the line of sight were perpendicular to its actual direction, i.e. the perturbation as intrinsically the same angular extent in the directions along the line of sight and perpendicular to it. This is illustrated schematically in figure (2).

Given this description we are now in a position to discuss the validity of our approximation. It consists in approximating the line of sight integral $`g(r)\mathrm{\Phi }_{DM,1}(r)f(\theta ,\phi )𝑑l`$ by $`g(R)\mathrm{\Phi }_{DM,1}(r)f(\theta ,\phi )𝑑l`$ where $`g`$ is any radial function. This approximation would be exact if $`g(r)`$ were constant in the relevant domain, i.e. if the line of sight had a constant $`r`$. As mentioned before this is the case in equation 48 if we assume isothermality. But the functions $`g(r)`$ we might deal with may scale roughly as $`r^2`$, as e.g. $`\rho _{g,0}^{}(r)/P_{g,0}(r)`$ in equation (44), thus it is far from being constant. The consequent error committed can be estimated by the quantity $`\mathrm{\Delta }rg^{}(r)`$ where $`\mathrm{\Delta }r`$ is the maximum $`r`$ discrepancy between the value assumed, $`g(R)`$, and the actual value as it is schematically illustrated in figure (3). In the worst case, $`g^{}(r)`$ scales as $`r`$. Then, using the obvious notations defined in this figure we get

$$(\mathrm{\Delta }r)_{max}=R(11/\mathrm{sin}(\theta \frac{\mathrm{\Delta }\theta }{2})).$$
(55)
Naturally this quantity is minimal for $`\theta 90^o`$ and diverges for $`\theta 0^o`$ when $`\mathrm{\Delta }\theta =0^o`$ : the error is minimal when the line of sight is nearly tangential ($`\theta 90^o`$) and so almost radial in this domain, and maximal when it is radial ($`\theta =0^o`$). This in principle is a very bad behavior, but the fact is that the closer $`\theta `$ is from $`0^o`$ the weaker the integrated perturbation is since it gets always more degenerate along the line of sight, i.e. the integrated perturbations tend to a radial behavior and will therefore be absorbed in the $`\mathrm{\Phi }_{DM,0}(r)`$ term. The extreme situation, i.e. when $`\theta =0^o`$ will trigger a mere radial image as long as the perturbation exhibits a kind of axial symmetry. This error is impossible to alleviate since we are dealing with a fully degenerate situation but will not flaw the method at all since the integrated perturbation will be null. This approximation will be validated numerically below.

### 3.3 How to obtain a X prediction ?

The previously derived map offers a great interest that we now aim at exploiting, namely the ability of precise X prediction. Indeed, for a given X spectral emissivity model, the X-ray spectral surface brightness is

$$S_X(E)=\frac{1}{4\pi (1+z)^4}n_e^2\mathrm{\Lambda }(E,T_e)𝑑l$$
(56)
where $`\mathrm{\Lambda }`$ is the spectral emissivity, $`z`$ is the redshift of the cluster and $`E`$ is the energy on which the observed band is centered. Hence we can write, assuming a satisfying knowledge of $`z`$ and $`\mathrm{\Lambda }`$ :

$`S_X(E)`$ $``$ $`{\displaystyle n_e^2T_e^{1/2}𝑑l}`$ (57)
$``$ $`{\displaystyle \rho _g^2T_g^{1/2}𝑑l}`$ (59)
$``$ $`{\displaystyle \rho _{g,0}^2T_{g,0}^{1/2}𝑑l}+2\epsilon {\displaystyle \rho _{g,0}T_{g,0}^{1/2}\rho _{g,1}f(\theta ,\phi )𝑑l}`$
$`+{\displaystyle \frac{1}{2}}\epsilon {\displaystyle \rho _{g,0}^2T_{g,0}^{1/2}T_{g,1}f(\theta ,\phi )𝑑l}`$
where we omitted to write the $`(r)`$s for clarity’s sake. If we now make use of the same approximation as used and discussed before, we can express directly this quantity in terms of observations $`y`$ and $`\varphi `$. We get indeed

$`S_X(E)`$ $``$ $`{\displaystyle \rho _{g,0}^2T_{g,0}^{1/2}𝑑l}`$ (61)
$`+\mathrm{\hspace{0.25em}\hspace{0.25em}\hspace{0.25em}2}\epsilon \rho _{g,0}(R)T_{g,0}^{1/2}(R){\displaystyle \rho _{g,1}f(\theta ,\phi )𝑑l}`$
$`+{\displaystyle \frac{1}{2}}\epsilon \rho _{g,0}^2(R)T_{g,0}^{1/2}(R){\displaystyle T_{g,1}f(\theta ,\phi )}`$
$``$ $`{\displaystyle \rho _{g,0}^2T_{g,0}^{1/2}𝑑l}+2\epsilon \rho _{g,0}(R)T_{g,0}^{1/2}(R)D_{g,1}(R,\phi )`$
$`+{\displaystyle \frac{1}{2}}\epsilon \rho _{g,0}^2(R)T_{g,0}^{1/2}(R)\zeta _{g,1}(R,\phi ).`$
Both the first order terms $`T_{g,0}`$ and $`\rho _{g,0}`$, and the second order corrections $`D_{g,1}`$ and $`\zeta _{g,1}`$ have been derived in the previous sections. We are thus able to generate self-consistently a X luminosity map from our previously derived maps. This is a very nice feature of this method. We will further discuss the approximation and its potential bias in the next section.

This derivation opens the possibility of comparing on the one hand SZ and WL observations with, on the other hand, precise X-ray measurements as done e.g. by XMM or CHANDRA. Note that in the instrumental bands of most of X-ray satellites the $`T_g`$ dependence is very weak and can be neglected. This can be easily taken into account by eliminating the $`T_g`$ dependence in the previous formula. Even if the interest of such a new comparison is obvious we will discuss it more carefully in the two following sections. In principle, one could also easily make some predictions concerning the density weighted X-ray temperature defined by the ratio $`n_g^2T_g𝑑l/n_g^2𝑑l`$ but the fact is that since the gas pressure and so the SZ effect tends to have a very weak gradient we are not able by principle to reproduce all the interesting features of this quantity, namely the presence of shocks.

## 4 Application on simulations

In order to demonstrate the ability of the method in a simplified context we used some outputs of the recently developed N-body + hydrodynamics code RAMSES simulating the evolution of a $`\mathrm{\Lambda }`$-CDM universe. The RAMSES code is based on Adaptative Mesh Refinement (AMR) technics in order to increase the spatial resolution locally using a tree of recursively nested cells of smaller and smaller size. It reaches a formal resolution of $`12\mathrm{kpch}^1`$ in the core of galaxy clusters (see Refregier and Teyssier 2000 and Teyssier 2001, in preparation, for details). We use here the structure of $`2`$ galaxy cluster extracted of the simulation to generate our needed observables, i.e. X-ray emission measure, SZ decrement and projected density (or projected gravitational potential).

The relevant observables, i.e. projected mass density, SZ decrement and for comparison purpose only the X-ray emission measure, of the 2 clusters are depicted using a logarithmic scaling in figure 4 and 5 (upper panels). This clusters have been extracted of the simulation at $`z=0.0`$ and thus tends to be more relaxed. They are ordinary clusters of virial mass (defined by $`\delta _{334}`$ in our particular cosmology) $`\mathrm{4.50\; 10}^{14}h^1\mathrm{M}_{\mathrm{}}`$ and $`\mathrm{4.15\; 10}^{14}h^1\mathrm{M}_{\mathrm{}}`$. Both exhibit rather regular shape, i.e. they have not undergo recently a major merge. The depicted boxes are respectively $`3.5h^1\mathrm{Mpc}`$ and $`4.0h^1\mathrm{Mpc}`$ wide. We smooth the outputs using a gaussian of width $`120h^1\mathrm{kpc}`$ thus degrading the resolution. We did not introduce any instrumental noise. This clusters are to a good approximation isothermal thus for the sake of simplicity we will assume that $`T_g`$ is constant making the discussion on $`T_{g,0}`$ and $`T_{g,1}`$ useless at this point. We apply the method previously described using perturbed spherical symmetry. We deduce by averaging over concentric annuli a zeroth order circular description of the gas density and then add to it some first order corrections. Note that since we assume isothermality SZ data give us straightforwardly a projected gas density modulo a temperature $`T_{g,0}`$ coefficient, thus we use the formulation of equation (48), exact in this context. This constant temperature is fixed using the hydrostatic equilibrium and the WL data.

In figure 4 and 5 (lower panels) we show the predicted X-ray emission measure to zeroth and first order as well as a map of relative errors. Note that to first order the shape of the emission measure is very well reproduced. The cross-correlation coefficients between the predicted and simulated X-ray emission measure are $`0.978`$ and $`0.986`$. Of course this is partly due to the assumed good quality of the assumed SZ data but nonetheless, it demonstrates the validity of our perturbative approach as well as of our approximation. The approximation performed in equation (61), i.e. the multiplication by the function $`\rho _{g,0}(R)`$ will naturally tends to cut out the perturbations at high $`R`$. This is the reason why the further perturbation are slightly less well reproduced and the relative errors tend to increase with $`R`$. Nevertheless, since the emission falls rapidly with $`R`$ as visible on the lower figures (note the logarithmic scaling) the total flux is well conserved, respectively to $`0.9\%`$ and $`9\%`$. This last number might illustrate that the large extent of the perturbations in the second case may limit our method. An ellipsoidal fit could have help decrease this value. Note that moreover the clump visible mainly in X-ray emission measure of figure 5 is not reproduce. This is natural because it does not appear through the SZ effect since the pression remains uniform throughout clumps. If resolved by WL, this substructure should anyway be treated separately, e.g. by considering the addition of a second very small structure. Note that the first cluster showed exhibits a spherical core elongated in the outter region thus it is not actually as ellipsoidal as it looks which may explain why our perturbed spherical symmetry works well.

## 5 Discussion

### 5.1 Hypothesis …and non hypothesis

Our approach makes several assumptions. Some general and robust hypothesis have been introduced and discussed in section 2.1. Note that we do not need to assume isothermality. Our key hypothesis consists in assuming the validity of a perturbative approach and in the choice of the nature of this perturbations, i.e. with a radial/angular part separation. Theoretical predictions, observations and simulations show that relaxed clusters are regular and globally spheroidal objects, which is what initially motivated our approach. Then in our demonstration on simulations, this turns out to be reasonable. Such an approach can not deal properly with sharp features as e.g. shocks waves due to infalling filaments. Then assuming the validity of the angular and radial separation, leads to the equality of this angular parts for all relevant physical quantities ($`P_g`$, $`T_g`$, $`\varphi _{DM}`$…), using to first order in $`\epsilon `$ the hydrostatic equilibrium and the equation of state. If this is not satisfied in practice then we could either question the validity of this separation or the physics of the cluster. Our experience with simulation shows that for reasonably relaxed clusters, i.e. not going through a major merge, the angular part of the perturbation is constant amongst observables. Thus it looks like the separation (and thus the equality of the angular perturbation) is a good hypothesis in general and its failure is a sign of non-relaxation, i.e. non-validity of our general physical hypothesis.

Then an important hypothesis lies in the validity of the approximation used. Note first that even if its form is general, its validity depends on the quantity which is assumed to be constant along the integral. In the case of the gas density obtained from the SZ map, it is an exact statement as soon as we assume the isothermality and since clusters in general are not too far from isothermality, this hypothesis is reasonable.

Now, some worth to remember “non hypothesis” are the isothermality and the sphericity (or ellipsoidality). This might be of importance. Indeed, in evaluating the Hubble constant from joint SZ and X-ray measurement it has been evaluated in Inagaki et al.  (1995); Roettiger et al.  (1997); Puy et al.  (2000) that, both the asphericity and the non-isothermality of the relevant cluster can yield some important bias (up to $`20\%`$). Even if this measure is not our concern here, it is interesting to note that this hypothesis are not required here.

### 5.2 The equivalent spheroidal symmetry case

So far, we have work and discussed the perturbed spherical symmetry case. If we turn to spheroidal symmetry the problem is very similar as long as we assume the knowledge of the inclination angle $`i`$ between the polar axis of the system and the line of sight. This is what we recall in appendix B which is directly inspired from Fabricant et al.  (1984): once the projection is nicely parametrised we get for the projected quantity , e.g. for the pressure :

$`y(\eta )`$ $`=`$ $`2{\displaystyle \frac{B_e}{R}}{\displaystyle _\eta ^{\mathrm{}}}{\displaystyle \frac{P_{g,0}(t)tdt}{(t^2\eta ^2)^{\frac{1}{2}}}}`$ (62)
$`P_{g,0}(t)`$ $`=`$ $`{\displaystyle \frac{1}{2\pi }}{\displaystyle \frac{R}{B_e}}{\displaystyle _t^{\mathrm{}}}P_{g,0}^{}(\eta ){\displaystyle \frac{d\eta }{\left(\eta ^2t^2\right)^{1/2}}}.`$ (63)
following the notations of appendix B. Since we are dealing with the same Abel integral we can proceed in two steps as we did before.

Even if the inclination angle is *a priori* not accessible directly through single observations it has been demonstrated that it is possible to evaluate it using the deprojection of an axially symmetric distribution of either X-ray/SZ maps or SZ/surface density maps Zaroubi et al.  (1998, 2000). Our approach in this work try to avoid to explicit the full 3-D structure rather than building it, and this is done in a simple self-consistent way therefore we will not get into the details of this procedure that will be discussed in a coming work (Doré et al. 2001, in preparation). Note also that axially symmetric configuration elongated along the line of sight may appear as spherical. This is a difficult bias to alleviate without any prior for the profile. In our case, our method will be biased in the sense that the deprojected profile will be wrong. Nevertheless, we might hope to reproduce properly the global quantities, like abundance of DM or gas and so to alleviate some well known systematics (see previous section), e.g. in measuring the baryon fraction.

## 6 Conclusion and outlook

It this paper we have presented and demonstrated the efficiency of an original method allowing to perform in a self-consistent manner the joint analysis of SZ and WL data. Using it on noise free simulation we demonstrated how well it can be used to make some x-ray surface brightness prediction, or equivalently emission measure. Our choice in this approach has been to hide somehow the deprojection by using some appropriate approximations. Thus we do not resolved fully the 3-D structure of clusters, but note that the work presented here is definitely a first step towards a full deprojection (Doré et al. 2001, in preparation). Some further refinements of the methods are under progress as well.

When applying the method to true data, the instrumental noise issue is an important matter of concern. Indeed, whereas the strong advantage of a parametric approach, e.g. using a $`\beta `$-model, is that it allows to adjust the relevant parameters, e.g. $`r_c`$ and $`\beta `$, on the projected quantities (the image) itself, which is rather robust to noise, it might be delicate to determine the profiles and its derivate by a direct deprojection. Nevertheless, our perturbative approach, as it first relies on a zeroth order quantity found by averaging over some annulus, a noise killing step (at least far from the center), and then work on some mere projected perturbation should be quite robust as well. Consequently we hope to apply it very soon on true data. Furthermore, in this context it should allow a better treatment of systematics (asphericity, non isothermality,…) plaguing any measure of the baryon fraction $`f_b`$ or the Hubble constant $`\mathrm{H}_0`$ using X-ray and SZ effect Inagaki et al.  (1995). These points will be discussed somewhere else (Doré et al. 2001, in preparation ).

## Acknowledgment

O.D. is grateful to G. Mamon, M. Bartelmann, S. Zaroubi and especially S. Dos Santos for valuable discussions. We thank J. Calrstrom et al. for allowing the use of some of their SZ images.

## Annexe : Deprojection in spheroidal symmetry

In this appendix we recall some useful results concerning spheroid projection derived by Fabricant, Gorenstein and Rybicki Fabricant et al.  (1984). In the context of spheroidal systems, cartesian coordinates system are the most convenient for projection. Thus, if the observer’s coordinate system $`(x,y,z)`$ is chosen such that the line of sight is along the $`z`$ axis and such that the polar axis of the spheroidal system $`z^{}`$ lies in the $`xz`$ plane at an inclination angle $`i`$ to the z-axis, then, in the cartesian coordinate system $`(x^{},y^{},z^{})`$ the general physical quantities relevant to our problem depends only on the parameter $`t`$ defined by

$`t^2`$ $`=`$ $`{\displaystyle \frac{x^2+y^2}{B_e^2}}+{\displaystyle \frac{z^2}{A_e^2}}`$ (64)
$`=`$ $`{\displaystyle \frac{(x\mathrm{cos}i+yz\mathrm{sin}i)^2+y^2}{B_e^2}}+{\displaystyle \frac{(z\mathrm{cos}ix\mathrm{sin}i)^2}{A_e^2}}.`$ (65)
If we project a physical quantity $`G(t)`$ on the observer sky plane $`xy`$ then,

$`I(x,y)`$ $`=`$ $`I(\eta )`$ (66)
$`=`$ $`{\displaystyle _{\mathrm{}}^+\mathrm{}}G(t)𝑑l`$ (67)
$`=`$ $`2{\displaystyle \frac{B_e}{R}}{\displaystyle _\eta ^{\mathrm{}}}{\displaystyle \frac{G(t)tdt}{(t^2\eta ^2)^{\frac{1}{2}}}}`$ (68)
where

$`\eta ^2`$ $``$ $`{\displaystyle \frac{x^2}{(RA_e)^2}}+{\displaystyle \frac{y^2}{(B_e)^2}}`$ (69)
$`\text{and }R`$ $``$ $`\sqrt{{\displaystyle \frac{B_e^2}{A_e^2}}\mathrm{cos}^2i+sin^2i}.`$ (70)
Of course this result shows that if we were to observe a spheroidal system we would map ellipses with an axial ratio equal to $`{\displaystyle \frac{B}{A}}={\displaystyle \frac{1}{R}}{\displaystyle \frac{B_e}{A_e}}`$. But the main result of this appendix is that we obtain at the end an Abel integral similar to the one obtained in the case of spherical system, where the radius as been replaced by the parameter $`t`$. This simple fact justifies the very analogous treatment developed in this paper for spherical and spheroidal systems.


text:
# Quantization of a generally covariant gauge system with two super Hamiltonian constraints

## I Introduction

Generally covariant theories such as Einstein’s theory of gravitation have the peculiar property of featuring a Hamiltonian that is constrained to vanish. This constraint is associated with the invariance of the action under reparametrizations. In the case of general relativity, the Arnowitt-Deser-Misner (ADM) Hamiltonian is the sum of four constraint functions, three of them are the supermomenta $`_a`$ (linear and homogeneous functions of the field momenta) and the other one is the super-Hamiltonian $``$ (a quadratic function of the field momenta). The supermomenta generate the change of the field under diffeomorphisms on the spatial hypersurfaces $`\mathrm{\Sigma }`$ of the foliated space-time manifold, and express the invariance of the action under these transformations. The super-Hamiltonian generates the evolution of the field under displacements that are normal to the hypersurface $`\mathrm{\Sigma }`$; this evolution is nothing but a physically irrelevant reparametrization of the dynamical trajectory of the system. Actually there are four constraint functions in each point of $`\mathrm{\Sigma }`$. So the algebra of constraints is rather complicated, because one should evaluate the Poisson brackets between constraints at different points. In particular, the Poisson brackets between two super-Hamiltonian take part in the algebra of constraints. This algebra was calculated by Dirac:

$`\{(x),(x^{})\}=^a(x)\delta _{,a}(x,x^{})+^a(x^{})\delta _{,a}(x,x^{}),`$ (3)
$`\{_a(x),(x^{})\}=(x)\delta _{,a}(x,x^{}),`$
$`\{_a(x),_b(x^{})\}=_b(x)\delta _{,a}(x,x^{})+_a(x^{})\delta _{,b}(x,x^{}),`$
and his first class character guarantees that the dynamical trajectories are consistent with the constraints

$$=0,_a=0.$$
(4)
According to Dirac, the quantization of a constrained system requires a factor ordering able to preserve the algebra of constraints at the level of operators (absence of anomalies). In the case of general relativity this issue remains unsolved . In order to avoid the regularization of operators, it is a common practice to study this question in finite-dimensional systems featuring constraints that resemble the algebra of general relativity. Although these types of systems have been widely studied (see, for example, Refs. ), it has been pointed out by Montesimos et al. that an important feature of the algebra (3)-(3) has not been sufficiently examined: usually people does not deal with Eq. (3) but only with (3),(3), i.e. frequently only one super-Hamiltonian is included (however, models with several commuting Hamiltonian constraints were considered in Refs. ). Actually, Montesinos et al. solved a simple model with two euclidean flat super-Hamiltonians plus a single supermomenta satisfying only Eqs. (3)-(3) with constant structure coefficients; thanks to the simple structure of the constraints, the ordering problem is trivially solved in this case. Our aim is to quantize a more complex system in order to deal with a non trival ordering of constraint operators that mimic the complete algebra (3)-(3).

It is well known that, besides the Dirac method, the Becchi-Rouet-Stora-Tyutin (BRST) formalism is a powerful tool to quantize a first class constrained system  (see also Ref.  for an account of the Batalin-Fradkin-Vilkovisky formalism). In recent works , we have taken advantage of its strength for obtaining the consistently ordered constraint operators belonging to generally covariant systems including only one super-Hamiltonian constraint. Here, we will show that the tools there developed are also useful in the treatment of a nontrivial system featuring two super-Hamiltonian constraints an $`m`$ supermomenta constraints.

We will start by defining the model. Then, the system will be quantized within the framework of the BRST formalism, where the nilpotency of the BRST generator must be proven in order to guarantee an anomaly free quantization. The Dirac constraint operators will be identified from the nilpotent BRST generator. In spite of this quantization will be performed for a system featuring an algebra with constant structure functions, the results will be extended to more general algebras by means of a unitary transformation. Finally, a physical inner product will be introduced to complete a consistent quantization procedure.

## II The model

Let us consider a system described by $`4n`$ canonical coordinates $`(q^i,p_i)`$, subjected to two super-Hamiltonian constraints

$$_1=\frac{1}{2}g^{i_1j_1}(q^{k_1})p_{i_1}p_{j_1}+\upsilon _1(q^{k_2}),$$
(5)
$$_2=\frac{1}{2}g^{i_2j_2}(q^{k_2})p_{i_2}p_{j_2}+\upsilon _2(q^{k_1}),$$
(6)
where $`i_1=1,\mathrm{},n`$ and $`i_2=n+1,\mathrm{},2n`$.

The metrics $`g^{i_1j_1}`$ and $`g^{i_2j_2}`$ are indefinite and non degenerated, and depend on the $`q^{i_1}`$’s and the $`q^{i_2}`$ s respectively. On the contrary, the potentials exhibit an opposite functional dependence: $`\upsilon _1=\upsilon _1(q^{i_2})`$ and $`\upsilon _2=\upsilon _2(q^{i_1})`$.

In addition, the system is also subjected to $`m`$ linearly independent supermomentum constraints

$$_a=\xi _a^ip_i,a=3,\mathrm{},m+2;i=(i_1,i_2),$$
(7)
where

$$\xi _a^{i_1}=\xi _a^{i_1}(q^{k_1}),\xi _a^{i_2}=\xi _a^{i_2}(q^{k_2}).$$
(8)
The special way the geometrical objects in the constraint functions depend on the coordinates, has been chosen for obtaining a constraint algebra which mimics Eqs. (3)-(3):

$$\{_1,_2\}=c_{12}^a_a,$$
(9)
$$\{_1,_a\}=c_{1a}^1_1,$$
(10)
$$\{_2,_a\}=c_{2a}^2_2,$$
(11)
$$\{_a,_b\}=C_{ab}^c_c.$$
(12)
It should be noticed that this system describes two interacting particles. In fact, it is not possible to decouple the system in two subsystems described by $`(q^{i_1},p_{i_1})`$ and $`(q^{i_2},p_{i_2})`$. As long as we know, models of this class in riemannian manifolds have not been studied yet in the literature.

We will start by quantizing an algebra with constant structure functions. Later, by taking into account the scaling properties of the super-Hamiltonians, we will extend the results to some algebras with non constant structure functions.

In order that the constraints (5), (6), and (7) effectively satisfy the algebra (9)-(12), metrics, vectors and potentials must fulfill certain relations. By substituting the constraints in Eq. (9) one obtains

$$g^{i_2k_2}\upsilon _{1,k_2}p_{i_2}g^{i_1k_1}\upsilon _{2,k_1}p_{i_1}=c_{12}^a(\xi _a^{i_1}p_{i_1}+\xi _a^{i_2}p_{i_2}),$$
(13)
then,It must be noticed that Eqs. (2.10) and (2.11) do not force $`c_{12}^a`$ to be a constant; they only imply that $`c_{12}^a={}_{}{}^{(1)}c_{12}^{a}(q^{i_1}){}_{}{}^{(2)}c_{12}^{a}(q^{i_2})`$, where $`{}_{}{}^{(1)}c_{12}^{a}\xi _a^{i_1}=0`$ and $`{}_{}{}^{(2)}c_{12}^{a}\xi _a^{i_2}=0`$.

$$c_{12}^a\xi _a^{i_2}=g^{i_2k_2}\upsilon _{1,k_2}$$
(14)
and

$$c_{12}^a\xi _a^{i_1}=g^{i_1k_1}\upsilon _{2,k_1}.$$
(15)
The substitution of the constraints in Eq. (10) yields

$$\frac{1}{2}(g_{,k_1}^{i_1j_1}\xi _a^{k_1}2g^{i_1k_1}\xi _{a,k_1}^{j_1})p_{i_1}p_{j_1}+\upsilon _{1,k_2}\xi _a^{k_2}=c_{1a}^1(g^{i_1j_1}p_{i_1}p_{j_1}+\upsilon _1)$$
(16)
then

$$g_{,k_1}^{i_1j_1}\xi _a^{k_1}2g^{i_1k_1}\xi _{a,k_1}^{j_1}=c_{1a}^1g^{i_1j_1}$$
(17)
and

$$\xi _a^{k_2}\upsilon _{1,k_2}=c_{1a}^1\upsilon _1.$$
(18)
In geometrical language, the relations (17) and (18) read

$`_{\stackrel{}{\xi }_a}\overline{\overline{g_1}}=c_{1a}^1\overline{\overline{g_1}},`$ (19)
$`_{\stackrel{}{\xi }_a}\upsilon _1=c_{1a}^1\upsilon _1.`$ (20)
Analogously, the substitution of the constraints in Eq. (11) leads to

$$\frac{1}{2}(g_{,k_2}^{i_2j_2}\xi _a^{k_2}2g^{i_2k_2}\xi _{a,k_2}^{j_2})p_{i_2}p_{j_2}+\upsilon _{2,k_1}\xi _a^{k_1}=c_{2a}^2(g^{i_2j_2}p_{i_2}p_{j_2}+\upsilon _2)$$
(21)
which means

$$g_{,k_2}^{i_2j_2}\xi _a^{k_2}2g^{i_2k_2}\xi _{a,k_2}^{j_2}=c_{2a}^2g^{i_2j_2}$$
(22)
and

$$\xi _a^{k_1}\upsilon _{2,k_1}=c_{2a}^2\upsilon _2$$
(23)
or, in geometrical language,

$`_{\stackrel{}{\xi }_a}\overline{\overline{g_2}}=c_{2a}^2\overline{\overline{g_2}},`$ (24)
$`_{\stackrel{}{\xi }_a}\upsilon _2=c_{2a}^2\upsilon _2.`$ (25)
Thus, the fulfillment of the algebra means that the supermomenta are conformal Killing vectors of the super-Hamiltonians.

Finally, by substituting the supermomenta in Eq. (12) one obtains that §§§Eq. (2.22) does not imply that the $`C_{ab}^c`$ s are constant, it only imposes the decomposition: $`C_{ab}^c={}_{}{}^{(1)}C_{ab}^{c}(q^{i_1}){}_{}{}^{(2)}C_{ab}^{c}(q^{i_2})`$ where $`{}_{}{}^{(1)}C_{ab}^{c}\xi _c^{i_1}=0`$ and $`{}_{}{}^{(2)}C_{ab}^{c}\xi _c^{i_2}=0`$.

$$(\xi _a^{i_1}\xi _{b,i_1}^{j_1}\xi _b^{i_1}\xi _{a,i_1}^{j_1})p_{j_1}+(\xi _a^{i_2}\xi _{b,i_2}^{j_2}\xi _b^{i_2}\xi _{a,i_2}^{j_2})p_{j_2}=C_{ab}^c(\xi _c^{i_1}p_{i_1}+\xi _c^{i_2}p_{i_2}).$$
(26)
## III The BRST formalism: classical and quantum generators

Our aim is to find the constraint operators that satisfy the algebra at the quantum level. The appropriate factor ordering can be found within the framework of the BRST formalism. For this, the original phase space is extended by including a canonically conjugate pair of fermionic ghosts $`(\eta ^a,𝒫_a)`$ for each constraint function. The central object is the BRST generator, a fermionic function $`\mathrm{\Omega }=\mathrm{\Omega }(q^i,p_j,\eta ^a,𝒫_b)`$ that captures all of identities satisfied by the set of first class constraints in a unique identity

$$\{\mathrm{\Omega },\mathrm{\Omega }\}=0.$$
(27)
The existence of $`\mathrm{\Omega }`$ is guaranteed at the classical level, and $`\mathrm{\Omega }`$ is unique up to canonical transformations in the extended phase space. It can be built by means of a recursive method . In the present case the result is (a closed algebras has “rank” equal to 1)

$$\mathrm{\Omega }=\eta ^1_1+\eta ^2_2+\eta ^a_a+\eta ^1\eta ^2c_{12}^a𝒫_a+\eta ^1\eta ^ac_{1a}^1𝒫_1+\eta ^2\eta ^ac_{2a}^2𝒫_2+\frac{1}{2}\eta ^a\eta ^bC_{ab}^c𝒫_c.$$
(28)
In order to quantize the extended system, the classical BRST generator must be realized as a Hermitian operator. The theory is free from BRST anomalies, if a Hermitian realization of $`\mathrm{\Omega }`$ can be found such that the classical property (27) becomes

$$[\widehat{\mathrm{\Omega }},\widehat{\mathrm{\Omega }}]=2\widehat{\mathrm{\Omega }}^2=0,$$
(29)
i.e., $`\widehat{\mathrm{\Omega }}`$ must be nilpotent. The BRST physical quantum states belong to the set of equivalence classes of BRST-closed states ($`\widehat{\mathrm{\Omega }}\psi =0`$) moduli BRST-exact ones ($`\psi =\widehat{\mathrm{\Omega }}\chi `$) (quantum BRST cohomology).

In order to get the Hermitian and nilpotent operator $`\widehat{\mathrm{\Omega }}`$ is helpful to write $`\mathrm{\Omega }`$ with the canonical and ghost momenta on an equal footing. So, let us adopt the notation

$$\eta ^{C_s}=(q^i,\eta ^\beta ),𝒫_{C_s}=(p_i,𝒫_\beta ),$$
(30)
where $`s=1,0`$, and $`\beta =(A;a)=(1,2;3,\mathrm{},m+2)`$. The index $`s`$ distinguishes original variables from the added ghost variables.

Also, $`\mathrm{\Omega }`$ is a sum of a term quadratic in the momenta

$$\mathrm{\Omega }^{quad}=\eta ^A_A\frac{1}{2}\underset{r,s=1}{\overset{0}{}}\mathrm{\Omega }^{A_rB_s}𝒫_{A_r}𝒫_{B_s}+\eta ^A\upsilon _A,$$
(31)
plus another term linear in the momenta

$$\mathrm{\Omega }^{linear}=\eta ^a_a+\eta ^1\eta ^2c_{12}^a𝒫_a+\eta ^1\eta ^ac_{1a}^1𝒫_1+\eta ^2\eta ^ac_{2a}^2𝒫_2+\frac{1}{2}\eta ^a\eta ^bC_{ab}^c𝒫_c\underset{s=1}{\overset{0}{}}\mathrm{\Omega }^{c_s}𝒫_{c_s}.$$
(32)
In general, it is a difficult step to get the Hermitian nilpotent BRST operator from its classical counterpart because a general method is not available (when not unrealizable: at the quantum level its mere existence is not guaranteed). Nevertheless, let us begin by proposing for the linear term the operator

$$\widehat{\mathrm{\Omega }}^{linear}=\underset{s=1}{\overset{0}{}}f^{\frac{1}{2}}\widehat{\mathrm{\Omega }}^{c_s}\widehat{𝒫}_{c_s}f^{\frac{1}{2}}$$
(33)
where $`f=f(q^i)`$ depends on all of original coordinates.

$`\widehat{\mathrm{\Omega }}^{linear}`$ will be Hermitian if $`f`$ satisfies

$$C_{a\beta }^\beta =f^1(f\xi _a^i)_{,i},$$
(34)
where $`C_{a\beta }^\beta =C_{ab}^b+c_{a1}^1+c_{a2}^2`$ (please remember that $`\beta `$ runs over all constraint labels). In Eq. (34), $`f`$ behaves as a volume in the gauge orbit of the supermomenta (see Ref. ).

The quadratic term can be caught in the Hermitian ordering:

$$\widehat{\mathrm{\Omega }}^{quad}=\frac{1}{2}\underset{r,s=1}{\overset{0}{}}f^{\frac{1}{2}}\widehat{𝒫}_{A_r}f\mathrm{\Omega }^{A_rB_s}\widehat{𝒫}_{B_s}f^{\frac{1}{2}}+\widehat{\eta }^A\upsilon _A.$$
(35)
Finally, the proposed Hermitian BRST operator can be rearranged in the $`\widehat{\eta }\widehat{𝒫}`$ order by repeatedly using the ghost (anti)commutation relations. After this procedure is completed, the classical structure of Eq. (28) will be reproduced at the quantum level, although it will be free of anomalies only if $`\widehat{\mathrm{\Omega }}`$ is nilpotent

$$\widehat{\mathrm{\Omega }}=\widehat{\eta }^1\widehat{}_1+\widehat{\eta }^2\widehat{}_2+\widehat{\eta }^a\widehat{}_a+\widehat{\eta }^1\widehat{\eta }^2c_{12}^a\widehat{𝒫}_a+\widehat{\eta }^1\widehat{\eta }^ac_{1a}^1\widehat{𝒫}_1+\widehat{\eta }^2\widehat{\eta }^ac_{2a}^2\widehat{𝒫}_2+\frac{1}{2}\widehat{\eta }^a\widehat{\eta }^bC_{ab}^c\widehat{𝒫}_c.$$
(36)
In this case, the result is

$$\widehat{}_1=\frac{1}{2}f^{\frac{1}{2}}\widehat{p}_{i_1}g^{i_1j_1}f\widehat{p}_{j_1}f^{\frac{1}{2}}+\upsilon _1,$$
(37)
$$\widehat{}_2=\frac{1}{2}f^{\frac{1}{2}}\widehat{p}_{i_2}g^{i_2j_2}f\widehat{p}_{j_2}f^{\frac{1}{2}}+\upsilon _2,$$
(38)
$$\widehat{}_a=f^{\frac{1}{2}}\xi _a^i\widehat{p}_if^{\frac{1}{2}}.$$
(39)
We still need to demand the nilpotency of the proposed $`\widehat{\mathrm{\Omega }}`$. We expect additional conditions over $`f(q^i)`$ since Eq. (34) does not completely fix it.

Proof of $`\widehat{\mathrm{\Omega }}^2=0`$. The proof is done by explicit calculation. Here we give an abridged demonstration; the full version can be found in Ref. . After the removing of the terms that cancel out due to the Jacobi identities such as

$`\{\{_1,_2\},_a\}+\{\{_2,_a\},_1\}+\{\{_a,_1\},_2\}=0`$ (41)
$`c_{12}^bC_{ba}^cc_{12}^c(c_{1a}^1+c_{2a}^2)=0,`$
or that are identically null \[$`(\eta ^a)^20`$, etc.\], one finally gets

$`[\widehat{\mathrm{\Omega }},\widehat{\mathrm{\Omega }}]=2\widehat{\eta }^1\widehat{\eta }^2([\widehat{}_1,\widehat{}_2]ic_{12}^b\widehat{}_b)`$ (44)
$`+2\widehat{\eta }^1\widehat{\eta }^a([\widehat{}_1,\widehat{}_a]ic_{1a}^1\widehat{}_1)`$
$`+2\widehat{\eta }^2\widehat{\eta }^a([\widehat{}_2,\widehat{}_a]ic_{2a}^2\widehat{}_2).`$
Therefore the nilpotency is apparent whenever the constraint operators $`\widehat{}_1,\widehat{}_2`$, and $`\widehat{}_a`$ realize the first class constraint algebra at the quantum level:

$$[\widehat{}_1,\widehat{}_2]=ic_{12}^a\widehat{}_a,$$
(45)
$$[\widehat{}_A,\widehat{}_a]=ic_{Aa}^{(A)}\widehat{}_{(A)}$$
(46)
(there is no sum over label $`A`$). Then, let us calculate explicitly Eq. (45):

$`[\widehat{}_1,\widehat{}_2]=[{\displaystyle \frac{1}{2}}f^{\frac{1}{2}}\widehat{p}_{i_1}g^{i_1j_1}f\widehat{p}_{j_1}f^{\frac{1}{2}}+\upsilon _1,{\displaystyle \frac{1}{2}}f^{\frac{1}{2}}\widehat{p}_{i_2}g^{i_2j_2}f\widehat{p}_{j_2}f^{\frac{1}{2}}+\upsilon _2]`$ (50)
$`=f^{\frac{1}{2}}[{\displaystyle \frac{1}{2}}f^1\widehat{p}_{i_1}g^{i_1j_1}f\widehat{p}_{j_1}+\upsilon _1,{\displaystyle \frac{1}{2}}f^1\widehat{p}_{i_2}g^{i_2j_2}f\widehat{p}_{j_2}+\upsilon _2]f^{\frac{1}{2}}`$
$`=f^{\frac{1}{2}}[{\displaystyle \frac{1}{2}}g^{i_1j_1}\widehat{p}_{i_1}\widehat{p}_{j_1}{\displaystyle \frac{i}{2}}f^1(g^{i_1j_1}f)_{,i_1}\widehat{p}_{j_1}+\upsilon _1,`$
$`,{\displaystyle \frac{1}{2}}g^{i_2j_2}\widehat{p}_{i_2}\widehat{p}_{j_2}{\displaystyle \frac{i}{2}}f^1(g^{i_2j_2}f)_{,i_2}\widehat{p}_{j_2}+\upsilon _2]f^{\frac{1}{2}}.`$
In order to satisfy Eq. (45) the commutator should not contain quadratic terms in the momenta, so

$$[\frac{i}{2}f^1(g^{i_1j_1}f)_{,i_1}\widehat{p}_{j_1},\frac{1}{2}g^{i_2j_2}\widehat{p}_{i_2}\widehat{p}_{j_2}]\text{and}[\frac{i}{2}f^1(g^{i_2j_2}f)_{,i_2}\widehat{p}_{j_2},\frac{1}{2}g^{i_1j_1}\widehat{p}_{i_1}\widehat{p}_{j_1}]$$

should be zero. This can be achieved by demanding that $`f(q^i)`$ factorizes as

$$f(q^i)=f_1(q^{i_1})f_2(q^{i_2}).$$
(51)
It should be remarked that the additional condition (51) is fully compatible with (34).

As a consequence of the property (51), the term

$$[\frac{i}{2}f^1(g^{i_1j_1}f)_{,i_1}\widehat{p}_{j_1},\frac{i}{2}f^1(g^{i_2j_2}f)_{,i_2}\widehat{p}_{j_2}]$$

is also zero.

Then it remains

$$[\widehat{}_1,\widehat{}_2]=ic_{12}^af^{\frac{1}{2}}\widehat{}_af^{\frac{1}{2}}f_1^1(g^{i_1j_1}\upsilon _{2,j_1}f_1)_{,i_1}+f_2^1(g^{i_2j_2}\upsilon _{1,j_2}f_2)_{,i_2},$$
(52)
where we used the classical relation (13) in order to rebuilt $`\widehat{}_a`$. The last two terms in Eq. (52) can be rewritten by using the classical relations (14) and (15)

$$f_1^1(c_{12}^a\xi _a^{i_1}f_1)_{,i_1}+f_2^1(c_{12}^a\xi _a^{i_2}f_2)_{,i_2}=c_{12}^af^1(\xi _a^if)_{,i}=c_{12}^a(C_{ab}^b+c_{a1}^1+c_{a2}^2),$$
(53)
to cancel out as a consequence of the Jacobi identity (41).

To complete the proof, let us now evaluate the Eq. (46):

$`[\widehat{}_A,\widehat{}_a]=[{\displaystyle \frac{1}{2}}f^{\frac{1}{2}}\widehat{p}_{i_A}g^{i_Aj_A}f\widehat{p}_{j_A}f^{\frac{1}{2}}+\upsilon _A,f^{\frac{1}{2}}\xi _a^k\widehat{p}_kf^{\frac{1}{2}}]`$ (57)
$`={\displaystyle \frac{1}{2}}f^{\frac{1}{2}}(\widehat{p}_{i_A}g^{i_Aj_A}f\widehat{p}_{j_A}\xi _a^k\widehat{p}_k\xi _a^k\widehat{p}_k\widehat{p}_{i_A}g^{i_Aj_A}f\widehat{p}_{j_A})f^{\frac{1}{2}}+i\xi _a^k\upsilon _{A,k}`$
$`={\displaystyle \frac{1}{2}}f^{\frac{1}{2}}\left(i\widehat{p}_{i_A}(g_{,k_A}^{i_Aj_A}\xi _a^{k_A}2g^{i_Ak_A}\xi _{a,k_A}^{j_A})\widehat{p}_{j_A}+fg^{i_Aj_A}[\xi _a^k(lnf)_{,k}]_{,i_A}+\xi _{a,i_Ak}^k\widehat{p}_{j_A}\right)f^{\frac{1}{2}}`$
$`+i\xi _a^k\upsilon _{A,k}.`$
After using Eqs. (17),(18) or (22), (23) one obtains

$`[\widehat{}_A,\widehat{}_a]=ic_{Aa}^A{\displaystyle \frac{1}{2}}f^{\frac{1}{2}}\widehat{p}_{i_A}g^{i_Aj_A}f\widehat{p}_{j_A}f^{\frac{1}{2}}+i\xi _a^k\upsilon _{A,k}`$ (59)
$`+{\displaystyle \frac{1}{2}}f^{\frac{1}{2}}\left(fg^{i_Aj_A}[\xi _a^k(lnf)_{,k}]_{,i_A}+\xi _{a,i_Ak}^k\widehat{p}_{j_A}\right)f^{\frac{1}{2}}.`$
The first two terms in (59) are $`ic_{Aa}^{(A)}\widehat{}_{(A)}`$ (there no sum over label $`A`$). Then we have to see whether the last term is zero. Actually,

$$[\xi _a^k(lnf)_{,k}]_{,i_A}+\xi _{a,i_Ak}^k=[\xi _a^k(lnf)_{,k}+\xi _{a,k}^k]_{,i_A}=[f^1(\xi _a^kf)_{,k}]_{,i_A}=C_{a\beta ,i_A}^\beta =0.$$
(60)
Thus the demonstration is completed. $`\widehat{\mathrm{\Omega }}`$ is Hermitian thanks to the choice (34) and is nilpotent due to the factorization (51).

## IV Unitary transformation and constraint operators

The former results can be generalized, whenever one takes into account transformations leaving the BRST system invariant. In fact, the essential properties of $`\widehat{\mathrm{\Omega }}`$ (Hermiticity and nilpotency) are not modified under a unitary transformation

$$\widehat{\mathrm{\Omega }}e^{i\widehat{C}}\widehat{\mathrm{\Omega }}e^{i\widehat{C}},$$
(61)
This transformation defines a new set of first class constraint operators. By choosing

$`\widehat{C}={\displaystyle \frac{1}{2}}[\widehat{\eta }^1F_1(q^i)\widehat{𝒫}_1\widehat{𝒫}_1F_1(q^i)\widehat{\eta }^1+\widehat{\eta }^2F_2(q^j)\widehat{𝒫}_2\widehat{𝒫}_2F_2(q^j)\widehat{\eta }^2]`$ (64)
$`=\widehat{\eta }^1F_1(q^i)\widehat{𝒫}_1+\widehat{\eta }^2F_2(q^j)\widehat{𝒫}_2+{\displaystyle \frac{i}{2}}F_1(q^i)+{\displaystyle \frac{i}{2}}F_2(q^i)`$
$`=\widehat{𝒫}_1F_1(q^i)\widehat{\eta }^1\widehat{𝒫}_2F_2(q^j)\widehat{\eta }^2{\displaystyle \frac{i}{2}}F_1(q^i){\displaystyle \frac{i}{2}}F_2(q^i),F_{1,2}(q)>0,`$
and using the identities

$$e^{i\widehat{\eta }^AF_A(q^i)\widehat{𝒫}_A}=1+i\widehat{\eta }^A(e^{F_A(q^i)}1)\widehat{𝒫}_A,$$
(65)
$$e^{i\widehat{𝒫}_AF_A(q^i)\widehat{\eta }^A}=1+i\widehat{𝒫}_A(e^{F_A(q^i)}1)\widehat{\eta }^A,$$
(66)
one obtains

$`\widehat{\mathrm{\Omega }}=\widehat{\eta }^1e^{\frac{F_1F_2}{2}}\widehat{}_1e^{\frac{F_1+F_2}{2}}+\widehat{\eta }^1\widehat{\eta }^2ie^{\frac{F_1+F_2}{2}}[\widehat{}_1,e^{F_2}]e^{\frac{F_1+F_2}{2}}\widehat{𝒫}_2`$ (71)
$`+\widehat{\eta }^2e^{\frac{F_2F_1}{2}}\widehat{}_2e^{\frac{F_1+F_2}{2}}+\widehat{\eta }^2\widehat{\eta }^1ie^{\frac{F_1+F_2}{2}}[\widehat{}_2,e^{F_1}]e^{\frac{F_1+F_2}{2}}\widehat{𝒫}_1`$
$`+{\displaystyle \frac{1}{2}}\widehat{\eta }^1\widehat{\eta }^2e^{F_1+F_2}\widehat{c}_{12}^a\widehat{𝒫}_a+\widehat{\eta }^ae^{\frac{F_1+F_2}{2}}\widehat{}_ae^{\frac{F_1+F_2}{2}}`$
$`+\widehat{\eta }^1\widehat{\eta }^a\left(\widehat{c}_{1a}^1ie^{\frac{F_1F_2}{2}}[\widehat{}_a,e^{F_1}]e^{\frac{F_1+F_2}{2}}\right)\widehat{𝒫}_1`$
$`+\widehat{\eta }^2\widehat{\eta }^a\left(\widehat{c}_{2a}^2ie^{\frac{F_1+F_2}{2}}[\widehat{}_a,e^{F_2}]e^{\frac{F_1+F_2}{2}}\right)\widehat{𝒫}_2+{\displaystyle \frac{1}{2}}\widehat{\eta }^a\widehat{\eta }^bC_{ab}^c\widehat{𝒫}_c`$
Thus one can identify constraint operators and structure functions in Eq. (71),

$$\widehat{H}_1=\frac{1}{2}e^{\frac{F_1F_2}{2}}f_1^{\frac{1}{2}}\widehat{p}_{i_1}g^{i_1j_1}f_1\widehat{p}_{j_1}f_1^{\frac{1}{2}}e^{\frac{F_1+F_2}{2}}+e^{F_1}\upsilon _1,$$
(72)
$$\widehat{H}_2=\frac{1}{2}e^{\frac{F_2F_1}{2}}f_2^{\frac{1}{2}}\widehat{p}_{i_2}g^{i_2j_2}f_2\widehat{p}_{j_2}f_2^{\frac{1}{2}}e^{\frac{F_1+F_2}{2}}+e^{F_2}\upsilon _2,$$
(73)
$$\widehat{}_a=e^{\frac{F_1+F_2}{2}}f^{\frac{1}{2}}\xi _a^i\widehat{p}_if^{\frac{1}{2}}e^{\frac{F_1+F_2}{2}}.$$
(74)
The resulting constraint operators, Eqs. (72),(73), correspond to scaled super-Hamiltonian constraints $`H_A=e^{F_A}_A`$. One can name $`G^{i_Aj_A}=e^{F_A}g^{i_Aj_A},V_A=e^{F_A}\upsilon _A`$ to write

$$\widehat{H}_1=\frac{1}{2}e^{\frac{F_1F_2}{2}}f_1^{\frac{1}{2}}\widehat{p}_{i_1}e^{F_1}G^{i_1j_1}f_1\widehat{p}_{j_1}f_1^{\frac{1}{2}}e^{\frac{F_1+F_2}{2}}+V_1,$$
(75)
$$\widehat{H}_2=\frac{1}{2}e^{\frac{F_2F_1}{2}}f_2^{\frac{1}{2}}\widehat{p}_{i_2}e^{F_2}G^{i_2j_2}f_2\widehat{p}_{j_2}f_2^{\frac{1}{2}}e^{\frac{F_1+F_2}{2}}+V_2,$$
(76)
with the corresponding structure functions

$`\widehat{C}_{AB}^B=ie^{\frac{F_1+F_2}{2}}[\widehat{}_A,e^{F_B}]e^{\frac{F_1+F_2}{2}}`$ (78)
$`={\displaystyle \frac{1}{2}}e^{\frac{F_A}{2}}f^{\frac{1}{2}}[\widehat{p}_{i_A}fg^{i_Aj_A}F_{B,j_A}+F_{B,j_A}fg^{i_Aj_A}\widehat{p}_{j_A}]f^{\frac{1}{2}}e^{\frac{F_A}{2}},`$
$$\widehat{C}_{12}^a=e^{\frac{F_1+F_2}{2}}\widehat{c}_{12}^a,$$
(79)
$$\widehat{C}_{Aa}^A=\widehat{c}_{Aa}^A+\xi _a^i(F_A)_{,i},$$
(80)
$$\widehat{C}_{ab}^c=\widehat{C}_{ab}^c.$$
(81)
All the operators and structure functions are ordered in such a way they satisfy,

$$[\widehat{H}_1,\widehat{H}_2]=\widehat{C}_{12}^1(q,p)\widehat{H}_1+\widehat{C}_{12}^2(q,p)\widehat{H}_2+\widehat{C}_{12}^a\widehat{}_a,$$
(82)
$$[\widehat{H}_A,\widehat{}_a]=\widehat{C}_{Aa}^{(A)}\widehat{H}_{(A)},$$
(83)
$$[\widehat{}_a,\widehat{}_b]=\widehat{C}_{ab}^c\widehat{}_c$$
(84)
(there is no sum over label $`A`$), i.e., the algebra is free from anomalies at the quantum level.

Of course, the scaling of constraints does not modify the dynamics of the system. However the constraint algebra looks different: the Eq. (83) now involves all constraints and the structure functions are no longer constant. Since the only consequence of the proposed unitary transformation is the scaling of the super-Hamiltonians, the commutators among supermomenta remain unchanged.

## V Physical inner product

To complete the quantization it is necessary to define a physical inner product where the spurious degrees of freedom are frozen by means of gauge fixing conditions (see also Ref. ). Here it is relevant to take into account the role played by the invariance transformations of the theory: (i) general coordinate transformation, (ii) linear combinations of supermomenta constraints, (iii) scaling of super-Hamiltonian constraints. The physical inner product for Dirac wave functions

$$(\phi _1,\phi _2)=𝑑q\left[\stackrel{m+2}{}\delta (\chi )\right]J\phi _1^{}(q)\phi _2(q)$$
(85)
(where $`J`$ is the Faddeev-Popov determinant and $`\chi `$ represents $`m+2`$ gauge conditions) must be invariant under any of these transformations. By regarding the behavior of the constraint operators under transformations (i-iii) one realizes that the Dirac wave functions should transform according to

$$\phi \phi ^{}=(detA)^{\frac{1}{2}}e^{\frac{F_1+F_2}{2}}\phi ,$$
(86)
$`A`$ being the matrix of the combination of linear constraints. Therefore the Faddeev-Popov determinant $`J`$ in the physical inner product should change in opposite way, in order that the inner product remains unchanged.

In Eq. (85) there are $`m`$ functions $`\chi `$ fixing $`m`$ coordinates associated to the supermomenta constraints, whose characteristics are very well known. So, let us pay attention to the remaining two gauge fixing functions, which come from the super-Hamiltonians and are involved with the scaling transformation. As it was shown in previous works, the scaling factor can be associated either with a positive definite potential (intrinsic time case ) or with the norm of a conformal Killing vector (extrinsic time case ). Up to this point it has not been necessary to make any assumption on the potentials of the super-Hamiltonians (except for their functional dependence). However, some care about the type of potential must be taken to fix the gauge. If the potentials are definite positive, the time is intrinsic and the scaling factor of the super-Hamiltonians will be the potentials themselves: $`F_A=lnV_A`$ and the associated gauge condition will be the one studied in Ref. . If the time is extrinsic, the scale factors are associated with the norm of a conformal Killing vector for each super-Hamiltonian: $`F_A=ln|\stackrel{}{\xi }_A|^2`$ and the inner product will be the one defined in Ref. .

However, there are $`m+2`$ constraints in the theory, so we need one more gauge condition. Since, this finite dimensional model mimics the constraint algebra of general relativity, the two subspaces $`(q^{i_1})`$ and $`(q^{i_2})`$ are interpreted as the field at two different points of the space-time, and so the super-Hamiltonians $`_1`$ and $`_2`$ are regarded as “the super-Hamitonian” evaluated in two diferent points in space-time. Then, we suppose that both super-Hamiltonian have the same type of time, and in fact, both contain the same time. So, we propose as the remaining gauge condition

$$\chi =\delta (t_1t_2).$$
(87)
This gauge condition can be naturally retrieved in the framework of the multiple time formalism . Therefore, in both cases -intrinsic and extrinsic time- the inner product is regularized by including Eq. (87) in the already known set of $`m+1`$ gauge conditions.

## VI Conclusions

In this work we managed to naturally extend the previous ordering findings of Refs.  to the interesting case of systems subject to more than one super-Hamiltonian constraint. As long as we know, this is the first work that treat ordering problems in such case and it also prove the power of the BRST formalism for providing the necessary tools for the treatment of invariance properties at the quantum level. By taking advantage of these invariance transformations of the theory, in particular the scaling invariance of the super-Hamiltonians, we were able to raise a nontrivial algebra between the super-Hamiltonians and to find the anomaly-free ordering, namely Eqs. (82)-(84). In this case, each scaling independently contribute with a term, involving also the momenta constraints. The role played by the invariance transformations, and how they modify the operator ordering, has been briefly discussed in Sec. V (see also Refs.  for a more in depth analysis).

When a finite dimensional system is proposed as a model for quantizing a general covariant theory (such as for example, general relativity), always remains the suspicion about the real competence of the model when it is compared with the infinite dimensional case. However, we can learn something about the ordering of the constraint operators that must be applied to both cases: the invariance transformation of the theory substantially modify the ordering. For example, if one admits that the super-Hamiltonian in the infinite dimensional case can be scaled, then the ordering for the supermomenta cannot be a simple functional derivation with respect to the canonical field variable, as usually appears in the literature, but it must “wear” the scaling factors in both sides of it.

###### Acknowledgements.

This research was supported by Universidad de Buenos Aires (Project No. TX 64) and Consejo Nacional de Investigaciones Científicas y Técnicas. D. M. S. thanks the members of the Department of Physics and Astronomy of the University of Pittsburgh, specially Roberto Gomez, for kindly providing facilities which made possible to finish this work.


text: # References