text:
# Complexity of decoupling and time-reversal for 𝑛 spins with pair-interactions: Arrow of time in quantum control

## Abstract

Well-known Nuclear Magnetic Resonance experiments show that the time evolution according to (truncated) dipole-dipole interactions between $`n`$ spins can be inverted by simple pulse sequences. Independent of $`n`$, the reversed evolution is only two times slower than the original one. Here we consider more general spin-spin couplings with long range. We prove that some are considerably more complex to invert since the number of required time steps and the slow-down of the reversed evolutions are necessarily of the order $`n`$. Furthermore, the spins have to be addressed separately. We show for which values of the coupling parameters the phase transition between simple and complex time-reversal schemes occurs.

Due to wide applications of Nuclear Magnetic Resonance in medicine, chemistry, and physics much effort has been undertaken to understand the dynamics of nuclear spins in solids and liquids . Nowadays, NMR plays an important role in the first experimental realizations of quantum computation .

Phenomenologically, the dynamical evolution of a single spin seems to be a relatively simple relaxation and dephasing process of an open quantum system, formally described by the Bloch equations . However, the fact that not all dephasing and relaxation processes are irreversible is important for practical purposes since refocusing techniques are applied with great success .

The simplest version of refocusing is possible when dephasing of the (uncoupled) spins is caused by spatial inhomogeneities of the strength of the magnetic field. In this case the time evolution can be reversed by sandwiching the natural dynamics by $`180^{}`$–rotations of all spins around an axis orthogonal to the magnetic field. More sophisticated versions of refocusing are necessary if dephasing and relaxation processes are caused by spin-spin interactions. An important example is the so-called dipole-dipole coupling .

We describe this interaction formally as follows. Let $`(^2)^n`$ be the Hilbert space of $`n`$ spins and $`\sigma _\alpha ^{(k)}`$ be the Pauli matrix $`\sigma _\alpha `$ acting on the $`k`$-th spin for $`1kn`$ and $`\alpha =x,y,z`$. Then the so-called truncated dipole-dipole Hamiltonian $`H_d`$ is given by

$`H_d:={\displaystyle \underset{k<l}{}}w_{kl}\left({\displaystyle \underset{\alpha }{}}\sigma _\alpha ^{(k)}\sigma _\alpha ^{(l)}3\sigma _z^{(k)}\sigma _z^{(l)}\right),`$

where $`w_{kl}`$ is the strength of the interaction between the spins $`k`$ and $`l`$.

We assume that all local unitary operations on the spins can be produced in arbitrarily small time (on a time scale governed by the strength of the coupling). This is often called the fast control limit. In practice this is done by external electro-magnetic pulses.

Consider a pulse implementing a rotation of each spin around the $`y`$-axis by $`90^{}`$, i. e., if the spin is in $`z`$-direction it is in $`x`$-direction after the rotation. Formally, one applies a unitary transformation on $`(^2)^n`$ of the form

$`v_y:=u_yu_y\mathrm{}u_y`$

with $`u_y^{}\sigma _zu_y=\sigma _x`$. Now we focus on the effect of the following scheme: (1) implement $`v_y`$ by the corresponding pulse, (2) wait the time $`t`$, and (3) implement $`v_y^{}`$, i. e., the inverse of $`v_y`$. The net effect of the scheme is an evolution

$`v_y^{}\mathrm{exp}(iH_dt)v_y=\mathrm{exp}(iv_y^{}H_dv_yt),`$

i. e., a dynamics according to the Hamiltonian $`v_y^{}H_dv_y`$. We choose a second rotation $`u_x`$ around the $`x`$-axis such that $`u_x^{}\sigma _zu_x=\sigma _y`$ and set $`v_x:=u_xu_x\mathrm{}u_x`$. It is decisive to note that

$`H_d=v_y^{}H_dv_y+v_x^{}H_dv_x.`$

This has the interesting consequence that for small times $`ϵ`$ one has

$`v_x^{}\mathrm{exp}(iH_dϵ)v_xv_y^{}\mathrm{exp}(iH_dϵ)v_y\mathrm{exp}(+iH_dϵ)`$

up to terms of second order of $`ϵ`$. This first order approximation is usually referred to as average Hamiltonian theory . The intuitive meaning is that the time evolution according to $`H_d`$ can be approximatively inverted by dividing the time $`t`$ into small intervals of length $`ϵ`$ and interspersing the natural time evolution by the above pulses. The inverted dynamics is two times slower than the natural one. We express it by saying that the time overhead of inverting the Hamiltonian $`H_d`$ is $`2`$. The task of time-reversal is a special case of the general problem of simulating Hamiltonians that has recently been considered by several authors .

Note that the decoupling and inversion schemes for the dipole-dipole Hamiltonian is simple with respect to three aspects: (a) it does not use selective pulses, i.e., in each time step the same unitary transformation is applied to all spins. (b) The number of steps of this inversion scheme does not increase with $`n`$ since it is always $`2`$, and (3) the time overhead does not increase with $`n`$ since it is also $`2`$.

The purpose of this article is to show that there are interactions that are considerably more complex to invert with respect to all three criteria. The most general spin-spin interaction on $`n`$ spins is given by

$`H_J:={\displaystyle \underset{k<l}{}}{\displaystyle \underset{\alpha ,\beta }{}}J_{kl;\alpha \beta }\sigma _\alpha ^{(k)}\sigma _\beta ^{(l)},`$

where $`J`$ is chosen to be a real symmetric $`3n\times 3n`$-matrix with zeros for $`k=l`$. Note that the symmetry of the coupling matrix $`J`$ does not imply any physical symmetry of the interaction. It is a consequence of our redundant notation that turns out to be very useful. The coupling matrix $`J`$ consists of $`3\times 3`$-blocks. The $`3\times 3`$-matrix $`J_{kl}`$ given by the block at position $`(k,l)`$ describes the coupling between the spins $`k`$ and $`l`$. We have $`J_{lk}=J_{kl}^T`$, i. e. the matrix describing the coupling between the spins $`l`$ and $`k`$ is just the transpose of the matrix describing the coupling between the spins $`k`$ and $`l`$. The blocks on the diagonal are zero matrices.

Let $`𝒞:=SU(2)SU(2)\mathrm{}SU(2)`$ be the control group of local operations on the spins. As already mentioned all operations in $`𝒞`$ are assumed to be implemented arbitrarily fast. The most general refocusing scheme (also called time-reversal or inversion) based on the first order approximation above is given by rotations $`v_1,v_2,\mathrm{},v_N𝒞`$, and (relative) times $`t_1,t_2,\mathrm{},t_N`$ such that the average Hamiltonian

$`\overline{H}:={\displaystyle \underset{j}{}}t_jv_j^{}H_Jv_j`$

is equal to $`H_J`$. Here $`N`$ is the number of time steps and $`\tau :=_{j=1}^Nt_j`$ is the time overhead. The time overhead gives the slow-down of the inverted Hamiltonian. In the following we consider both the number of time steps and the time overhead as complexity measures. In contrast, if the sequences are chosen in such a way that the average Hamiltonian is zero one has a decoupling scheme. Each inversion scheme with $`N`$ steps defines a decoupling scheme with $`N+1`$ steps in a straightforward way (note that the notion of time overhead does not make sense for decoupling). Conversely, each decoupling scheme can be converted to an inversion scheme as follows. The equation

$`{\displaystyle \underset{j=0}{\overset{N}{}}}t_jv_j^{}H_Jv_j=0`$

implies

$`{\displaystyle \underset{j=1}{\overset{N}{}}}(t_j/t_0)(v_jv_0^{})^{}H_J(v_jv_0^{})^{}=H_J`$

by elementary calculation. Therefore we shall restrict our attention to inversion.

The condition on inversion schemes can be rewritten by expressing the action of the control operations on the coupling matrix $`J`$ as follows. Note that any unitary operation $`uSU(2)`$ corresponds to a rotation on the $`3`$-dimensional Bloch sphere via the relation

$`u^{}\left({\displaystyle \underset{\alpha }{}}c_\alpha \sigma _\alpha \right)u={\displaystyle \underset{\alpha }{}}\stackrel{~}{c}_\alpha \sigma _\alpha ,`$

where the vector $`\stackrel{~}{c}=(\stackrel{~}{c}_x,\stackrel{~}{c}_y,\stackrel{~}{c}_z)`$ is obtained by applying a rotation $`USO(3)`$ on the vector $`c=(c_x,c_y,c_z)`$. It is straightforward to verify that conjugation of $`H_J`$ by $`v:=u^{(1)}u^{(2)}\mathrm{}u^{(n)}`$ corresponds to conjugation of $`J`$ by a block diagonal matrix of the form

$`V:=U^{(1)}U^{(2)}\mathrm{}U^{(n)}{\displaystyle \underset{k=1}{\overset{n}{}}}SO(3).`$

The condition for time inversion is hence given by

$$J=\underset{j}{}t_jV_jJV_j^T,$$
(1)
where the orthogonal matrix $`V_j`$ corresponds to the unitary transformation $`v_j`$ for $`j=1,\mathrm{},N`$. Time reversal and decoupling schemes that apply for general coupling $`J`$ have been presented in . These schemes can even be used when the coupling $`J`$ is unknown. Here we focus on optimality criteria for schemes referring to specific couplings and show that there are interactions that cannot be inverted significantly more efficiently than an unknown interaction.

We restrict our attention to interactions with an additional symmetry, namely a Hamiltonian of the following form

$`H:={\displaystyle \underset{k<l}{}}w_{kl}{\displaystyle \underset{\alpha ,\beta }{}}A_{\alpha \beta }\sigma _\alpha ^{(k)}\sigma _\beta ^{(l)}.`$

The matrix $`W:=(w_{kl})`$ is a real symmetric $`n\times n`$-matrix with zeros on the diagonal. It describes the coupling strengths and the signs of the interactions between all spins. The matrix $`A=(A_{\alpha \beta })`$ is a real symmetric $`3\times 3`$-matrix characterizing the type of the coupling. This means that all spins interact with each other via the same interaction and that only the coupling strength and the sign varies. It is important that in this special case the coupling matrix $`J`$ can be expressed as a tensor product of $`W`$ and $`A`$, i. e. $`J=WA`$.

The matrix $`W`$ characterizes the coupling topology of the spins. In the language of graph theory the matrix $`W`$ is the adjacency matrix of a weighted graph. Therefore many results on graph spectra (eigenvalues of the adjacency matrices) can be used to derive lower and upper bounds on the number of time steps and the time overhead.

To discuss the complexity aspects we distinguish between the following three cases:

1. $`A`$ is traceless.

All spins can be subjected to the same transformations in each time step, the number of time steps and the time overhead are at most $`2`$.
2. $`A`$ has negative and positive eigenvalues but $`tr(A)0`$.

The spins have to be addressed separately, the number of time steps necessarily grows for increasing $`n`$. But the time overhead does not depend on $`n`$. It depends only on the eigenvalues of $`A`$.
3. $`A`$ is either positive or negative semidefinite, i. e., the non-zero eigenvalues have the same sign.

Then the spins have to be addressed separately, the number of time steps is at least $`n1`$, and the time overhead is also at least $`n1`$.

To prove these statements we assume w.l.o.g. that the interaction between all pairs is of the form

$`a_x\sigma _x\sigma _x+a_y\sigma _y\sigma _y+a_z\sigma _z\sigma _z,`$

where $`a_x,a_y,a_z`$ are the eigenvalues of $`A`$. This can always be achieved by rotating the reference frame .

Case 1. Assume $`A`$ to be traceless. Let $`S`$ be a rotation in the Bloch sphere that realizes the cyclic permutation of the axis according to $`xyzx`$. Then we have

$`{\displaystyle \underset{j=0}{\overset{2}{}}}V^jJV_{}^{j}{}_{}{}^{T}=\mathrm{𝟎},`$

where $`V:=SS\mathrm{}S`$ and $`V^j`$ is the $`j`$th power of $`V`$. Here $`S`$ is a rotation in the Bloch sphere that realizes the cyclic permutation of the axis according to $`xyzx`$. The $`j`$th power of $`V`$ is denoted by $`V^j`$. Hence the Hamiltonian is inverted by a sequence of length $`2`$. Due to the equation

$`J=VJV^T+V^2JV_{}^{2}{}_{}{}^{T}`$

the time overhead of this inversion scheme is $`2`$.

Case 2. The fact that the spins have to be addressed separately has also been noted in . It can be seen as follows. If all spins are subjected to the same transformation each $`3\times 3`$–block of $`J`$ is conjugated by an element of $`SO(3)`$. Such a conjugation preserves the trace of each block matrix. This trace is even preserved under positive linear combination, i.e. the trace of each resulting $`3\times 3`$ matrix has the same sign as the trace of the original matrix. Therefore a sequence applying the same transformation to all spins can never change the sign of the trace in any block. But this is required for obtaining the inverse entries.

We show that the number of time steps cannot be constant, i. e. independent of $`n`$. Assume w.l.o.g. that $`tr(A)>0`$. Choose a partition of $`SO(3)=S_1S_2\mathrm{}S_p`$ into equivalence classes such that for any pair $`O,\stackrel{~}{O}S_i`$ ($`i=1,\mathrm{},p`$) implies $`tr(OA\stackrel{~}{O}^T)0`$ ($`p`$ depends on $`A`$). This is always possible since the partition can be chosen in such a way that all orthogonal transformations in the same class are sufficiently close together. Then $`OA\stackrel{~}{O}^T`$ almost corresponds to a conjugation of $`A`$ by an orthogonal matrix and its trace cannot deviate too much from $`tr(A)`$ due to the continuity of the trace.

Let $`O_j^{(k)}SO(3)`$ be the operation performed on the $`k`$th spin in time step $`j`$ (described as action on the coupling matrix). In each time step at least $`n/p`$ spins are transformed by operations of the same class***Here $`x`$ denotes the smallest integer $`m`$ such that $`xm`$.. After $`N`$ time steps at least $`n/p^N`$ spins have been transformed by operations of the same class in each step. Let $`k`$ and $`l`$ be such a pair that has been transformed by operations of the same equivalence class in each time step. Then the trace of the average coupling between the spins $`k`$ and $`l`$ after the sequence is given by

$$tr(\underset{j}{}t_jO_j^{(k)}AO_{j}^{(l)}{}_{}{}^{T})0.$$
(2)
If the sequence is an inversion scheme then we must have

$`{\displaystyle \underset{j}{}}t_jO_j^{(k)}AO_{j}^{(l)}{}_{}{}^{T}=A.`$

But this condition contradicts the inequality (2). Therefore such a pair $`k`$ and $`l`$ must not exist. This gives the lower bound $`N\mathrm{log}n/\mathrm{log}p`$ on the number of time steps $`N`$.

We show that the time overhead does not depend on $`n`$. We assume w.l.o.g. $`a_x>0>a_z`$ for the eigenvalues of $`A`$. First we describe a partial decoupling scheme selecting for instance the $`\sigma _z\sigma _z`$ coupling terms while switching off the $`xx`$ and $`yy`$ terms. Following  such a decoupling can be achieved by certain sequences of local conjugations by the unitary $`i\sigma _z`$. Each time step of this scheme is described by a column of a Hadamard matrix. The entries determine which spins are conjugated by $`i\sigma _z`$ transformations. The idea is that the $`xx`$ and $`yy`$ terms acquire in exactly half of the time steps a minus sign. Note that this scheme does not weaken the $`zz`$ coupling. In other words, although this scheme requires a number of time steps of the order $`n`$ the time overhead for the simulation of the interaction $`w_{kl}a_z\sigma _z^{(k)}\sigma _z^{(l)}`$ by the original one is $`1`$.

There are two cases depending on whether $`a_y0`$ or $`a_y<0`$. Consider the case $`a_y0`$ w.l.o.g. First switch off the $`xx`$ and $`yy`$ interactions. Apply on each spin a conjugation by $`u_y`$ (defined as in the beginning). This simulates the interaction

$`{\displaystyle \underset{k<l}{}}w_{kl}a_z\sigma _x^{(k)}\sigma _x^{(l)}.`$

By applying this interaction for the (relative) time $`a_x/a_z`$ we have reversed the $`xx`$-components in the original Hamiltonian. Since $`a_y`$ is not negative we can reverse similarly the $`yy`$-components with time overhead $`a_y/a_z`$. To invert the $`zz`$-components we use the $`xx`$-components. The time overhead is $`a_z/a_x`$. In summary we see that the time overhead is independent of $`n`$. It only depends on the eigenvalues of $`A`$.

Case 3. Let $`t_1,t_2,\mathrm{},t_N`$ be the (relative) times and $`V_1,V_2,\mathrm{},V_N`$ with $`V_j_{k=1}^nSO(3)`$ be the operations of an inversion scheme, i. e.,

$$\underset{j}{}t_jV_j(WA)V_j^T=WA.$$
(3)
First note that an inversion scheme for a given Hamiltonian $`H_J`$ can also be applied to a rescaled Hamiltonian obtained by changing the strength and the sign of an arbitrary spin pair interaction. If all the coefficients $`w_{kl}`$ are non-vanishing for $`kl`$ then we can equivalently describe a time inversions scheme for the coupling matrix $`J:=\stackrel{~}{W}A`$ where all non-diagonal entries of $`\stackrel{~}{W}`$ are $`1`$ and the diagonal consist of zeros. Hence the fact that interactions decrease with distance of the spins is irrelevant for considerations of the complexity of time inversion (as long as the interaction cannot be neglected). Therefore we assume w.l.o.g. all non-diagonal entries of $`W`$ to be $`1`$.

For each time step let $`V_j`$ be the block diagonal matrix $`V_j=U_j^{(1)}\mathrm{}U_j^{(n)}`$. Set $`R:=_{j=1}^Nt_jV_j(\mathrm{𝟏}A)V_j^T`$. We add $`R`$ to both sides of eq. (3) and obtain

$$\underset{j=1}{\overset{N}{}}t_jV_j\left((W+\mathrm{𝟏})A\right)V_j^T=WA+R.$$
(4)
The rank of the matrix $`(W+\mathrm{𝟏})`$ is $`1`$ since all its entries are $`1`$. Consequently, the rank of the left hand side of equation (4) is at most $`Nr`$, where $`r`$ is the rank of $`A`$. Note that $`R`$ is a semi-positive matrix. Therefore the rank of the right hand side is at least the number of positive eigenvalues of $`WA`$. The matrix $`W`$ has $`n1`$ positive eigenvalues (by adding the identity to $`W`$ one easily verifies that $`W`$ has the eigenvalues $`n1,1,1,\mathrm{},1`$). Hence $`WA`$ has $`(n1)r`$ positive eigenvalues. We conclude that the rank of the right hand side is at least $`(n1)r`$. This proves the statement $`Nn1`$.

We show now that the time overhead also grows with $`n`$. Let $`\lambda _{\mathrm{min}}`$ and $`\lambda _{\mathrm{max}}`$ be the minimal and maximal eigenvalues of $`J`$, respectively. Note that $`\lambda _{\mathrm{min}}`$ is negative since $`J`$ is traceless. Let $`t_1,t_2,\mathrm{},t_N`$ be the (relative) lengths of the time steps. Then we have

$`{\displaystyle \underset{j}{}}t_j\lambda _{\mathrm{max}}\lambda _{\mathrm{min}},`$

since the maximal eigenvalue of the sum of matrices in eq. (1) is at most the sum of their maximal eigenvalues  and $`\lambda _{\mathrm{min}}`$ is the maximal eigenvalue of $`J`$. We conclude that the time overhead $`\tau :=_jt_j`$ for the inverse evolution is at least

$`\tau \lambda _{\mathrm{max}}/\lambda _{\mathrm{min}}.`$

Since the maximal and minimal eigenvalues of $`W`$ are $`(n1)`$ and $`1`$, respectively, the maximal and minimal eigenvalues of $`J`$ are $`(n1)a_{\mathrm{max}}`$ and $`a_{\mathrm{max}}`$ where $`a_{\mathrm{max}}`$ is the maximal eigenvalue of $`A`$. We conclude that the minimal time overhead is at least $`n1`$.$`\mathrm{}`$

The truncated dipole-dipole coupling is an example for case 1. The corresponding matrix $`A`$ is given by $`A:=diag(1,1,2)`$. An important example for case 3 is the strong scalar coupling where $`A`$ is the identity matrix. By combining both types one can easily find examples for case 2 .

Besides the practical applications of refocusing the lower bounds on the inversion complexity may be seen as a small contribution to the old debate on the origin of phenomenological irreversibility: Whereas the considered dynamical evolutions are strictly reversible the complexity of a real implementation of the time inversion is growing with the number of interacting particles. Note that the relevance of complexity theory for the definition of physical entropy has been advocated for by several authors . Of course our complexity bounds are restricted to the assumptions of average Hamiltonian theory. Inversion schemes that are not based on this approach and include higher order terms of the Magnus expansion may be less complex. However, our results present an example how to develop a complexity theory of time reversal in physics.


format: text

text:
# On the metal abundance in the core of M87

## 1 Introduction

Our proximity to the Virgo cluster allows us to study its core on physical scales unachievable in other similar systems. M87, the giant elliptical galaxy at the center of Virgo, is not surprisingly amongst the most studied extra-galactic objects in the sky. At X-ray wavelengths it has been observed by all major observatories. Deprojection of the X-ray images obtained with the IPC instrument on Einstein (Fabian et al. 1984) showed that the gas in M87 has cooling times shorter than the Hubble time. Further evidence of cooling gas has come from the spectrometers on-board Einstein, which showed that the gas in M87 was not single-temperature (e.g. Canizares et al. 1982). Subsequent observations with ROSAT (Nulsen & Böhringer 1995), ASCA (Matsumoto et al. 1996) and BeppoSAX (Guainazzi & Molendi 2000), showed evidence of a temperature decrement with decreasing radius. The ASCA and BeppoSAX data also showed evidence of a metal abundance increase towards the center of M87 on scales of a few arcmin, which correspond to a few tens of kpc.

A first set of results from the XMM-Newton observation of M87 have been recently published in Böhringer et al. (2001) and Belsole et al. (2001). One of the most striking results presented is a very marked drop in the abundance profiles of $`Si`$ and $`Fe`$ within 1 arcmin of the center of M87. Böhringer et al. (2001) suggest that resonant scattering may be responsible for this drop. Mathews et al. (2001), using detailed radiation transfer calculations, remark that the resonant scattering effect is insufficient to explain the observed profile and that some sort of continuous opacity is required.

In a recent paper (Molendi & Pizzolato 2001) we have shown that the multi-phase appearance of the EPIC spectra of M87 within 2 arcmin from the core is not ascribable to the presence of a truly multi-phase gas but is rather the consequence of having, within the same line of sight, emission from gas at different radii with different temperatures. In this Letter we revisit the XMM-Newton observation of M87; our goal is to remeasure the metal abundance profile by performing a more advanced spectral analysis of the EPIC data with models which take into account the projection effects leading to the multi-phase appearance of the gas in M87.

The remainder of this Letter is organized as follows. In section 2 we describe the data preparation. In section 3 we discuss the spectral models. In section 4 we present and discuss the results of our analysis.

## 2 Observations and Data Preparation

We use XMM-Newton EPIC data from the PV observation of M87/Virgo. Details on the observation, as well as results from a first analysis of this object, have already been published in Böhringer et al. (2001) and Belsole et al. (2001).

We have obtained calibrated event files for the MOS1, MOS2 and PN cameras with SASv5.0. Data were manually screened to remove any remaining bright pixels or hot column; for the PN camera we also exclude regions contaminated by out of time events. Periods in which the background is increased by soft-proton flares have been excluded using an intensity filter; we rejected all events accumulated when the count rate exceeds 15 cts/100s in the $`[1012]`$ keV band for the MOS cameras and 25 for the PN. We have accumulated spectra in 4 concentric annular regions centered on the emission peak. We extend our analysis out to 3 arcmin from the emission peak. The bounding radii for the regions are: 0′′ and 30′′ for region 1; 30′′ and 1 for region 2; 1 and 2 for region 3 and 2 and 3 for region 4. We have removed point sources, and the substructures which are clearly visible from the X-ray image (e.g. Belsole et al. 2001) except in the innermost region. Here we have kept the nucleus and knot A, the reason being that on angular scales this small it is not possible to exclude completely their emission. As discussed in section 3 we prefer to fit the spectrum of this region with a model which includes a power-law component to fit the two point-like sources; note that we include only one power-law component because the two sources have similar spectra (Böhringer et al. 2001).

Spectra have been accumulated for MOS1 and MOS2 independently. The Lockman hole observations have been used for the background. Background spectra have been accumulated from the same detector regions as the source spectra. The line of sight galactic absorbing column depth towards the background field is somewhat smaller, $`N_\mathrm{H}6\times 10^{19}`$ cm-2, than that found on the line of sight of M87. This has a negligible impact on the soft part of our M87 spectra because the $`N_\mathrm{H}`$ variation will lead to a maximum difference of 10% in the background intensity, which is in turn never more than 2% of the source intensity in the 0.5-1.0 keV energy range.

The vignetting correction has been applied to the spectra rather than to the effective areas, similarly to what has been done by other authors who have analyzed EPIC data (Arnaud et al. 2001). Spectral fits were performed in the 0.5-4.0 keV band. Data below 0.5 keV were excluded to avoid residual calibration problems in the MOS response matrices at soft energies. Data above 4 keV were excluded because above this energy the spectra show a substantial contamination from hotter gas emitted further out in the cluster, on the same line of sight.

MOS and PN spectra are analyzed separately. As discussed in Molendi (2001), there are residual cross-calibration uncertainties between the two instruments which, as we shall see in the next sections, sometimes lead to different measurements of spectral parameters.

## 3 Spectral Modeling

All spectral fitting has been performed using version 11.0.1 of the XSPEC package. All models discussed below include a multiplicative component to account for the galactic absorption on the line of sight of M87. The column density is always fixed at a value of 1.8$`\times 10^{20}`$ cm-2, which is derived from 21cm measurements (Lieu et al. 1996). We have compared our data with three different spectral models.

Model I is a single temperature model (vmekal model in XSPEC), which allows to fit separately individual metals abundances. This model has 13 free parameters: the temperature $`T`$, the normalization and the abundance of $`O`$, $`Ne`$, $`Na`$, $`Mg`$, $`Al`$, $`Si`$, $`S`$, $`Ar`$, $`Ca`$, $`Fe`$ and $`Ni`$, which are all expressed in solar units.

Model II features two single temperature components (vmekal + vmekal in XSPEC) This model has two free parameters more than model I: the temperature and normalization of the second thermal component. The metal abundance of each element of the second thermal component is linked to the same parameter of the first thermal component. This model has been used in the past (e.g. Makishima et al. 2001 and refs. therein) as an alternative to cooling-flow models.

Model III includes a single temperature plus a multi-phase component (vmekal + vmcflow in XSPEC). This model has two free parameters more than model I: the minimum temperature, $`T_{\mathrm{min}}`$, and the normalization of the multi-phase component. The other parameters of the multi-phase component are not free: the maximum temperature, $`T_{\mathrm{max}}`$, and the metal abundance of each element are linked respectively to the temperature, $`T`$, and the metal abundances of the same elements of the single-phase component. In a recent paper (Molendi & Pizzolato 2001) we have shown that such a model fits adequately the EPIC spectra of three clusters which feature a temperature decrement in their core. The vmcflow component we use has originally been written to fit spectra from multi-phase gas, i.e. gas that at a given physical radius is characterized by a distribution of temperatures rather than by a single temperature. Here we are using it to describe a different physical scenario where, at a given physical radius the gas is all at one temperature, and the multi-phase appearance of the spectrum comes from having emission from many different physical radii all within the same line of sight.

For the spectrum accumulated in the innermost region models II and III also include a power-law component to model the emission of the nucleus and of knot A.

## 4 Results

In Figure 1 we plot the $`\mathrm{\Delta }\chi ^2`$ between models I and II (filled circles) and models I and III (open circle) as a function of radius. The solid (dashed) horizontal line indicates the $`\mathrm{\Delta }\chi ^2`$ value for which the statistical improvement of the model II (III) fit with respect to model I is significant at the 99% level according to the F-test. The results plotted in the top panel refer to the MOS spectra while those shown in the bottom panel are for the PN spectra. From Figure 1 we see that, for all bins and for both MOS and PN, the improvement of models II and III with respect to model I is significant at more than the 99% level, implying that the spectra are characterized by more than just one temperature.

In the case of the innermost bin we have verified that including only the second thermal component for models II and III, provides a significantly better fit with respect to model I and that the further inclusion of the power-law component is also statistically significant. Thus for this bin we have evidence in favor of a second thermal component and of a power-law component.

In Figures 2 and 3 we report respectively the MOS and PN radial abundance profiles of $`O`$ (top panel), $`Si`$ (middle panel) and $`Fe`$ (bottom panel), while in Figures 4 and 5 we show those for $`Mg`$ (top panel), $`S`$ (middle panel) and $`Ar`$ (bottom panel). Full symbols indicate the measurements with model I, empty squares and triangles those with models II and III respectively.

Figures 2, 3, 4 and 5 show that, for the two innermost bins, all the metal abundances measured with model I, either with MOS or PN, are systematically smaller than those measured with models II or III. The difference between the measurements becomes larger for smaller radii as the significance of the improvement of models II and III with respect to model I (see Figure 1) increases.

We note that, while for $`Si`$, $`Fe`$, $`S`$ and $`Ar`$ the measurements obtained with PN and MOS are consistent, for $`O`$ and $`Mg`$ abundances derived with the MOS are systematically higher than those derived with the PN. The reason lies in the current cross-calibration uncertainties between the spectral response of the two instruments (Molendi 2001). However, the shape of the abundance profiles derived with models II and III is the same for MOS and PN: neither shows evidence of a dramatic metal abundance drop such as is observed with model I.

Clearly single temperature models substantially underestimate the metal abundance in multi-temperature spectra. This is not a new result, indeed something very similar has been observed in the ASCA spectra of galaxies and galaxy groups (e.g. Buote 2000a,b) and goes under the name of “$`Fe`$-bias”. The interested reader can find a detailed discussion of the “$`Fe`$-bia” in Appendix A of Buote (2000a). The residuals of a single temperature model applied to a multi temperature spectrum around the Fe-L line complex (see fig. 5 of Buote 2000b) are also quite illustrative.

Figure 1 tells us that models II and III behave very similarly in terms of their improvement with respect to model I. In the outermost radial bin, where we find the smallest evidence for emission from more than a single temperature, the abundances measured with models II and III do not differ substantially from those derived from model I. As we move inwards the evidence for multi-temperature gas becomes stronger, the $`Fe`$-bias becomes more important, and the abundance measurements with models II and III differ more and more from those obtained with model I. In the innermost bin, the power-law emission from the nuclear sources contributes to raising the continuum with respect to the line emission thereby driving the abundances obtained with model I further down.

Since the metal abundances in the innermost region are sensitive to the normalization of the power-law component we have performed the following test to verify that this parameter is sufficiently well determined by the fitting procedure. We have accumulated a spectrum for the innermost bin excising the emission from the nucleus and knot A, using an exclusion radius of 7′′ around the sources. We recall that for such a small exclusion radius the accumulated spectrum will still be contaminated by the emission of the point-sources. We have then applied models II and III to this spectrum. We find that the best fitting values for the metallicities are in agreement with those derived by fitting the spectrum where the two sources have not been eliminated.

Close inspection of the $`O`$, $`Si`$, $`Fe`$ and $`S`$ profiles reveals a small decrement in the innermost bin, the difference between the measurement in this region and the one in region 2 , which has the largest abundance, is significant at less than 2$`\sigma `$ for all elements in both MOS and PN with the exception of $`Si`$ for the PN, which is significant at about 2.7$`\sigma `$.

In summary, the metal abundance profiles derived by fitting single temperature models are very different from the ones that are obtained by fitting the more appropriate spectral models II or III. The former all show a dramatic abundance decrease towards the center at radii smaller than one arcmin. The latter show a flattening of the abundance profile within 1 arcmin and a hint of a decrease in the innermost bin. In the light of these results the evidence for absorption in the core of M87 is greatly reduced. Observations at higher angular resolution with the Chandra observatory will likely settle the issue definitively.

###### Acknowledgements.

We thank the many people who have contributed to building, calibrating and operating the EPIC instrument on-board XMM-Newton. S. Ghizzardi is thanked for useful discussions.


text:
# The formation and evolution of binary systems. III. Low-mass binaries in the Praesepe cluster.Based on observations obtained at the Canada-France-Hawaii Telescope

## 1 Introduction

Binary and multiple systems provide a fossil record of the star formation process. In the last decade, studies of various galactic populations (field stars, open clusters, star forming regions) have led to the conclusion that most solar-type stars occur in binary systems rather than in isolation (e.g., Duquennoy & Mayor 1991, Mathieu et al. 2000). Hence, the most common output of protostellar collapse appears to be the formation of multiple systems. Beyond this indisputable observational result, the way stellar systems form remains an important issue (Bodenheimer et al. 2000), and so a robust determination of the detailed properties of binary stars, e.g., the distributions in their orbital periods, mass-ratios, and orbital eccentricities, is critically needed to guide theoretical models (Clarke et al. 2001, Ghez 2001).

Equally important in order to get clues to the star formation process is to determine whether the properties of binaries depend upon the environment in which they form, and whether these properties evolve over time or, to the contrary, remain stable during pre-main sequence and main sequence evolution: that is, are the statistical properties of binary populations universal and do they unambiguously reflect the processes which gave them birth, or do they vary both over time and from place to place in the solar neighbourhood? In order to address these issues, multiple systems have to be sampled and characterized in various types of environments and in stellar populations that have reached different stages of evolution.

A number of studies have been devoted to these issues. Large scale searches for binaries have been completed among low-mass field dwarfs (Duquennoy & Mayor 1991, Fisher & Marcy 1992, Tokovinin 1992), T Tauri stars in young stellar associations (e.g., Leinert et al. 1993, Ghez et al. 1993, Simon et al. 1995) and low-mass stars in young clusters (e.g., Bouvier et al. 1997, Duchêne et al. 1998, Patience et al. 1998, Mermilliod & Mayor 1999). The results indicate that in most surveyed regions the low-mass stars exhibit a similar fraction of close visual binaries with semi-major axes between a few tens to about 1000 AU (see Duchêne 1999 for a summary). There are, however, two notable exceptions. First, the Taurus star-forming association appears to contain about twice as many binaries as are present in the field (Leinert et al. 1993, Ghez et al. 1993). Second, there is a marked deficit of wide binaries, with semi-major axes in the range from 1000 to 5000 AU, in the Orion cluster as compared to the field binary population (Scally et al. 1999).

One of the difficulties in interpreting these results is that the various studies were performed with different techniques, thus introducing different observing biases (see Duchêne 1999). Motivated by the need to obtain homogeneous data sets for close visual binaries, especially in cluster environments whose study has been somewhat neglected compared to star forming regions and field populations, we have started high angular resolution surveys of low-mass stars in nearby young clusters. For this project, we have selected a number of clusters having ages in the range from 2 Myr to 700 Myr, in order to investigate any evolution of the binary content with time. The first 2 papers of this series reported on the results obtained for the $``$100 Myr Pleiades cluster (Bouvier et al. 1997, Paper I), the $``$2 Myr IC 348 cluster (Duchêne et al. 1999, Paper II), while preliminary results for the $``$80 Myr Alpha Per cluster have been reported by Eislöffel et al. (2001).

We report here the results obtained on the binary population of the 700 Myr old Praesepe cluster from high angular resolution observations of 149 low mass cluster members. Section 2 describes our data aquisition and analysis techniques, which are similar to those used in Papers I and II. The binary frequency among the low-mass stars in Praesepe and the binary properties of the cluster are derived in Section 3. Combining these data with those available in the literature for other young clusters, star forming associations, and the field, we discuss in Section 4 the implications of the trends (or lack of) observed in binary frequency as a function of time and environmental conditions for the binary formation process.

## 2 Observations and data analysis

The sample was primarily drawn from the proper motion study of Praesepe candidate members by Jones & Stauffer (1991, hereafter JS91). From this catalogue, we first selected 70 candidates in a B-V range between 0.52 and 1.4 which had a proper motion membership probability larger than 90%. These candidates were observed in February and December, 1997. For another observing session in January 1998, 53 additional candidates were selected from JS91 in the same B-V range, which had a proper motion membership probability larger than 75% and visual photometry consistent with membership in a (V,B-V) color-magnitude diagram. We completed our sample with 26 additional candidates having both radial velocities and photometry consistent with Praesepe membership, which were selected from Mermilliod’s Open Clusters Database (WEBDA, Mermilliod 1999) and amongst Praesepe halo candidates (Mermilliod et al., 1990).

A total of 149 Praesepe G and K dwarfs were thus observed in 1997 and 1998 at the Canada-France-Hawaii Telescope with the adaptive optics system PUEO (Rigaut et al. 1998). The IR camera was MONICA (Nadeau et al. 1994) in Feb. 1997 and KIR (Doyon et al. 1998) in Dec. 1997 and Jan. 1998, which provide a field of view of 9″ and 36″, respectively. We used the same aquisition procedure for all runs. The targets were observed successively in four quadrants of the camera, for a total integration time of typically 60s; the 4 exposures were subsequently registered and added to produce a final image that was properly sky subtracted and flat-field corrected. The primaries were first observed in either the H or K band, depending on the atmospheric turbulence conditions, in order to optimize the adaptive correction. Whenever a binary was detected in real-time, it was observed in at least one other filter in order to subsequently check the membership of the companion in a color-magnitude diagram. During each run, the images produced by the adaptive optics system were diffraction-limited in the H and K bands, providing a spatial resolution of 0.09″ and 0.13″ FWHM, respectively.

Aperture photometry was performed using IRAF/APPHOT and calibrated using several UKIRT Faint Standards observed during the 3 runs. The photometric accuracy is of order of 0.05 magnitudes in the JHK bands. Table 1 lists the near infrared magnitudes of non resolved primaries while those of Praesepe binaries are listed in Table 2. For the binaries, differential photometry was obtained by fitting a template PSF simultaneously to the primary and the companion within IRAF/DAOPHOT. The PSF was provided by unresolved Praesepe stars that were observed just before and/or just after the exposure on the binary. Several PSF templates were used for each binary, thus providing an estimate of the photometric error. Typically, the differential photometry is accurate to within 0.02 magnitudes, but the error may be larger for binaries with a separation close to the resolution limit or for companions close to the detection limit. The magnitude difference between the companion and the primary was combined with the aperture photometry of the system to provide the magnitude of each component.

The separation and position angle were derived from the photocenter of the components in the image as provided by the PSF fitting algorithm. The plate scale and orientation of the detector were calibrated by observing IDS astrometric standards (Van Dessel & Sinachopoulos, 1993) during each run. The rms error is typically 5 mas on the separation and 0.1° on the position angle.

## 3 Results

We detected 26 binaries and one triple system having separations less than 7″ among 149 Praesepe G and K dwarfs. Table 2 provides their names and cross-identifications, visual photometry, and indicates whether they were previously known as either photometric or spectroscopic binaries.

### 3.1 Comparison with previous work

Of the 26 binaries reported here, eight were previously known as spectroscopic binary (SB) systems (Mermilliod & Mayor 1999): three long-period SBs (BDA 297, 322, 533) are resolved here with a projected separation of 23 AU or less, while three others have too short a period to be resolved (BDA 287, 540, 1184) and so the companion reported here makes them triple systems. The two remaining systems, BDA 365 and BDA 495, are known to be spectroscopic triples (Mermilliod et al. 1994). They both consist of a long-period binary with one of the components itself being a short-period binary. Both objects are spatially resolved here, and, given their angular separations and spectroscopic orbits, it is likely that we resolved the widest binary pair of these systems.

Figure 1 shows a (V, B-V) color-magnitude diagram (CMD) for the whole sample of G and K Praesepe members observed in this survey. The different symbols indicate previously known spectroscopic binaries, visual binaries detected here by adaptive optics, and “single” stars. Not surprisingly, the sample contains a number of (presumably short-period) SBs located on the main and binary sequences of the cluster that are not resolved by adaptive optics. More interestingly, about 15 objects that are not detected as binaries, either from spectroscopy or adaptive optics, lie more than 0.5 mag above the main sequence of the cluster and must therefore be nearly equal-mass binary systems. In order to have escaped detection, these systems must have an orbital period in the range from about 3 103 to 3 104 days, i.e., a semi-major axis in the range from about 4 to 15 AU. Finally, although many of the binaries that are resolved in the present study are displaced above the main sequence as expected, about a third of them are located on the cluster main sequence and were not previously detected as binaries through photometry. BDA 1184 (the triangle at V = 11.85, B-V = 0.79) is the most extreme case: it appears to be a single star but is in fact a triple system, a short-period spectroscopic binary and a wider companion at 1$`\stackrel{}{.}`$28. This simply means that many faint and red companions are hidden in visible light and are too light-weight to show up in radial-velocity observations.

A (K, V-K) CMD, also shown in Figure 1, is more appropriate to detect very red companions photometrically. In this diagram, BDA 1184 shows a vertical displacement of 0.42 mag. Another example, BDA 2418, which lies right on the single star sequence in the (V, B-V) plane, already shows a displacement of 0.15 mag in (I, V-I) and 0.35 mag in (K, V-K). However, the star BDA 287, another triple system, is still located on the single star sequence even in (K, V-K), because the K magnitude difference is 5.13.

To conclude this section, the common use of (V, B-V) colour-magnitude diagram to study binarity in open clusters from photometric data is not a good choice or strategy. (I, V-I) planes are certainly better and JHK observations provide still much more information on faint red companions. Observations through the $`BVIK`$ filters would therefore permit the detection of a wealth of new binary candidates, which would raise the binary frequency to values that are probably more realistic. However, direct observations are still needed to detect systems with large magnitude differences or with a multiplicity of order higher than 2.

### 3.2 Binary frequency and orbital period distribution

Table 2 also lists the astrometric and infrared properties of the systems and their components. Although BDA 322 is clearly elongated on the images, its separation is too small for deriving precise astrometric and photometric properties, and so we list only approximate values in Table 2. The masses of the individual binary components and the resulting mass ratios, $`q=M_2/M_1`$, were derived from the JHK photometry of the components using the mass-magnitude relationships from Baraffe et al. (1998) models, assuming an age of 0.7 Gyr and a distance modulus for Praesepe of $`(mM)=6.28`$ ($`d=180`$pc, Robichon et al. 1999). For the primaries or secondaries that are themselves unresolved spectroscopic binaries, this method would overestimate their mass by up to about 20% for equal-brightness components.

In order to ascertain photometric membership, we plotted the primaries and secondaries in various JHK color-magnitude diagrams, comparing their location in these diagrams with the 700 Myr isochrone from Baraffe et al. (1998). The primaries and secondaries are all consistent with being Praesepe members within the photometric errors. However, the third component of the BDA 401 system, BDA 401C, lies far away from the isochrone and is probably a field object. We therefore ignore it in the following discussion and consider BDA 401 to be a double. Finally, a few binaries were observed at only one wavelength, but these systems are so tight that their binary nature is not in serious doubt.

That we find only one chance projection at a distance of less than 7″ in a sample of 149 stars is consistent with the 2MASS Point Source Catalogue Statistics, which predict about 1200 objects per square degree down to a magnitude of K=15 in the direction of Praesepe. This translates into $``$0.015 objects within a 7″ radius and leads to an estimate of 2 chance projections in the present sample.

In order to estimate the binary frequency among the Praesepe G and K dwarfs, a correction factor has to be applied to the number of systems actually detected to account for the detection limit of our survey. The largest magnitude difference we are able to detect between the secondary and the primary is shown in Figure 2 as a function of angular separation. At separations close to the diffraction limit, the detection of faint companions is limited by the constrast against the bright primary while, at large separations, the detection is background limited. We therefore proceed to derive a correction factor in various separation bins as follows.

In each bin, we convert the maximum reachable contrast between the secondary and the primary, $`\mathrm{\Delta }m_{max}`$, into a minimum mass ratio $`q_{min}`$ using the mass-magnitude relationship from the Baraffe et al. (1998) models. We assume that the mass-ratio distribution of the Praesepe binaries is the same as that of the field G dwarfs derived by Duquennoy & Mayor (1991, DM91). This assumption is consistent with the mass ratio distribution derived for our binary sample. While the overall q–distribution is rather flat (see Table 2), restricting the analysis to the separation range where we can detect all companions down to $`q=0.1`$ (sep$``$0.8″, see Table 3) so that we can compare to DM91, we count 10 companions. Mass ratios in this small subsample range from 0.11 to 0.92, with a mean of 0.38. In DM91 survey, for binaries with periods longer than 104 days, the average mass ratio is 0.40. Both results are very similar, though ours admitedly relies on only a few systems. The fraction of missed companions is then found by integrating DM91’s mass-ratio distribution between q=0.1 and q=$`q_{min}`$. We finally apply this correction factor to the number of detected binaries to obtain the total number of systems with mass-ratios larger than 0.1 in this separation range. The 1$`\sigma `$ uncertainties on BF in each $`\mathrm{log}\text{P}`$ bin correspond to Poisson noise, i.e., $`\sigma =\sqrt{N_d}\times (1+N_u/N_d)`$, where $`N_d`$ and $`N_u`$ are the number of detected and missed systems, respectively (see Table 3).

In order to establish the distribution of orbital periods, i.e., the frequency of binary systems in each $`\mathrm{log}\text{P}`$ bin, we convert the intervals of projected separation into bins of orbital periods. Angular separations $`\rho `$ are statistically corrected for projection effects to yield the semi-major axis according to: $`\mathrm{log}a=\mathrm{log}(\rho \times d)+0.1`$ (DM91), where $`d=180`$ pc is the distance to the Praesepe cluster. Kepler’s 3rd law with an average mass of 1.3 $`M_{}`$ for the system yields the corresponding range of orbital periods (see Table 3).

The last lines of Table 3 list and also Figure 3 illustrates the derived frequency of Praesepe low-mass binaries in each $`\mathrm{log}\text{P}`$ bin between 4.4 and 6.9 (P in days), as compared to the frequency of G dwarf binaries over the same separation range111As in our previous paper on the Pleiades (Bouvier et al. 1997), we define the binary frequency, B.F., as the number of binary orbits divided by the number of primaries in the sample. This is equivalent to the companion star fraction, $`csf=(B+2T)/(S+B+T)`$, where S, B, and T are single, binary, and triple systems, respectively.. Even though the statistical uncertainties in each bin are somewhat large, especially at the shortest orbital periods where the incompleteness correction is significant, the overall binary frequency appears to be very similar among the solar-type Praesepe stars and the field G dwarfs, amounting to 25.3$`\pm `$5.5% and 23.8%, respectively, for the 4.4-6.9 range in $`\mathrm{log}\text{P}`$. Restricting the comparison to the range of orbital periods not affected by detection biases ($`\mathrm{log}\text{P}`$ from 5.4 to 6.9), the BF amounts to 10.7$`\pm `$2.6% for Praesepe stars compared to 13.4% for field dwarfs. We thus conclude that the frequency of long period binaries among Praesepe G and K stars is indistinguishable from that measured among the G field dwarfs by DM91. This conclusion seems to apply to short-period systems as well ($`\mathrm{log}\text{P}4.0`$) as Mermilliod’s & Mayor’s (1999) investigation of spectroscopic binaries in the Praesepe cluster yields BF$`=25\pm 5\%`$, as compared to 21% for field dwarfs.

Praesepe has nearly the same age as the Hyades cluster. Unfortunately, comparison of the BF between the two clusters is limited by the fact that there is only a slight overlap in the separation ranges covered by Patience et al.’s (1998) survey in the Hyades and ours in Praesepe. In the 15–50 AU common range of semi-major axes, we detect 8 companions, i.e. an observed companion fraction of 5.3$`\pm `$1.9%. In the same separation range, Patience et al. (1998) detected 9 companions, out of which 3 would not have been detected in our survey given their flux ratios. This amounts to a companion star fraction of 3.7$`\pm `$1.5%, less than 1$`\sigma `$ lower than our Praesepe estimate. Within this restricted range of orbital periods, we thus do not find any significant difference between Hyades and Praesepe binary fraction, though this conclusion is obviously based on small number statistics.

## 4 Discussion

In this section, we first briefly review current observational results on the statistics of binary frequency in young clusters and associations, including results from both this and previous papers of this series and other published work. The comparative analysis of the properties of binary populations in different environments and at various evolutionary stages provides constraints on binary formation models, which are discussed below with the emphasis on the difficulties that are encountered by the various formation scenarios. The discussion is restricted to visual (a $``$ 20-1000 AU), low-mass binaries ($``$ 0.5-1.1$`\text{ M}_{}`$) that are resolved by high angular resolution techniques.

### 4.1 Binary statistics in clusters and associations

The three young clusters we surveyed for low-mass binaries, IC 348 ($``$2 Myr, Duchêne et al. 1999), the Pleiades (120 Myr, Bouvier et al. 1997), and Praesepe (700 Myr), have been observed with the same instrumentation, thus resulting in similar detection biases. Moreover, the results have been analyzed in a consistent way, in particular with regard to the incompleteness corrections. The three clusters are found to exhibit identical binary frequencies over the $`\mathrm{log}\text{P}`$range $``$4.5-7.0, $`<`$BF$`>`$ $``$ 0.25 $`\pm `$ 0.05, which is also consistent with the BF measured for solar-type field dwarfs in the same range of orbital periods (BF $``$ 0.24, DM91). Similar results have been obtained for other clusters, Orion ($``$2 Myr, Prosser et al. 1994, Petr et al. 1998, Simon et al. 1999), Alpha Persei ($``$ 80 Myr, Eislöffel et al. 1999), and the Hyades ($``$600 Myr, Patience et al. 1998) — all of which exhibit a low-mass BF in a restricted $`\mathrm{log}\text{P}`$range that is consistent with the field dwarf BF when the results are analyzed in a uniform way (see Duchêne 1999).

Since the ages of these clusters cover a huge range from 2 to 700 Myr, the results suggest that low-mass binaries are formed at a very early stage of cluster evolution (at an age $``$ 1 Myr) and that the binary fraction does not evolve much thereafter until the cluster eventually dissolves into the field ($``$ 1 Gyr). The lack of a significant evolution of the binary population during the secular dynamical evolution of a cluster is consistent with recent numerical simulations (e.g. Kroupa 2000).

Note, however, that this conclusion may not hold for more massive and/or spectroscopic binaries. Abt & Willmarth (1999) have reported marginal evidence for a BF that increases with a cluster’s age, from the Orion Nebula Cluster to Praesepe, for spectroscopic binaries with A-type primaries. They interpret this trend as the possible signature of binary formation by capture in evolving clusters and/or preferential escape of single stars during the secular dynamical evolution of clusters (de la Fuente Marcos 1997). No such trend is seen for low-mass wide binaries.

Another general result is that pairs of nearly coeval clusters (e.g., IC348 and Orion, Alpha Persei and the Pleiades, Praesepe and the Hyades) not only exhibit similar fractions of low-mass visual binaries over the separation range probed by adaptive optics and speckle techniques ($``$20-1000 AU), but the distribution of orbital periods is consistent as well, with admittedly large uncertainties in the shape of the $`\mathrm{log}\text{P}`$ distribution (see Fig. 3). These similarities suggest that either all these clusters were formed under very similar conditions and have evolved in the same way, or else the binary content of clusters and their properties depend only weakly on initial conditions.

In marked contrast with the results obtained for cluster binaries, the binary frequency in low-density star forming associations, most notably the Taurus-Auriga cloud, is higher by a factor of about 2 than that observed in both clusters and field solar-type stars (Leinert et al. 1993, Ghez et al. 1993). Ghez (2001) also reported a significant difference between the $`\mathrm{log}\text{P}`$distributions of binaries in clusters and those in associations, the latter harbouring a larger fraction of wider binaries than the former.

The different properties of binary populations in clusters and associations can be a signature of different formation mechanisms in these environments (e.g., Durisen & Sterzik 1995). Alternatively, if one assumes a universal formation mechanism that yields the same initial BF in clusters and associations, the observed differences could reflect dynamical processes acting very early-on, such as the rapid disruption of wide primordial binaries in clusters through gravitational encounters (e.g. Kroupa et al. 1995). We discuss these two possibilities in turn below.

### 4.2 A universal mechanism for binary formation?

Among the various possible ways of forming low-mass binaries, tidal capture has been shown to be inefficient even in the densest protostellar clusters (Clarke & Pringle 1991, Kroupa 1995, Clarke 2001) and fission of massive protostars or protostellar disks, which could conceivably yield the tightest binaries, seems to be prevented by the development of bar-like instabilities (Durisen et al. 1986, Bate 1998). Therefore, multiple fragmentation during protostellar collapse appears today to be the most promising mechanism for creating wide multiple systems (Bodenheimer 2001).

Recent collapse calculations indicate that the likely output of multiple fragmentation is the formation of small-N protostellar aggregates, where N$``$ 3-10 (e.g., Burkert, Bate & Bodenheimer 1997, Klessen & Burkert 2000). These aggregates experience rapid dynamical decay and eventually leave a bound binary system, while other fragments are dynamically ejected mostly as single remnants (e.g., McDonald & Clarke 1993, 1995; Sterzik & Durisen 1998). Since few-body interactions occur on a small scale within protostellar aggregates (r $``$ a few 100 AU) and on a very short time scale ($`10^4`$ yr), the resulting primordial binary fraction is not expected to depend strongly upon the global properties of the star forming region.

The fraction of primordial binaries that results from the dynamical decay of protostellar aggregates is usually identified with the high BF observed in loose associations like Taurus. Then, the lower BF measured in clusters is thought to result from the rapid disruption of primordial binaries, through destructive gravitational encounters that occur on a time scale of less than 1 Myr (e.g., Kroupa, Petr, & MacCaughrean 1999). Since the rate of gravitational encounters scales with the local stellar density, this scenario conceivably accounts for the observed trend of lower binary fractions in denser star forming regions (Patience & Duchêne 2001), and it is further supported by the paucity of wide binaries observed in the ONC (Scally et al. 1999) and, more generally, in young open clusters (Ghez 2001).

This mode of binary formation is not exempt from difficulties, however. One issue is whether the Taurus binaries can be regarded as representative of a universal population of primordial binaries. The frequency of primordial binaries that is expected from the decay of small-N aggregates is of order of BFp $``$ 1 / (N-1), i.e., at most 50% for N=3. This is significantly lower than the BF measured in the Taurus association, which amounts to $``$80% for stars in the mass range $`0.31.2\text{ M}_{}`$, with little dependence on the primary mass (Leinert et al. 1993). This discrepancy could be solved if single fragments that were dynamically ejected from protostellar aggregates had escaped from their birth place. With typical ejection velocities of 3-4 $`\text{ km}\text{ s}^1`$(Sterzik & Durisen 1995), they would be located a few parsecs away from their birth site at an age of 2 Myr, i.e., a few degrees away from the Taurus stellar groups (Gomez et al. 1995). A widely distributed population of X-ray emitting T Tauri stars has been detected with ROSAT over the Taurus cloud (Wichmann et al. 1996, 2000, Frink et al. 1997), but these stars do not seem to be preferentially single (Köhler & Leinert 1998), as would be expected if they were escapers.

An intriguing possibility is that the ejected fragments are very low mass, indeed substellar, objects (Sterzik & Durisen 1999, Clarke & Reipurth 2001) that might so far have escaped detection in the Taurus cloud. Although the search for brown dwarfs in Taurus has been somewhat disappointing (Luhman 2000), it has only concentrated on very limited areas centered on the small Taurus stellar groups. A widely distributed population of (single) substellar objects over the Taurus cloud could reconcile the high BF frequency measured for Taurus stars with the lower BF expected from the decay of small-N aggregates, which includes both stellar and substellar fragments. In support of this hypothesis, we note that current determinations of the substellar IMF do indicate that isolated brown dwarfs are numerous in clusters (e.g. Luhman et al. 2000, Moraux et al. 2001) and appear to be preferentially single objects (Martín et al. 1999). If originally ejected from small-N aggregates, substellar fragments may be more easily retained in the deep potential well of dense clusters than in loose associations (de la Fuente Marcos & de la Fuente Marcos 2000), which might explain why they have not been found in the central regions of Taurus.

Another aspect of the models that is challenged by the observations is whether the lower BF of clusters compared to associations can be understood as the mere result of the disruption of primordial binaries. Models that describe the dynamical evolution of primordial binaries in young clusters, starting from an initial distribution similar to the one observed in Taurus, show that the destruction rate of wide primordial binaries ($`\mathrm{log}\text{P}`$5-7) is a sensitive function of the initial stellar density (e.g., Kroupa 1995; Kroupa, Aarseth & Hurley 2000). Yet, all clusters studied so far appear to harbour the same BF to within a few percent in this $`\mathrm{log}\text{P}`$range. The lack of dispersion in the BF measured for clusters is then surprising, given that it is unlikely all clusters surveyed so far have formed with precisely similar densities. For instance, the stellar density in the Trapezium cluster is of order of 5 104 pc-3 (McCaughrean & Stauffer 1994) whereas, at a similar age, it is about 5 103 pc-3 in IC 348 (Herbig 1998). If gravitational encounters leading to the disruption of primordial binaries are the dominant mechanism that yields a lower BF in clusters, one would expect to observe somewhat different binary fractions between clusters themselves.

Hence, while scenarios of binary formation and evolution that assume an initially large fraction of primordial binaries in all star forming regions, followed by a rapid erosion of the binary population in dense clusters, have recently become quite popular, it remains to be seen whether the difficulties outlined above can be solved.

### 4.3 Do local conditions impact on the binary formation process?

As an alternative to a universal formation mechanism, it is probably too early yet to rule out binary and, indeed, single star formation as the direct outcome of cloud collapse and fragmentation, without going through the transient episode of small-N protostellar aggregates. Unfortunately, the theory and simulations of fragmentation are not yet predictive enough, and the final product of protostellar collapse can depend sensitively on initial conditions (see Bodenheimer et al. 2000 for a review), e.g., the radial density profile of the parental cloud (Burkert et al. 1997), its temperature (Sterzik and Durisen 1995), turbulence (Klein 2001), the magnetic field (Boss 2001), etc. The large scale environment, such as cloud-cloud collisions, or other external impulsive processes, such as supernova blasts, may also impact on the fragmentation process (Whitworth 2001). Hence, one might expect that different initial conditions in star forming regions lead to significant variations in the properties of the young stellar populations they harbour.

Since the initial conditions that led to star formation in a given molecular cloud are usually poorly known, it is somewhat difficult to constrain this alternative mode of binary and single stars formation with current observations. As noted above, however, one of the striking results of the recent binary surveys is the quasi-universality of the BF in clusters, which all appear to harbour the same fraction of solar-type, wide (sep $``$ 20 AU) binaries to within the statistical uncertainties. This suggests that the fragmentation process, if directly responsible for the formation of binary systems, might not be as sensitive to local conditions as numerical simulations tend to indicate. On the other hand, while the results for cluster binaries are homogeneous and similar to those obtained for the field binary population, the much larger BF observed in the Taurus cloud would seem to suggest that gross variations in the local conditions do impact on the fragmentation process.

In this respect, it is interesting to note that binary frequency is not the only difference that exists between the stellar populations of the Orion cluster and Taurus association. Significant differences have also been found in the distribution of stellar angular momentum among their low-mass T Tauri stars (Clarke & Bouvier 2000), and in the distribution of their stellar masses, with Taurus harbouring apparently both fewer high-mass stars and fewer very-low mass objects than Orion (Hartmann& Kenyon 1995, Luhman 2000). It is tempting to think that the differences observed in the fundamental properties (mass and angular momentum distributions, binary frequency) of the Orion and Taurus populations are causally related and point to a common origin that reflects intrinsically different modes of star formation in clusters and in associations (e.g., Myers 1998, Williams et al. 2000, Motte & André 2001).

## 5 Conclusion

From an adaptive optics imaging survey of 149 G and K-type primaries of the Praesepe cluster, we find that solar-type cluster members harbour the same proportion of close visual binaries as do G-type field dwarfs. Long lived open clusters, such as Praesepe, probably started their evolution as extremely dense protostellar clusters. Yet, only about 10% of the field population is thought to result from the dissipation of such rich clusters. At the other extreme, Taurus-like regions of distributed star formation have very low star forming efficiencies. Hence, as recently advocated by Adams & Myers (2001), most field stars must have been born in stellar groups which dissipate in a few million years, corresponding to initial conditions somewhat intermediate between dense protostellar clusters destined to become young open clusters and loose associations. The very similar binary fraction measured for solar-type stars in young open clusters (Praesepe, Pleiades, Alpha Per) and in the field thus suggest that the formation and evolution of low-mass binaries is not very sensitive to local conditions.

The main limitation of studies like the present one which aim at constraining the star formation process through the investigation of young binaries is that they have been mostly concerned with low-mass systems so far and somewhat neglected higher mass binaries (see, however, Preibisch et al. 1999, Garcìa & Mermilliod 2001, Duchêne et al. 2001). If multiple fragmentation of collapsing clouds is the dominant mode for the formation of multiple stellar systems, the mass distribution of fragments in small protostellar groups may largely determine the resulting binary frequency for a given primary mass. For instance, we argued above that the high BF observed for T Tauri stars in Taurus might merely be the result of neglecting a putative population of single brown dwarfs distributed over the cloud. In regions where high mass stars are formed, such as in Orion, more single ejected fragments would be of solar mass or so, thus resulting in a lower binary fraction among low mass stars. The investigation of such causal relationships between the fundamental properties of young stars, e.g., between binary fraction and the mass function, requires the consideration of the whole stellar population of the star forming region with a complete census of multiple systems at all primary masses, which is not available today for any star forming region.

A promising new way to better understand the formation of multiple systems is to investigate, in different environments, extremely young stellar objects still embedded in their natal cloud at the end of the protostellar collapse. The high degree of multiplicity of such Class 0 and Class I “protostellar” sources starts to be revealed from high angular resolution studies in the millimeter range (Looney, Mundy & Welch 2000). The advent of adaptive optics systems equipped with near-IR wavefront sensors on large telescopes will now open the way to large scale surveys of embedded protobinaries with a tenfold increase in angular resolution compared to current millimeter studies, reaching separations as small as a few astronomical units. Such studies will provide unprecedented details on the fragmentation process at a very early stage of evolution of young systems, before any significant dynamical evolution of protostellar systems has occurred.

###### Acknowledgements.

We thank Jean-Luc Beuzit and Olivier Lai for reobserving in Nov.-Dec. 1999 some suspected binaries with the same instrumentation and Isabelle Baraffe for computing and providing a 700 Myr isochrone from her models of low mass stars. We acknowledge useful discussions with Cathie Clarke and Pavel Kroupa on cluster dynamics and with Frédérique Motte on prestellar cores.


text:
# Characterization of harmonic and subharmonic functions via mean–value properties

## 1. Introduction and results

Let $`h`$ be a harmonic function in the closed ball centred at a point $`a`$ and with radius $`r`$, $`\overline{B}(a,r)`$. Then it is well–known that $`h`$ satisfies the mean–value properties

(1)
$$h(a)=_a(r):=\frac{1}{|B(a,r)|}_{B(a,r)}h(x)𝑑x$$

and

(2)
$$h(a)=𝒮_a(r):=\frac{1}{|B(a,r)|}_{B(a,r)}h(x)𝑑s.$$

In both cases, there exist converse results stating that if the equality is satisfied, then $`h`$ must be harmonic. In the first case, for instance, if we assume that $`h`$ is locally integrable on a domain $`\mathrm{\Omega }`$ and satisfies (1) whenever $`\overline{B}(a,r)`$ is contained in $`\mathrm{\Omega }`$, then we have that $`h`$ is harmonic in $`\mathrm{\Omega }`$. A similar result holds for (2) provided we assume now that $`h`$ is continuous in $`\mathrm{\Omega }`$ — for these classical results, see \[ABR, GT\], for instance.

A second, lesser known, set of converse properties to these mean value results assumes that the equalities are satisfied in a set $`\mathrm{\Omega }`$ for all integrable harmonic functions defined in $`\mathrm{\Omega }`$. Then, it is possible to show that this implies that $`\mathrm{\Omega }`$ is a ball \[AG, Be, ES, Ko, Ku, PS\]. In particular, \[Be\]\[Ko\] and \[PS\] consider the case where the surface and volume averages are equal.

Motivated by this, and since we could not find any reference in the literature to a converse result of the first type when it is assumed that the two averages coincide, in this paper we consider several relations between the integral means of harmonic and subharmonic functions. We begin by considering the question of what can be said when $`h`$ satisfies

(3)
$$_a(r)=𝒮_a(r)$$

for balls in a set $`\mathrm{\Omega }`$. We prove that, under suitable assumptions, this last equality also implies that $`h`$ is harmonic. More precisely, we have the following

###### Theorem 1.1.

Let $`h`$ be a continuous function in $`\mathrm{\Omega }`$. Then $`h`$ is harmonic in $`\mathrm{\Omega }`$ if and only if (3) is satisfied for all balls whose closure is contained in $`\mathrm{\Omega }`$.

###### Remark 1.

Note that in the equality (3) it is not assumed a priori that the value of the averages remains fixed when $`a`$ is a given point but the radius changes.

###### Remark 2.

Continuity is essential, since otherwise any function coinciding with a harmonic function on all but a set of finite points, for instance, would clearly satisfy (3) but not be harmonic.

From the results in \[Be\]\[Ko\] and \[PS\] we have that the equality of the averages (3) for all harmonic functions is possible only for balls. M. Rao has shown that in fact the averages of harmonic functions on general domains must satisfy a one–sided inequality of Harnack type \[R\]. Here we improve on this result establishing also a lower bound by using the methods of \[Be\] and \[PS\] to obtain the following

###### Theorem 1.2.

Let $`\mathrm{\Omega }`$ be a domain satisfying the uniform exterior sphere condition. Then there exist constants $`0<c_11c_2<\mathrm{}`$, depending only on $`\mathrm{\Omega }`$, such that

$$\frac{c_1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h𝑑s\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h𝑑x\frac{c_2}{|\mathrm{\Omega }|}_\mathrm{\Omega }h𝑑s$$

for all non–negative harmonic functions $`h`$ defined in $`\mathrm{\Omega }`$ and which take continuous values on the boundary. If either of the constants $`c_1`$ or $`c_2`$ can be taken to be equal to one, then $`\mathrm{\Omega }`$ is a ball.

The proof of the second part of the result also uses a variation of a well–known result of Serrin’s \[S\], which we believe to be interesting in its own right – see Lemma 2.1 for the details.

Another natural question to ask is whether a result similar to Theorem 1.1 holds for subharmonic functions. Here we call a real–valued function $`u`$ defined in $`\mathrm{\Omega }`$ subharmonic, if $`u`$ is continuous and for every ball $`B`$ whose closure is contained in $`\mathrm{\Omega }`$ and every harmonic function $`h`$ in $`B`$ satisfying $`hu`$ on $`B`$ we have that $`hu`$ in $`B`$ \[GT\]. We have that subharmonic functions satisfy

(4)
$$_a(r)𝒮_a(r),$$

so it makes sense to ask if this inequality characterizes subharmonic functions. In this respect we have the following

###### Theorem 1.3.

Let $`u`$ be a continuous function in $`\mathrm{\Omega }`$. Then $`u`$ is subharmonic in $`\mathrm{\Omega }`$ if and only if (4) holds for all balls whose closure is contained in $`\mathrm{\Omega }`$.

From the proof of Theorem 1.2 it follows that the second inequality in that result also holds for $`C^2`$ subharmonic functions and this can easily be improved to $`C^0`$ sub-harmonic functions. Also for particular domains or in dimension $`1`$ one easily exhibits sequences of positive subharmonic functions for which the volume average is bounded and the surface average becomes unbounded. This can be improved for relatively general domains provided the geometry of the boundary allows a sufficiently rich supply of harmonic subfunctions, e.g., if the domain satisfies the *uniform exterior sphere condition*.

###### Theorem 1.4.

Let $`\mathrm{\Omega }`$ be as above, and $`u`$ be a non–negative subharmonic function in $`\mathrm{\Omega }`$ and continuous in $`\overline{\mathrm{\Omega }}`$. Then there exists a constant $`c_2`$ as in Theorem 1.2 such that

$$\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }u𝑑x\frac{c_2}{|\mathrm{\Omega }|}_\mathrm{\Omega }u𝑑s.$$

On the other hand, there exists a sequence of positive subharmonic functions $`u`$ for which the volume average in $`\mathrm{\Omega }`$ remains bounded, while the surface average is unbounded.

Finally, we consider a property of subharmonic functions proved by Beardon \[B\], which states that if $`u`$ is subharmonic in a set $`\mathrm{\Omega }`$, then

(5)
$$𝒮_a(\kappa r)_a(r)$$

for all balls whose closure is contained in $`\mathrm{\Omega }`$ and where

(6)
$$\kappa =\{\begin{array}{cc}\frac{1}{2},\hfill & n=1\hfill \\ e^{1/2},\hfill & n=2\hfill \\ \left(\frac{2}{n}\right)^{1/(n2)},\hfill & n3.\hfill \end{array}$$

We prove that this property again characterizes subharmonic functions. More precisely, we have

###### Theorem 1.5.

Let $`u`$ be a continuous function defined on a domain $`\mathrm{\Omega }`$ and assume that condition (5) is satisfied for all balls whose closure is contained in $`\mathrm{\Omega }`$ and for some $`\kappa `$ satisfying

$$0<\kappa \kappa _1:=\frac{n}{4}+\frac{1}{2}\sqrt{\frac{n^2}{4}+2n}.$$

Then $`u`$ is subharmonic in $`\mathrm{\Omega }`$.

Since the upper bound for $`\kappa `$ in the theorem is larger than (or equal to, in the case where $`n`$ is one) the values of $`\kappa `$ given in (6), this yields a converse of Beardon’s property.

## 2. Harmonic functions

### 2.1. Proof of Theorem 1.1

Assume that (3) holds for a given fixed point $`a`$ in $`\mathrm{\Omega }`$, and all $`r`$ in $`(0,R)`$, for some positive number $`R`$ such that $`\overline{B}(a,R)`$ is contained in $`\mathrm{\Omega }`$. Denoting by $`\omega _n`$ the volume of the unit ball in $`^n`$ we have that

$$\frac{1}{\omega _nr^n}_{B(a,r)}h(x)𝑑x=\frac{1}{n\omega _nr^{n1}}_{B(a,r)}h(x)𝑑s$$

and thus

$$n_0^r\left[_{B(a,t)}h(x)𝑑s\right]𝑑t=r_{B(a,r)}h(x)𝑑s.$$

Letting

$$\phi (r)=_{B(a,r)}h(x)𝑑s,$$

we have that

$$n_0^r\phi (t)𝑑t=r\phi (r).$$

Since $`h`$ is continuous, we can differentiate this with respect to $`r`$ which gives that $`\phi `$ satisfies the differential equation

$$r\phi ^{}(r)+(1n)\phi (r)=0.$$

Hence $`\phi (r)=cr^{n1}`$ for some constant $`c`$ and

$$\frac{c}{n\omega _n}=\frac{\varphi (r)}{n\omega _nr^{n1}}=𝒮_a(r),$$

from which it follows that, for $`a`$ and $`r`$ as above, the average $`_a(r)`$ is independent of $`r`$. Since $`h`$ is continuous, the value taken by this average must be attained by $`h`$ for some point in the ball $`B(a,r)`$. As this holds for arbitrarilly small $`r`$, we have that this value must coincide with $`h(a)`$. Applying now the converse result to the mean–value equality (1) completes the proof of Theorem 1.1.

### 2.2. The integral Harnack inequality

Let $`v`$ be the solution of the equation

(7)
$$\{\begin{array}{cc}\mathrm{\Delta }v+1=0,\hfill & x\mathrm{\Omega }\hfill \\ v=0,\hfill & x\mathrm{\Omega }.\hfill \end{array}$$

We have that

$$\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h𝑑x=\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h\mathrm{\Delta }v𝑑x=\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h\left(\frac{v}{\nu }\right)𝑑s.$$

Since there exist constants $`c_1`$ and $`c_2`$ such that

$$0<c_1<\frac{v}{\nu }<c_2<\mathrm{},$$

the inequalities follow.

If $`c_1`$ can be taken to be equal to $`1`$, we have that

(8)
$$_\mathrm{\Omega }h\left(\frac{1}{|\mathrm{\Omega }|}+\frac{1}{|\mathrm{\Omega }|}\frac{v}{\nu }\right)𝑑s0.$$

Let $`c`$ be a positive real number such that

$$h_0=c+\frac{1}{|\mathrm{\Omega }|}+\frac{1}{|\mathrm{\Omega }|}\frac{v}{\nu }$$

is positive, and choose now $`h`$ to be equal to $`h_0`$ on the boundary. Since

$$_\mathrm{\Omega }\frac{1}{|\mathrm{\Omega }|}+\frac{1}{|\mathrm{\Omega }|}\frac{v}{\nu }ds=0,$$

replacing $`h`$ in (8) gives

$$_\mathrm{\Omega }\left(\frac{1}{|\mathrm{\Omega }|}+\frac{1}{|\mathrm{\Omega }|}\frac{v}{\nu }\right)^2𝑑s0,$$

from which it follows that $`v/\nu `$ is constant on the boundary. By Serrin’s result \[S\], $`\mathrm{\Omega }`$ must be a ball.

If $`c_2`$ equals $`1`$, we now have

$$\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h\left(\frac{v}{\nu }\right)𝑑s\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h𝑑s$$

and taking $`h=v/\nu `$ on the boundary yields

$$\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }\left(\frac{v}{\nu }\right)^2𝑑s\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }\left(\frac{v}{\nu }\right)𝑑s=\frac{|\mathrm{\Omega }|}{|\mathrm{\Omega }|}.$$

The result now follows from the following lemma, which is a variation on Serrin’s result \[S\].

###### Lemma 2.1.

Let $`v`$ be the solution of equation (7). Then

$$\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }\left(\frac{v}{\nu }\right)^2𝑑S\left(\frac{|\mathrm{\Omega }|}{|\mathrm{\Omega }|}\right)^2,$$

with equality if and only if $`\mathrm{\Omega }`$ is a ball.

###### Proof.

By the Cauchy–Schwarz inequality we have that

(9)
$$_\mathrm{\Omega }\left(\frac{v}{\nu }\right)𝑑S|\mathrm{\Omega }|^{1/2}\left[_\mathrm{\Omega }\left(\frac{v}{\nu }\right)^2𝑑S\right]^{1/2}.$$

Thus

$$\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }\left(\frac{v}{\nu }\right)^2𝑑S\left(\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }\frac{v}{\nu }𝑑S\right)^2=\left(\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }\mathrm{\Delta }v𝑑x\right)^2=\left(\frac{|\mathrm{\Omega }|}{|\mathrm{\Omega }|}\right)^2.$$

On the other hand, equality holds in (9) if and only if $`v/\nu `$ is constant on the boundary and thus, by Serrin’s result, $`\mathrm{\Omega }`$ must be a ball. ∎

## 3. Subharmonic functions

### 3.1. Proof of Theorems 1.3 and 1.4

We begin by proving the following

###### Lemma 3.1.

If $`u`$ is a continuous function in $`\mathrm{\Omega }`$ such that $`_a(r)𝒮_a(r)`$ holds for all balls whose closure is contained in $`\mathrm{\Omega }`$, then $`u`$ satisfies the maximum principle in $`\mathrm{\Omega }`$.

###### Proof.

Let $`M`$ be the maximum of $`u`$ in $`\overline{\mathrm{\Omega }}`$, and assume that this is attained at an interior point. Define

$$𝒜=\{x\mathrm{\Omega }:u(x)=M\}$$

From the continuity of $`u`$ we have that $`𝒜`$ is closed in $`\mathrm{\Omega }`$ and that $`𝒮`$ is a continuous function of the radius. Let $`x_0`$ be a point in $`𝒜`$, and take $`r_1`$ such that $`\overline{B(x_0,r_1)}`$ is in $`\mathrm{\Omega }`$. We have that either there exists $`\overline{r}`$ in $`(0,r_1)`$ such that $`𝒮_{x_0}(\overline{r})<M`$ or not. In the second case, $`𝒮_{x_0}(r)=M`$ for all $`r`$ in $`(0,r_1)`$, which immediately yields that $`x_0`$ is an interior point of $`𝒜`$.

In the first case, let

$$r_{}=inf\{r(0,\overline{r}):𝒮_{x_0}(r)=S_{x_0}(\overline{r})\}.$$

Again by continuity, we have that $`r_{}`$ is positive and also that $`𝒮_{x_0}(r_{})=𝒮_{x_0}(\overline{r})<M`$. Hence

$$\begin{array}{ccc}\omega _nr_{}^n_{x_0}(r_{})\hfill & =\hfill & _{B(x_0,r_{})}v𝑑x\hfill \\ & =\hfill & n\omega _n_0^r_{}r^{n1}𝒮_{x_0}(r)𝑑r\hfill \\ & >\hfill & \omega _nr_{}^n𝒮_{x_0}(r_{})\hfill \end{array}$$

and thus $`_{x_0}(r_{})>𝒮_{x_0}(r_{})`$, contradicting the hypothesis. ∎

Proof of Theorem 1.3: Assume that $`w`$ satisfies the inequality $`_a(r)𝒮_a(r)`$ but that it is not subharmonic. Then, there exists a ball $`B`$ and a function $`h`$ harmonic in $`\overline{B}`$ such that $`\overline{B}`$ is contained in $`\mathrm{\Omega }`$ and $`hw`$ on $`B`$, but $`h(x_0)<w(x_0)`$ for some $`x_0`$ in $`B`$. Since the function $`wh`$ also satisfies the inequality but not the maximum principle, we have a contradiction and hence $`w`$ must be subharmonic. ∎

Proof of Theorem 1.4: Let $`v`$ be subharmonic. Let $`h`$ be the harmonic solution of the Dirichlet problem with boundary values specified by $`v`$. Then $`vh`$ is subharmonic and the maximum principle shows that $`vh`$. Hence Theorem 1.2 applied to $`h`$ implies

$$\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }v𝑑x\frac{1}{|\mathrm{\Omega }|}_\mathrm{\Omega }h𝑑x\frac{c_2}{|\mathrm{\Omega }|}_\mathrm{\Omega }h𝑑s=\frac{c_2}{|\mathrm{\Omega }|}_\mathrm{\Omega }v𝑑s.$$

For the second part of the statement consider first the trivial one dimensional case. Assume (without loss of generality) $`\mathrm{\Omega }=(0,1)`$ and let, for $`k`$ a positive integer, $`v_k(x)=\mathrm{max}\{1,kk^2x,k+k^2(x1)\}`$. Then the surface average is $`k`$ and the volume average is bounded.

To generalize the one dimensional case we will need for each element of the sequence a maximum of a finite number of harmonic functions which is large in a neighborhood of $`\mathrm{\Omega }`$ and small elsewhere. To deal with a fairly general geometry of $`\mathrm{\Omega }`$ requires using in the construction harmonic functions other than affine functions and will require a covering argument to control the oscillation of these on $`\mathrm{\Omega }`$. For domains satisfying the uniform exterior sphere condition translating and rescaling fundamental solutions will be enough. Details follow below.

Let $`k>2`$ be an integer and $`n>1`$. For each $`y\mathrm{\Omega }`$, let $`\overline{B}_{\delta _k}(\overline{y}^k)`$ denote a closed ball such that $`\overline{B}_{\delta _k}(\overline{y}^k)\overline{\mathrm{\Omega }}=\{y\}`$ and $`\delta _k>0`$. Here $`(\delta _k)`$ is a sequence of positive numbers smaller than the radius in the uniform exterior sphere condition and verifying $`\delta _k0`$.

Denote the fundamental solution of the Laplace operator by $`\mathrm{\Psi }(|x|)`$. Define, for each $`k>2`$ and $`y\mathrm{\Omega }`$

$$v_y^k(x)=k^2\mathrm{\Psi }(|x\overline{y}^k|/\delta _k)+k^2\mathrm{\Psi }(|y\overline{y}^k|/\delta _k)+k.$$

Clearly $`v_y^k(y)=k`$ and $`v_y^k<k`$ on $`\mathrm{\Omega }\{y\}`$. By continuity and compactness there is a finite number (dependent on $`k`$) of points $`y_j^k\mathrm{\Omega }`$ and corresponding balls $`B_{r_{jk}}(y_j^k)`$ forming a covering of $`\mathrm{\Omega }`$ such that each $`v_{y_j}^k>k1`$ in $`B_{r_{jk}}(y_j^k)`$.

Let $`v^k=\mathrm{max}\{\mathrm{max}_j\{v_{y_j}^k\},1\}`$. Then $`v^k`$ is subharmonic, $`kv^k1`$ in $`\overline{\mathrm{\Omega }}`$ and $`v^k>k1`$ on $`\mathrm{\Omega }`$.

Let $`\mathrm{\Omega }_\delta =\{x\mathrm{\Omega }:\mathrm{dist}(x,\mathrm{\Omega })\delta \}`$. If $`z\mathrm{\Omega }_{\delta _k}`$ then $`|z\overline{y}_j^k|2\delta _k`$ for all $`\overline{y}_j^k`$ and $`v_{y_j}^k(z)k^2\mathrm{\Psi }(2)+k^2\mathrm{\Psi }(1)+k`$ implying that, for sufficiently large $`k`$, $`v^k=1`$ in $`\mathrm{\Omega }_{\delta _k}`$. This allows, for sufficiently large $`k`$, the estimates

$`{\displaystyle _\mathrm{\Omega }}v^k𝑑s`$ $`>|\mathrm{\Omega }|(k1),`$
$`{\displaystyle _\mathrm{\Omega }}v^k𝑑x`$ $`={\displaystyle _{\mathrm{\Omega }_{\delta _k}}}1+{\displaystyle _{\mathrm{\Omega }\mathrm{\Omega }_{\delta _k}}}v^k𝑑x|\mathrm{\Omega }|+k|\mathrm{\Omega }\mathrm{\Omega }_{\delta _k}|.`$
The result follows by choosing the sequence $`(\delta _k)_{k>2}`$ in such a way that it implies bundedness of the last expression. ∎

### 3.2. Beardon’s property

The converse of Beardon’s property follows from three facts: a standard mollification argument that shows the result holds if it holds for smooth functions, a computation for smooth functions that shows that Beardon’s inequalities holding for balls centered at a point $`x_0`$ imply the Laplacian is non–negative at $`x_0`$, provided $`\kappa <\kappa _1`$, and a well known device to extend the latter result to $`\kappa \kappa _1`$. We present these as three lemmas.

###### Lemma 3.2.

Assume $`u`$ is a continuous function for which (1.5) holds and is not subharmonic. Then a smooth subharmonic function exists for which the same properties hold.

###### Proof.

Let $`u`$ be such a continuous function. Then, in some point $`\overline{x}`$ in the interior of some ball $`B\mathrm{\Omega }`$, $`u`$ is bigger then the solution $`h`$ of Dirichlet’s problem for the the Laplace equation in $`B`$ having $`u`$ as boundary data. Denote $`\gamma u(\overline{x})h(\overline{x})>0`$.

Let $`\delta =\mathrm{dist}(B,\mathrm{\Omega })`$ and $`0<ϵ<\delta /4`$ and consider the usual mollifiers $`\rho _ϵ`$ supported in $`B_ϵ(0)`$, $`\phi `$ a continuous cut-off function which is $`0`$ in $`\mathrm{\Omega }\mathrm{\Omega }_{\delta /4}`$ and $`1`$ in $`\mathrm{\Omega }_{\delta /2}`$ and extend $`u`$ by $`0`$ in the complement of $`\mathrm{\Omega }`$.

Then Beardon’s inequality holds for $`\rho _ϵ(u\phi )`$ in $`\mathrm{\Omega }_{3\delta /4}`$. To check this statement let $`B_r(a)\mathrm{\Omega }_{3\delta /4}`$ and notice that

$$\begin{array}{cc}& \frac{1}{r^{n1}\kappa ^{n1}\omega _n}_{B_{kr}(a)}\rho _ϵ(u\phi )𝑑S\hfill \\ & =\frac{1}{r^{n1}\kappa ^{n1}\omega _n}_{B_{kr}(a)}\left(_{B_ϵ(0)}\rho _ϵ(y)u(xy)\phi (xy)𝑑y\right)𝑑S(x)\hfill \\ & =\frac{1}{r^{n1}\kappa ^{n1}\omega _n}_{B_ϵ(0)}\left(\rho _ϵ(y)_{B_{kr}(a)}u(xy)𝑑S(x)\right)𝑑y\hfill \\ & \frac{n}{r^n\omega _n}_{B_ϵ(0)}\left(\rho _ϵ(y)_{B_r(a)}u(xy)𝑑x\right)𝑑y\hfill \\ & =\frac{n}{r^n\omega _n}_{B_r(a)}\left(_{B_ϵ(0)}\rho _ϵ(y)u(xy)𝑑y\right)𝑑x\hfill \\ & =\frac{n}{r^n\omega _n}_{B_r(a)}\rho _ϵ(u\phi )𝑑x.\hfill \end{array}$$

As $`\rho _ϵ(u\phi )u`$ uniformly on $`B\{\overline{x}\}`$ as $`ϵ0`$, we have, for sufficiently small $`ϵ>0`$, $`\rho _ϵ(u\phi )\gamma /2<h`$ on $`B`$ and $`\rho _ϵ(u\phi )(\overline{x})\gamma /2>h(\overline{x})`$.

Hence $`\rho _ϵ(u\phi )\gamma /2`$ has all the desired properties in $`\mathrm{\Omega }_{3\delta /4}`$. ∎

###### Lemma 3.3.

Let $`uC^2(B_R(a))`$ and assume (1.5) holds for $`u`$ for some $`\kappa `$ satisfying $`0<\kappa <\kappa _1`$ and all $`r`$ satisfying $`0<r<R`$. Then $`\mathrm{\Delta }u(a)0`$.

###### Proof.

Assume $`a=0`$. Beardon’s inequality can be rewritten as

$$\begin{array}{cc}\hfill 0& \frac{n}{r^n}_0^r\left(_{B_t(0)}u(x)𝑑S(x)\right)𝑑t\frac{1}{r^{n1}\kappa ^{n1}}_{B_{\kappa r}(0)}u(x)𝑑S(x)\hfill \\ & =\frac{n}{r^n}_0^r\left(t^{n1}_{B_1(0)}u(tx)𝑑S(x)\right)𝑑t_{B_1(0)}u(\kappa rx)𝑑S(x).\hfill \end{array}$$

Let now

$$\psi (r)=n_0^r\left(t^{n1}_{B_1(0)}u(tx)𝑑S(x)\right)𝑑tr^n_{B_1(0)}u(\kappa rx)𝑑S(x).$$

As $`\psi 0`$ we have that

$$\underset{¯}{\mathrm{lim}}_{r0}\frac{\psi (r)}{r^{n+2}}0,$$

which in turn equals

$$\underset{r0}{lim}\frac{\psi ^{}(r)}{(n+2)r^{n+1}},$$

provided the latter limit exists. Now

$$\begin{array}{cc}\hfill \psi ^{}(r)& =nr^{n1}_{B_1(0)}u(tx)𝑑S(x)nr^{n1}_{B_1(0)}u(\kappa rx)𝑑S(x)\hfill \\ & r^n_{B_1(0)}\kappa u(\kappa rx)x𝑑S(x)\hfill \\ & =nr^{n1}_{B_1(0)}\left[_{\kappa r}^r\frac{d}{dt}\left(u(tx)\right)𝑑t\right]𝑑S(x)r^n_{B_1(0)}\kappa u(\kappa rx)x𝑑S(x)\hfill \\ & =nr^{n1}_{B_1(0)}\left[_{\kappa r}^ru(tx)x𝑑t\right]𝑑S(x)r^n_{B_1(0)}\kappa u(\kappa rx)x𝑑S(x)\hfill \\ & =nr^{n1}_{\kappa r}^r\left[_{B_1(0)}u(tx)x𝑑S(x)\right]𝑑tr^n_{B_1(0)}\kappa u(\kappa rx)x𝑑S(x)\hfill \\ & =nr^{n1}_{\kappa r}^r\left[_{B_1(0)}t\mathrm{\Delta }u(tx)𝑑x\right]𝑑tr^{n+1}_{B_1(0)}\kappa ^2\mathrm{\Delta }u(\kappa rx)𝑑x,\hfill \end{array}$$

where the last step follows by applying the divergence theorem. Thus,

$$\begin{array}{cc}\hfill \underset{r0}{lim}\frac{\psi ^{}(r)}{r^{n+1}}& =\underset{r0}{lim}\left[\frac{n}{r^2}_{\kappa r}^r\left[_{B_1(0)}t\mathrm{\Delta }u(tx)𝑑x\right]𝑑t\kappa ^2_{B_1(0)}\mathrm{\Delta }u(\kappa rx)𝑑x\right]\hfill \\ & =\underset{r0}{lim}\left[\frac{n}{r^2}_{\kappa r}^r\left[_{B_1(0)}t\mathrm{\Delta }u(tx)𝑑x\right]𝑑t\right]\kappa ^2_{B_1(0)}\mathrm{\Delta }u(0)𝑑x.\hfill \end{array}$$

Since

$$\underset{r0}{lim}\left[\frac{n}{r^2}_{\kappa r}^r\left[_{B_1(0)}t\mathrm{\Delta }u(tx)𝑑x\right]𝑑t\right]=\frac{n}{2}_{B_1(0)}\left(1\kappa \right)\mathrm{\Delta }u(0)𝑑x,$$

we finally obtain that

$$\left[\frac{n}{2}(1\kappa )\kappa ^2\right]\mathrm{\Delta }u(0)0,$$

and thus if $`\kappa (0,\kappa _1)`$, it follows that $`\mathrm{\Delta }u(0)0`$. ∎

The next lemma is only relevant in establishing the converse of Beardon’s property in the case $`n=1`$.

###### Lemma 3.4.

The previous lemma holds if $`\kappa =\kappa _1`$.

###### Proof.

Again assume $`a=0`$. Assume Beardon’s inequality holds for some smooth function $`u`$ with $`\kappa =\kappa _1`$. Let $`v(x)=|x|^2`$. Then

$$\frac{1}{\omega _n|\kappa |^{n1}r^{n1}}_{B_{\kappa r}}v(x)𝑑S=\kappa ^2r^2$$

and

$$\frac{n}{\omega _nr^n}_{B_r}v(x)𝑑x=_0^r\left(_{B_\rho }|\rho |^2𝑑S\right)𝑑\rho =\frac{n}{n+2}r^2.$$

Hence Beardon’s property holds for $`v`$ with $`\kappa =\kappa _0`$, hence it holds for $`u+ϵv`$ for any $`ϵ>0`$ and with $`\kappa =\kappa _0`$ and consequently $`\mathrm{\Delta }u(0)+2nϵ=\mathrm{\Delta }(u+ϵv)(0)0`$ for any $`ϵ>0`$. But then $`\mathrm{\Delta }u(0)0`$. ∎

From the proof of Lemma 3.3 we have that if $`_a(r)𝒮_a(\kappa r)`$ for $`\kappa >\kappa _1`$ then

$$\left[\frac{n}{2}(1\kappa )\kappa ^2\right]\mathrm{\Delta }u(a)0,$$

and we again conclude that $`\mathrm{\Delta }u(a)0`$. This of course already followed from Theorem 1.3 (of which this argument gives a different proof). On the other hand, we have that for $`u(x)=|x|^{2p}`$ ($`p`$ integer)

$$_0(r)=\frac{nr^{2p}}{2p+n},\text{ while }𝒮_0(r)=\kappa ^{2p}r^{2p},$$

giving that $`_0(r)𝒮_0(\kappa r)`$ provided

$$\kappa \left(\frac{n}{2p+n}\right)^{1/(2p)},$$

which converges to one as $`p\mathrm{}`$. Thus we see that the best possible value of $`\kappa `$ for this property of subharmonic functions is one.