format: text

text:
# Radio Halos and Relics in Clusters of Galaxies and Detection Statistics

## 1. VLA observations

In a recent paper, Giovannini, Tordi, & Feretti (1999) presented a sample of halo and relic cluster radio sources found in the NRAO VLA Sky Survey (NVSS, Condon et al., 1998). They used as cluster sample the X-Ray-brightest Abell-type clusters (XBACs) presented by Ebeling et al. (1996). The complete sample contains 177 clusters; some more clusters are present but do not belong to the complete sample because of their too low galactic latitude or too high redshift (larger than 0.2; see Giovannini, Tordi, & Feretti (1999) for a detailed discussion).

Six clusters where a central halo and/or a peripheral relic source was found for the first time, were observed by us with $``$ 5 hrs long observations each, with the VLA in the C and D configurations. They are: A115, A520, A773, A1664, A2254 and A2744.

In all clusters the presence of the extended source was confirmed. The new results are in good agreement with NVSS data, confirming the very high reliability of the NVSS survey also for low brightness extended sources. The size of diffuse sources is larger in the new maps, as expected because of the better sensitivity reached in pointed observations.

Adding to the list given by Giovannini, Tordi, & Feretti (1999) all clusters where the presence of a halo and/or a relic source was reported in literature, we have a sample of $``$ 40 cluster of galaxies with a diffuse source.

## 2. Frequency of halos and relics with the X-Ray luminosity

The number of new diffuse halos detected in the NVSS is rather high, especially considering that these sources are characterized by a steep spectrum and therefore that they are better imaged at frequencies lower than 1.4 GHz and also taking into account the limited surface brightness sensitivity.

We used these data to derive the occurrence of a diffuse source in clusters of galaxies of different X-ray luminosity. In Table 1 we report the results obtained for the complete sample (Col. 4) and the total sample (Col. 5). The X-ray luminosity is computed in the Rosat band (0.1 - 2.4 Kev).

The percentage of clusters with halos and relics increases with X-Ray luminosity. There is evidence that this increase is more relevant for central halo sources than for peripheral relic sources. The clusters hosting radio halos and/or relics have an X-Ray luminosity significantly higher than the other clusters (confidence level $`>`$ 99.9% estimated with a Kolmogorov-Smirnov test).

## 3. Frequency of halos and relics with the redshift

To investigate if the frequency of diffuse sources depends on the redshift, we have considered only the 19 clusters with an X-ray luminosity $`>`$ 9 $`\times `$ 1044 erg s-1 which are always visible in the complete sample (z $``$ 0.2). We find that 55% of clusters with 0.15 $`<`$ z $``$ 0.2 have a diffuse source, to be compared with 33% of clusters with 0.044 $``$ z $``$ 0.1. The difference is not statistically significant because the too small number of clusters; a deeper and larger sample is necessary to investigate this point.

## References

Condon, J. J., Cotton, W. D., Greisen, E. W., Yin, Q. F., Perley, R.A., Taylor, G. B., & Broderick, J. J.: 1998, AJ, 115, 1693

Ebeling, H., Voges, W., Bohringer, H., Edge, A. C., Huchra, J. P., & Briel, U. G.: 1996, MNRAS, 283, 1103

Giovannini, G., Tordi, M., & Feretti L.: 1999, New Astron., 4, 141


text:
# Removing radio interference from contaminated astronomical spectra using an independent reference signal and closure relations

## 1 Introduction

The growing level of radio frequency interference (RFI) is a recognized problem for research in radio astronomy. Fortunately, the technological advances that are giving rise to the increasing background of radiation – through increased telecommunications and wide-spread use of high speed electronics – are also providing some of the tools necessary for separating astronomical signals from undesirable RFI contamination. New radio telescopes will necessarily have RFI suppression, excision and cancellation algorithms intrinsic to their designs. No one technical solution will make radio observations immune to interference; successful mitigation is most likely to be a hierarchical or progressive approach throughout the telescope, combining new instrumentation and algorithms for signal conditioning and processing Ekers & Bell (2000); Ellingson (1999).

Techniques from the communication industry that are finding application in radio astronomy experiments include (1) adaptive beam forming with array telescopes that steer nulls of the instrument reception pattern in the directions of sources of RFI Hampson et al. (1998); Ellingson & Hampson (2000); Leshem et al. (2000); Smolders & Hampson (2000); Kewley et al. (2000), (2) parametric signal modeling techniques, where the RFI signal is received and decoded to obtain a high signal-to-noise ratio (SNR) reference signal for subtraction from the astronomical data Ellingson et al. (2000); Leshem & van der Veen 1999a ; Leshem & van der Veen 1999b , and (3) adaptive filtering using a reference horn to obtain a high SNR copy of the RFI for real-time cancellation from the signal path ahead of the standard radio astronomical backend processors Barnbaum & Bradley (1998).

This paper describes an intuitive but powerful RFI cancellation technique that is suitable for radio spectroscopy where time-averages are recorded. The method requires computation of cross power spectra between the RFI contaminated astronomical signals and high signal-to-noise ratio RFI “reference signals” obtained from a receiving system that senses the RFI but not the astronomical signal. The correction term that removes the unwanted RFI is computed from closure relations obeyed by the RFI signal. The test applications reported here derived the reference signal either from a separate horn antenna aimed at the RFI source or from a second feed horn at the focus of the Parkes telescope, as illustrated in Fig 1. For these experiments, we recorded digitally sampled baseband signals from two polarizations for both the reference and astronomy feeds, and then we performed the cross correlations in software off-line. However, the method could use correlation spectrometers of the sort already in use at radio observatories. With minor design enhancements, future generation correlators could incorporate this technique with the additional benefit of the faster electronics, greater digital precision, and higher data rates that technological advance promises.

There are a number of advantages to performing the RFI in a “post-correlation” stage. Foremost is that the RFI subtraction remains an option in the data reduction path, rather than a commitment made online and permanently. Furthermore, the correlation method is effectively a coherent subtraction, since the correlation functions retain the information describing relative phase between the RFI entering in the astronomy data stream and the RFI entering the reference antenna. As we show in this paper, this means that the RFI noise power is largely subtracted, leaving only the usual components of system noise.

This paper provides a description of the post-correlation RFI cancellation technique and illustrates its success with data from the Parkes telescope. A mathematical overview shows (1) why unknown multipaths do not cause the algorithm to break down, (2) how to simply construct a suitable RFI reference spectrum, and (3) how to build an inverse filter to obtain immunity to low signal levels at frequencies that suffer destructive interference by multipathing in the reference horn signal path.

## 2 Mathematical description of the method

In our mathematical model, we make the assumption that the RFI source emits a single signal $`i(t)`$. (At the RFI source, the signal from a single power amplifier feeds an antenna of unknown, but irrelevant, polarization.) The RFI that appears differently in the recorded data channels at Parkes has experienced scatterings with different path lengths and amplitudes, so that the received signals are linear sums of time delayed versions of the original broadcast $`i(t)`$.

In fact, the model is applicable to multiple interferers within the spectrometer band, provided they do not overlap in frequency.

Consider for the moment a single interferer, which propagates through the four signal paths $`s_1`$, $`s_2`$, $`s_3`$, and $`s_4`$ that will be processed: there are two astronomical channels, $`s_1(t)`$ and $`s_2(t)`$, which convey the voltages $`a_A(t)`$ and $`a_B(t)`$ from the celestial sources for the two independent polarizations from the Parkes Telescope receiver along with contamination from the RFI signal $`i(t)`$. Radiation from astronomical sources may be partially polarized, causing $`a_A(t)`$ and $`a_B(t)`$ to be correlated to some degree. The two reference channels carry $`s_3(t)`$ and $`s_4(t)`$, containing the representations of $`i(t)`$ but negligible signal from the celestial sources. For example, the measured signal in channel 1 comprises channel noise $`n_1(t)`$, plus the convolutions of the impulse responses for the multiple scattering paths for the interference signal and the true astronomy signal: $`s_1(t)=H_A(t)\mathrm{}a_A(t)+H_1(t)\mathrm{}i(t)+n_1(t)`$. Here $`H_A(t)`$ and $`H_1(t)`$ are impulse responses for the astronomy signal path and interference respectively, and $`\mathrm{}`$ is the convolution operator.

For many purposes, an intuitive picture of the multipathing results from considering the scattering sites to be achromatic mirror-like scatterers each with relative effective areas $`G_{A,j}`$ and $`G_{1,k}`$, and attaching the path delay to each separate version of the astronomy and RFI signals, $`s_1(t)=_jG_{A,j}a_A(t\tau _{1,j})+_kG_{1,k}i(t\tau _{1,k})+n_1(t)`$. The time delays $`\tau `$ are determined by the different path lengths $`L_{1,j}`$ to give $`\tau _{1,j}=L_{1,j}/c`$, where $`c`$ is the speed of light.

All the signal paths are vulnerable to corruption by stochastic noise. The noise terms $`n_1(t)`$, $`n_2(t)`$, $`n_3(t)`$, and $`n_4(t)`$ should ideally be uncorrelated among the different data paths. Unfortunately, in real astronomical systems, there is likely to be low level coupling between the two orthogonal polarizations of a feed horn or common stray radiation pickup from spillover that will make a weakly correlated noise floor in some of the cross power spectra. This will form a systematic limitation to the accuracy of the subtraction.

In this experiment, the goal is to explore the usefulness of cross-correlation spectra to correct for the effects of RFI in time-averaged spectra. These spectra are products of scaled sums of the of the Fourier transforms of the astronomy signals, the RFI and noise. In the tests with real data in section 4, we compute estimates of the complex spectra

$`S_1(f)`$ $`=`$ $`g_AA_A+g_1I+N_1`$
$`S_2(f)`$ $`=`$ $`g_BA_B+g_2I+N_2`$
$`S_3(f)`$ $`=`$ $`g_3I+N_3`$
$`S_4(f)`$ $`=`$ $`g_4I+N_4`$ (1)
from Fourier transforms of finite length time series of discrete samples of the four signals. The transforms contain contributions from the celestial sources $`A_A(f)`$ and $`A_B(f)`$, the RFI $`I(f)`$, and the noise in each channel $`N_i(f)`$, modulated by the associated complex voltage gains, which are the Fourier transforms of the impulse responses $`H(t)`$. The gains for each channel separate into dependencies on (1) the path delay $`L/c`$, which appears in a frequency dependent phase term, according to the Shift Theorem of Fourier transforms, and (2) a possible additional frequency dependence $`g(f)`$ of each delay path:

$$g_1(f)=g_{1,k}(f)e^{i2\pi fL_{1,k}/c}$$
(2)
These complex gain and delay factors are sufficiently general to include complicated scatterers and propagation through dispersive and lossy media.

The real power spectra for the four data channels have the following form, once terms that average toward zero are omitted and the complex gains are assumed to be constant over the time span for which the spectra are computed:

$`P_1(f)`$ $`=`$ $`<S_1S_1^{}>`$
$`=`$ $`|g_A|^2<|A_A|^2>`$
$`+|g_1|^2<|I|^2>+<|N_1|^2>`$
$`P_2(f)`$ $`=`$ $`<S_2S_2^{}>`$
$`=`$ $`|g_B|^2<|A_B|^2>`$
$`+|g_2|^2<|I|^2>+<|N_2|^2>`$
$`P_3(f)`$ $`=`$ $`<S_3S_3^{}>`$
$`=`$ $`|g_3|^2<|I|^2>+<|N_3|^2>`$
$`P_4(f)`$ $`=`$ $`<S_4S_4^{}>`$ (3)
$`=`$ $`|g_4|^2<|I|^2>+<|N_4|^2>`$
We use the superscript symbol  to represent complex conjugation and the $`<\mathrm{}>`$ notation to signify averages over an integration time $`t_{int}`$; in the tests we describe in section 4, we find the method is effective for $`t_{int}`$ as long as $``$1 sec. We adopt a normalization where the power levels $`<|A_A|^2><|I|^2><|N_1|^2>`$ so that, for example, in data channel 1 the signal to noise ratio $`SNR_1|g_A|^2`$, and the interference to noise ratio $`INR_1|g_1|^2`$.

The goal of the RFI cancellation will be to form estimates of the $`|g_1|^2<|I|^2>`$ and $`|g_2|^2<|I|^2>`$ terms and then subtract them from $`P_1(f)`$ and $`P_2(f)`$, while leaving the astronomical signal (and noise) behind.

This discussion has assumed that the complex gain terms are constant over the integration time $`t_{int}`$, allowing us to separate the interference from the gain in expressions such as:

$$<g_1Ig_1^{}I^{}>=|g_1|^2<|I|^2>$$
(4)
In anticipation of the discussion of section 4, we note that this assumption will fail for extended integration times since the scattering paths that lead the RFI signal to the telescope feed will change, resulting in loss of precision in the cancellation scheme and and leading to substantial residuals in the corrected spectra.

The complex cross power spectra for all combinations of the four data channels have the following form when the leading contributions are retained:

$`C_{12}(f)`$ $`=`$ $`<S_1S_2^{}>`$
$`=`$ $`g_Ag_B^{}<A_AA_B^{}>`$
$`+g_1g_2^{}<|I|^2>`$
$`+g_1<IN_2^{}>+g_2^{}<N_1I^{}>`$
$`+<N_1N_2^{}>`$
$`C_{ij}(f)`$ $`=`$ $`<S_iS_j^{}>`$ (5)
$`=`$ $`g_ig_j^{}<|I|^2>`$
$`+g_i<IN_j^{}>+g_j^{}<N_iI^{}>`$
$`+<N_iN_j^{}>`$
$`\mathrm{for}\mathrm{i}\mathrm{j},\mathrm{j}>2`$
In order to cancel the dominant RFI terms in the power spectra Eqns 3, we need to compute quantities of the form $`|g_1|^2<|I|^2>`$, which can be subtracted from the measured $`<|S_1|^2>`$. One possibility would be a combination of auto and cross power spectra of the form

$`|g_1|^2<|I|^2>`$ $`=`$ $`{\displaystyle \frac{g_1g_3^{}g_1^{}g_3}{g_3^{}g_3}}<|I|^2>`$
$``$ $`{\displaystyle \frac{|C_{13}|^2}{<|S_3|^2>}}`$
$``$ $`{\displaystyle \frac{|C_{13}|^2}{|g_1|^2<|I|^2>+<|N_3|^2>}}`$
However, the problem encountered with this approach is that $`<|S_3|^2>`$ is biased by the term $`<|N_3|^2>`$, which averages over time to the total power spectrum of the noise in data channel 3. This might be compensated by calibrating the noise spectrum, in order to improve the estimate of $`|g_3|^2<|I|^2>=<|S_3|^2><|N_3|^2>`$ to use in the denominator of Eqn 2.

An alternate combination that avoids this bias forms the estimate of $`|g_1|^2<|I|^2>`$ from three cross power spectra, in which the noise terms have the form $`<N_iN_j^{}>`$ and $`<IN_j^{}>`$ which average toward zero as the integration time increases:

$`|g_1|^2<|I|^2>`$ $`=`$ $`{\displaystyle \frac{g_1g_3^{}g_1^{}g_4}{g_3^{}g_4}}<|I|^2>`$ (7)
$``$ $`{\displaystyle \frac{C_{13}C_{14}^{}}{C_{34}^{}}}`$
$`|g_2|^2<|I|^2>`$ $`=`$ $`{\displaystyle \frac{g_2g_3^{}g_2^{}g_4}{g_3^{}g_4}}<|I|^2>`$ (8)
$``$ $`{\displaystyle \frac{C_{23}C_{24}^{}}{C_{34}^{}}}`$
$`g_1g_2^{}<|I|^2>`$ $`=`$ $`{\displaystyle \frac{g_1g_4^{}g_2^{}g_3}{g_3g_4^{}}}<|I|^2>`$ (9)
$``$ $`{\displaystyle \frac{C_{14}C_{23}^{}}{C_{34}}}`$
The expressions involving the measured cross power spectra are “approximate,” since the cross power spectra result from finite integrations, and the noise terms will limit the precision of the cancellation.

The occurrence of the $`C_{34}`$ term in the denominator for Eqns 7-9 indicates there will be a problem in implementing a correction scheme in frequency ranges where $`C_{34}`$ becomes small or zero. In many situations when the signal to noise ratio for the spectra in the numerator is very high, the cross power spectra in the numerator will also be small or zero whenever $`C_{34}`$ is small, so that divergence will be canceled. A simple means to avoid division by small numbers in the presence of noise in this kind of situation is to create a filter of the following form:

$`|g_1|^2<|I|^2>`$ $``$ $`{\displaystyle \frac{C_{13}C_{14}^{}C_{34}}{\psi (f)+C_{34}C_{34}^{}}}`$ (10)
$`|g_2|^2<|I|^2>`$ $``$ $`{\displaystyle \frac{C_{23}C_{24}^{}C_{34}}{\psi (f)+C_{34}C_{34}^{}}}`$ (11)
$`g_1g_2^{}<|I|^2>`$ $``$ $`{\displaystyle \frac{C_{14}C_{23}^{}C_{34}^{}}{\psi (f)+C_{34}C_{34}^{}}}`$ (12)
where $`\psi (f)`$ is the square of the power spectrum of the noise present in $`C_{34}(f)`$. Whenever $`\psi (f)`$ becomes small compared to $`C_{34}(f)`$, the expressions Eqns 10-12 revert to Eqns 7-9. When the noise exceeds the signal power in $`C_{34}`$, the computed correction tends to zero. In practice, during the test described here, a constant $`\psi _0`$ was used in place of $`\psi (f)`$. Alternatively, division by zero can be avoided by testing the amplitude of $`C_{34}(f)`$ for significance above a noise threshold and setting the correction to zero when the significance criterion is not met.

These corrections for the the auto correlation spectra (Eqns 7 and 8) are expected to be real valued. Therefore a logical test for the accuracy of the correction is that phases computed for the frequencies of strong RFI signal should be close to zero. In fact, this is a statement of phase closure. Note that the denominators of Eqns 10-12 are purely real, and the numerators, such as $`C_{13}C_{14}^{}C_{34}`$, form logical triangles for computing closure phases.

An amplitude closure relation

$$\frac{C_{13}C_{24}}{C_{23}C_{14}}=\frac{g_1g_3^{}g_2g_4^{}}{g_2g_3^{}g_1g_4^{}}=1$$
(13)
can also be constructed and tested. It too will suffer from divergence of the quotient in frequency ranges where the cross power spectra are noisy and $`C_{23}`$ and $`C_{14}`$ have small amplitude. Here we avoid including $`C_{12}`$, which typically has significant cross correlated power in addition to the RFI signal, due either to polarized celestial flux or cross talk between the channels.

## 3 Noise and the Accuracy of the Cancellation

In this section we assess the importance of the Interference-to-Noise ratio. First we expand the autocorrelation spectra, keeping all cross terms, including those that average toward zero. Then $`P_1(f)`$ becomes

$`P_1(f)`$ $`=`$ $`|g_A|^2<|A_A|^2>+|g_1|^2<|I|^2>+2\mathrm{R}\mathrm{e}\left[g_1g_A^{}<IA^{}>+g_1<IN_1^{}>\right]`$ (14)
$`+2\mathrm{R}\mathrm{e}\left[g_A<AN_1^{}>\right]+<|N_1|^2>`$
The complex correction spectra were described by equations 7, 8 and 9. Including the cross terms and noting that when

$$g_3\mathrm{and}\mathrm{g}_4\mathrm{g}_1\mathrm{and}\mathrm{g}_2$$
(15)
and the interference power to noise power ratios

$`INR_3`$ $`=`$ $`{\displaystyle \frac{|g_3|^2<|I|^2>}{<|N_3|^2>}}|g_3|^21`$
$`INR_4`$ $`=`$ $`{\displaystyle \frac{|g_4|^2<|I|^2>}{<|N_4|^2>}}|g_4|^21`$ (16)
then the correction $`CX_1`$ for $`P_1(f)`$ becomes:

$`CX_1`$ $``$ $`{\displaystyle \frac{C_{13}C_{14}^{}}{C_{34}^{}}}`$ (17)
$``$ $`|g_1|^2<|I|^2>+2\mathrm{R}\mathrm{e}\left[g_1g_A^{}<IA^{}>+g_1<IN_1^{}>\right]`$
$`+{\displaystyle \frac{g_Ag_1^{}}{g_3^{}}}<AN_3^{}>+{\displaystyle \frac{g_A^{}g_1}{g_4^{}}}<A^{}N_4>`$
$`+{\displaystyle \frac{g_1^{}}{g_3^{}}}<N_1N_3^{}>+{\displaystyle \frac{g_1}{g_4}}<N_1^{}N_4>{\displaystyle \frac{|g_1|^2}{g_3^{}g_4}}<N_3^{}N_4>`$
The terms in Eqn 17 involving $`I`$ all appear in the power spectrum of Eqn 14, so that application of this correction $`CX_1`$ to $`P_1(f)`$ leads to a result with no residual contamination by the RFI:

$`P_1(f)CX_1`$ $`=`$ $`|g_A|^2<|A_A|^2>+<|N_1|^2>`$ (18)
$`+2\mathrm{R}\mathrm{e}\left[g_A<AN_1^{}>\right]{\displaystyle \frac{g_Ag_1^{}}{g_3^{}}}<AN_3^{}>{\displaystyle \frac{g_A^{}g_1}{g_4^{}}}<A^{}N_4>`$
$`{\displaystyle \frac{g_1^{}}{g_3^{}}}<N_1N_3^{}>{\displaystyle \frac{g_1}{g_4}}<N_1^{}N_4>+{\displaystyle \frac{|g_1|^2}{g_3^{}g_4}}<N_3^{}N_4>`$
The complete cancellation of the RFI terms is consistent with the concept that post correlation subtraction is equivalent to coherent subtraction of the RFI electric field $`i(t)`$ in the time domain, which should leave no trace of the RFI signal nor an increase in noise power. Provided the $`INR`$ for the reference horn channels is substantially greater than the $`INR`$ for the Parkes feed channels, the noise terms due to $`N_3`$ and $`N_4`$ will be smaller by a factor of order $`\beta \sqrt{INR_3/INR_1}|g_3|/|g_1|`$ than the normal noise contributions arising from $`N_1`$ plus the astronomical signal power. The terms in Eqn 18 other than $`|g_A|^2<|A_A|^2>`$ and $`<|N_1|^2>`$, such as $`(g_A^{}g_1/g_4^{})<A^{}N_4>|g_A|\beta ^1t^{1/2}`$ and $`(g_1^{}/g_3^{})<N_1N_3^{}>\beta ^1t^{1/2}`$, average toward zero (in inverse proportion to the square root of the integration time), provided the noise in the signal channels is independent. The next higher order terms, which are not included in Eqn 18, have dependencies such as $`<IN_1^{}><N_1^{}N_4>/g_4<|I|^2>t^1`$ and $`g_1^{}<IN_3^{}><IN_1>/g_3^{}<|I|^2>\beta ^1t^1`$, which converge toward zero faster than the dominant noise terms.

## 4 The test data

The astronomical dataset used in testing these algorithms is a dual linear polarization data stream from the CSIRO ATNF 64m telescope at Parkes in Australia. One configuration has two polarizations from the central beam of the the Parkes multibeam receiver Staveley-Smith et al. (1996) and two polarizations from a reference horn aimed at an interfering source Bell et al. (2000). A second configuration uses both polarizations from two beams of the Parkes Multibeam system, which are directed at slightly different areas of the sky. The center frequency of the datasets was 1499 MHz in each case.

The datasets we used are labeled SRT00501, SRT00502 and SRT00601 Bell et al. (2000). The main interfering source is a NSW Government digital point-to-point microwave link. Examples of the time-averaged spectrum for the RFI signal are shown in Fig 2. Further details are available in the ACA databases Australian\_Communications\_Authority (1998); Sarkissian (1999). The 4 signals were downconverted to base band and passed through 5 MHz low pass filters. Each signal was then digitized with 2-bit precision at a 20MHz sampling rate to achieve a factor 2 oversampling and recorded using the CPSR recorder van Straten et al. (2000).

The data processing for these tests simulates a radio astronomy backend by computing power spectra and cross power spectra in software. The sampled voltages are treated in 8192 sample blocks (410$`\mu `$sec durations). Fourier transformation of each block yields a spectrum of 4096 independent complex coefficients.

Examples of 25 second averages of the power spectra are shown in Fig 3. The RFI spectrum is significantly different in the spectra for all four data channels. Since the data are 2$`\times `$ oversampled, we kept the total power and cross power spectra of the lower 2048 complex Fourier coefficients. The RFI subtraction steps were performed on these spectra, as illustrated below, and subsequently the spectra were block averaged by 4 to keep power spectra of length 512 spectral channels for further calibration and display at 9.8 kHz resolution.

Examples of the cross power spectra are shown in Fig 4. Both $`C_{12}`$ and $`C_{34}`$ have broadband correlated power, which is clear in the integrated spectra both as a significant non-zero amplitude and as a well defined trend in phase across the band. The phase gradients that are clearly visible in the $`C_{13}(f),C_{14}(f),C_{23}(f)`$, and $`C_{24}(f)`$ spectra indicate differential path delays. The delay for $`C_{13}(f),C_{14}(f)`$ and for the channel numbers above 300 in $`C_{23}(f)`$, and $`C_{24}(f)`$ amounts to 15 time steps (0.75$`\mu `$sec). The steeper phase gradient for channels below 300 in $`C_{23}(f)`$, and $`C_{24}(f)`$ imply path delays of 60 time steps $`=3\mu `$sec or a path length of $``$900m. Apparently, the frequency dependent reception pattern of the far sidelobes of the Parkes dish, coupled with complicated scatterers in the field, can lead to quite complex spectral dependence for the RFI. A strength of our method is the simplicity with which it handles this complex multipathing.

Examples of the complex correction spectra defined in Eqns 7-9 are shown in Fig 5. The spectra for Correction A and Correction B have zero phase over the spectral range where the signal is significantly non-zero, confirming that phase closure applies.

Fig 6 shows averages of the amplitude closure quantity from Eqn 13, plotted as amplitude and phase across the spectrum. There is deterioration in the closure relation when the RFI signal is low, as expected.

Fig 2 compares three time-averaged power spectra representing the RFI signal sensed by the reference horn: power spectra for each of the two polarizations and one cross power spectrum. For use in correction schemes, the cross power spectrum $`C_{34}(f)`$ has the advantage that it is not contaminated by the positive bias of the receiver noise total power that is seen in the autocorrelation power spectra. In practice, the low level correlated signal in $`C_{34}(f)`$, due to actual cross correlated broadband power or due to correlated quantization noise, may limit the accuracy of the cross power spectrum as an estimate of the RFI.

Fig 7 compares uncorrected and corrected spectra for two scans. For the purposes of display, the spectra have been crudely passband calibrated by dividing the spectra by gain templates formed from the scan average of the total power spectra of scan SRT00601, which was recorded while the sky frequency for the Parkes data channels was tuned off the RFI frequency. Since this gain template is in common to the processing for both scans SRT00501 and SRT00502, some of the common structure in the spectra in Fig 7 results from this common gain template.

The time dependence of the RFI sensed by the Parkes telescope is displayed in image format in Fig 8, side by side with the spectra after the RFI subtraction. These spectra received the same processing as described for the averages in Fig 7, with the additional step of subtracting a third order polynomial spectral baseline from each time step. This additional step was done to remove some faint variations in the total power level that occurred from one integration to the next.

The cross polarized spectrum from the Parkes A$`\times `$B is shown in Fig 9. The figure includes the raw spectrum and the corrected spectrum after subtraction of the correction shown in Fig 5. Remaining in the corrected spectrum, there is a slow modulation of the power across the 5 MHz band, probably indicating that this power represents broad band noise signal that is scattered within the Parkes telescope structure with the same delay path lengths associated with traditional standing waves and path lengths of a few hundred feet.

The integration time over which the RFI corrections are applied is a critical parameter when the signal to noise ratio of the reference signal data path is low. If the integration time is too short, then some of the derived correction will be noise, and this will be folded into the resulting spectrum. On the other hand, if the integration time is long, the impulse response functions that couple the astronomy and RFI signals to the receiver will vary, and the mathematics in Eqn 4 and Eqns 7-9 will break down. Fig 10 shows a series of tests with a range of integration times ($`t_{int}`$ 8 msec to 8 sec) for application of the RFI subtraction. In all cases, the plots give the grand averages of the entire 25 seconds of the scan after application of the algorithm on the shorter data segments. The parameter $`\psi _o`$ has an appropriate value in each case to produce a stable noise level and representation of spectral features in the frequency range away from the RFI contamination. For these data, the algorithm was most effective for integration times of only a $``$ 1 sec or less. When treated on 8 sec averages, substantial RFI remains unsubtracted.

## 5 Discussion of Limitations

The wide range of delays for the scattered RFI led to problems in our initial experiments, which simulated an FX correlator with 2048 sample transforms. The RFI illustrated in Fig 4 has two strong components that are delayed by 15 and 60 time steps. A 60 time step delay is $``$3% of the time block being processed, so that the cross correlation is not being performed on fully overlapped data streams. Increasing the window to 8192 was adequate to reduce the residuals to a level that was barely visible above the noise in the corrected spectra. An additional but less significant improvement resulted from applying an constant delay offset of 35 time steps to the reference signal at the input to the “FX correlator” to make it closer to the average of the principal delays in the Parkes data channels. A traditional time-domain lag correlation spectrometer would not encounter this problem, provided a sufficient number of lags are allocated to fully cover the range of delays experienced by the RFI.

In these tests at Parkes, there is a possibility that some uncancelled signal may be present due to a second transmitter operating at these same frequencies. The reference horn was pointed at the stronger, nearby transmitter, while the second transmitter, located at approximately three times the distance, can be scattered into the Parkes Telescope signal paths without being sensed by the reference horn.

The coarse digitization of the RFI reference signals will eventually form a limitation to the precision of the subtraction. Crude quantization generates an artificial noise floor throughout the spectrum, effectively scattering power out of the narrow band RFI. Since both polarizations from the reference horn are recorded at relatively high signal to noise ratio, the quantization noise is also correlated between the two data channels, so that there is corruption of the cross power spectrum as well as for the autocorrelation spectra. (The implication is that the term $`<N_3N_4>`$ in Eqn 5 will not average to zero with increased integration.)

## 6 The Toxicity Test

A crucial requirement of an RFI subtraction algorithm is that it must leave the astronomical signal of interest unaltered. To test the current method, we added simulated galaxy signals to the two Parkes input data streams. Independent gaussian noise was filtered with a double-horned galaxy profile and the instrumental IF filter passbands and then injected into the data pipeline just before the correlation stage. Fig 11 shows a comparison between the RFI corrected spectra both with and without the galaxy. To highlight the difference between the injected signal and the output after RFI subtraction, the bottom panel shows (1) the difference between the output with and without the added signal and (2) the ratio of this difference to the synthetic galaxy profile added to the input. The RMS deviation of the ratio about unity is 0.005 for polarization A and 0.006 for polarization B. No systematic deviations are seen across the band, other than rises in noise level at the edges where the galaxy profile is approaching zero at the edge of the profile. The conclusion is that this method does no systematic harm to the astronomical signals.

## 7 The Parkes Two Feed Experiment

The mathematical description (Eqns. 3 to 12) can be equally well applied to a case that uses a second feed from the Parkes Telescope as the “reference horn.” Both feeds may receive astronomical signals, but since the feeds point different directions in the sky, these are independent astronomical signals, which will not correlate and therefore will not be subtracted from each other by this algorithm. The two feeds do sense the same RFI signal $`i(t)`$, although through different scattering paths. This is sufficient commonality that the cross power spectral approach should permit each feed to serve as the reference antenna for the other. Similar experiments have been reported Sault (1997); Kewley et al. (1999).

Fig. 12 shows a comparison between the autocorrelation spectra measured for the Parkes two feed experiment, before and after RFI subtraction. There is noticeable difference among the four channels in the effectiveness of the RFI subtraction. This probably results from hte bias breated by the broadband polarized flux or correlated noise in the cross polarized spectra $`3\times 4`$ that is being used as the “reference” for data channels 1 and 2. As shown in Fig. 13, this noise floor is higher in feed 2 ($`INR`$ 35:1) than in feed 1 ($`INR`$ 100:1), causing the rfi template spectra derived from the cross power spectra to be less faithful, which in turn leads to larger error in the correction spectra to be subtracted from feed 1. For comparison, the reference horn spectrum $`C_{34}(f)`$ in Fig 4 has $`INR`$ 1000:1

## 8 Conclusions

RFI subtraction can be performed using cross power spectra between the astronomy data channels and RFI “reference” channels. In principle, the reference channel can also be an astronomy channel provide it carries an astronomy signal that is uncorrelated with the astronomy in the channel that is being corrected.

The tests made at Parkes demonstrate that a specifically designed reference sensor provided a higher signal-to-noise ratio reference signal – and consequently cleaner cancellation – than that obtained from a second horn feed at the Parkes Telescope focus, whose principal function is to illuminate the Parkes dish.

A refinement will be to implement this scheme using two reference antennas that are spatially separated (as in the lower diagram of Fig 1) in order to avoid correlated noise contributions, while still obtaining as clear and stable path to the RFI source as possible. The cross power spectra from the two spatially separated antennas would form an optimal RFI reference spectrum $`C_{34}`$ for use in the Eqns 10-12. To avoid problems with differential delay causing loss of coherence in the reference signal, the spectrometer would need to operate with spectral resolution $`\mathrm{\Delta }f=(\mathrm{\Delta }\varphi /2\pi )c/L<<c/L`$ where $`\mathrm{\Delta }\varphi `$ is the allowable phase rotation across a spectrometer channel and $`L`$ is the spatial separation of the sensors. The 2.4 kHz spectrometer resolution emulated in software for the study reported here would allow spatial separations of up to 1200m, if $`\mathrm{\Delta }\varphi `$ is required to be less than 0.02$`\pi `$ radians.

There are a number of advantages to performing this type of “post-correlation” RFI subtraction:

1. Provided the required correlation products are recorded (i.e. the online system is capable of recording correlation functions with a sufficiently large number of delay lags), the RFI subtraction can be performed offline, where it remains an option in the data reduction path, rather than a commitment made online and permanently.

2. The method is not vulnerable to the effects of sporadic RFI, which hurt many algorithms that have an initialization period while they acquire the RFI signal and optimize their cancellation parameters.

3. Nor is the result influenced by changes in beam shape during adaptive nulling.

4. The correlation method is effectively a coherent subtraction, since the correlation functions retain the information describing relative phase between the RFI entering in the astronomy data stream and the RFI entering the reference antenna. We showed in section 3, this means that the RFI noise power is largely subtracted, leaving only system noise.

5. Generalization of the method to an array of telescopes is straightforward but demands additional correlator capacity. If there are two reference signal sensors, labeled “$`x`$” and “$`y`$,” that sense negligible astronomical signal, then their cross power spectrum $`C_{xy}(f)`$ containing a high INR signal can be used to correct any other power spectrum $`C_{ij}(f)`$ through the closure relation $`C_{ix}C_{jy}^{}/C_{xy}^{}`$. The $`i`$, $`j`$ indices can denote orthogonal or parallel polarizations drawn from any combinations of antennas in the array or auto correlation, when $`i=j`$.

6. A modification of the method can be applied to pulsar data streams in which a digital correlator replaces the narrow band filter bank used in compensating for pulse dispersion. One possible implementation would construct cross-power “coupling spectra” $`X(f)`$ that are valid for the time interval $`t_{int}`$ during which the “$`g`$ factors” of Eqn 1 are stable. The coupling spectrum is

$$X_{13}(f)=\frac{g_1^{}g_4}{g_3^{}g_4}\frac{C_{14}^{}}{C_{34}^{}}$$

where $`C_{14}`$ and $`C_{34}`$ are averaged for up to 1 second as appropriate. Then the correction $`CX_1`$ can be computed and applied to cancel RFI in measurements of $`P_1(f)`$ on shorter timescales:

$$CX_1=|g_1|^2<|I|^2>=X_{13}g_1g_3^{}<|I|^2>X_{13}<S_1S_3^{}>$$

where $`<S_1S_3^{}>=C_{13}`$. Alternatively one could construct the corrected time sequence $`s_1(t)`$ by subtracting the correction

$$SX_1=g_1I=X_{13}^{}g_3IX_{13}^{}S_3$$

from $`S_1(f)`$ and inverse Fourier transforming to obtain an RFI cancelled version of $`s_1(t)`$. Computationally efficient schemes could be implemented that include coherent dedispersion Hankins (1974) in the same transform operations as the RFI cancellation.

7. The method can be generalized to removal of solar radiation whose multi-path scattering effects give rise to the spectral “standing wave” problem. The important difference that the Sun generates a dual polarized signal with a variable polarized component; this will necessitate a processing path more akin to the sub-space decomposition Ellingson (1999); Leshem et al. (2000), in order to identify orthogonal components of the solar “rfi” signals.

8. The present generation of digital correlators, which typically operate with single bit to 9 level precision, could implement this method for use in moderate levels of RFI for testing and astronomical observation in the near future.

The disadvantages of the method are: (1) The data rates will be high, since the method requires a full-scale cross correlator, preferably with multi-bit precision to accept large SNR RFI reference signals, that must dump spectra after relatively short integrations of less than $``$1 sec. Of course, these data rates are lower than recording the full baseband, but they are substantially higher than a real-time adaptive filter approach, which would allow long spectral integrations once the RFI has been canceled. (2) When the interference is much stronger than the astronomical signal, a 1 or 2 bit sampler is “captured” so that the data stream consists primarily of $`\pm `$2, and the zero crossings are determined by the phase of the interferer. Under these conditions, the astronomical data will be largely lost. These applications will require correlators with greater digital precision.

The case where multiple interferers occupy the same frequency band will require a greater number of sensors and more correlator capacity devoted to processing the the RFI signals, as laid out in the analyses of Sault (1997) and Ellingson (1999) .

The authors are grateful to W. van Stratten, M. Bailes, S. Anderson, S. Ellingson, R. Sault, P. Perillat, R. Ekers, J. Bunton, L. Kewley, M. Smith and P. Sackett for helpful comments and discussion. F.B. is grateful to the ATNF in Epping, NSW, the Department of Astronomy at OSU, Columbus, OH, and to the IAS, Princeton, NJ for their hospitality while this work was done.


text:
# Time Operator for the Quantum Harmonic Oscillator: Resolution of an Apparent Paradox

## 1  Introduction

In the space spanned by the eigenstates of the coordinates operator $`q`$ and the momentum operator $`p`$, we consider the quantum harmonic oscillator described by the following Hamiltonian

$$H=\frac{1}{2m}p^2+\frac{1}{2}m\omega ^2q^2$$
(1)
and we consider a time operator $`\chi `$, which is an operator conjugate to the Hamiltonian in that space:

$$[\chi ,H]=i\mathrm{}$$
(2)
where $`\chi =\chi (q,p)`$ is a function of the operators$`(q,p)`$, but is not an explicit function of the time $`t`$. For an eigenstate of $`H`$, the( continuous) eigenvalue of $`\chi `$ is an angle variable whose period is energy-dependent and whose domain is $`(\mathrm{},\mathrm{})`$. It is convenient to introduce a phase operator $`\varphi `$ by

$$\omega \chi =\frac{\pi }{2}\mathit{1}\varphi +G(H)$$
(3)
where $`G(H)`$ is some function of $`H`$. For an eigenstate of $`H`$, the (continuous) eigenvalue of $`\varphi `$, which we denote by $`\phi `$ and call the phase variable, is also an angle whose domain is $`(\mathrm{},\mathrm{})`$; but its period is the constant $`2\pi `$.

Aharanov and Bohm  discussed the time operator for the case of a free particle. A time operator for the harmonic oscillator was derived by Bender and Dunne - and was re-derived using similar considerations by Lewis, Lawrence, and Harris , who were motivated by its application as a starting point of a perturbation theory for the construction of invariant quantum operators

Susskind and Glogower  showed that a time operator for the harmonic oscillator does not exist in the Hilbert space of energy eigenstates, and they proposed a pair of alternate operators for use in this Hilbert space. Cannata and Ferrari  also proposed an alternate operator fro use in the Hilbert space. In a Comment on the paper by Lewis et al, Smith and Vaccaro  disputed the existence of a time operator for the harmonic oscillator in any space by pointing at an apparent inconsistency. They incorrectly asserted that, because the eigenvalue spectrum of $`H`$ is discrete, the commutation relation (2) does not have a representation in the Hilbert space of the eigenstates of $`H`$. Although Lewis et al denied the validity of this apparent inconsistency in a Reply , they did not address the issue in detail. The existence of time operators that have been derived , together with the apparent inconsistency , constitutes an apparent paradox. This apparent inconsistency is addressed in this letter. It is removed by correctly taking into account the non-hermeticity of the Hamiltonian in a space that contains the phase variable itself in addition to periodic functions of the phase variable. As a result, the commutation relation (2) is seen to have a representation in the Hilbert space, and the apparent paradox is thus resolved. The present authors believe that a time operator for the harmonic oscillator can also be derived simply and directly from the classical time function by adapting an Moyal quantization . This derivation is left for the future publication.

The apparent inconsistency raised by Smith and Vaccarots1 is based on the assertion that the relation

$$<m|[B,H]|n>=<m|B|n>(nm)\mathrm{}\omega $$
(4)
holds for any operator $`B`$, where

$$H|n>=(n+\frac{1}{2})\mathrm{}\omega |n>$$
(5)
Relation (4) is correct if operation by $`B`$ on any eigenstate of $`H`$ generates a state that lies in the space spanned by the eigenstates of $`h`$. However $`B`$ does not have that property if it is conjugate to $`H`$, as in the case with the time operator $`\chi `$. Then, relation (4) is incorrect. This can be understood clearly by working in the phase variable representation, i.e., in the representation of the eigenstates of $`\phi `$; in that representation both $`\phi `$ and $`\chi `$ are diagonal. The action of $`\phi `$ on one of its eigenstates is simply multiplication by the eigenvalue; and according to (3), the action of $`\chi `$ on an eigenstate of $`\phi `$ is multiplication by a linear function of $`\phi `$. Both $`\phi `$ and $`\chi `$ are hermitian in the full space of the eigenstates of $`\phi `$ and in every subspace thereof. The reason why (4) is incorrect when applied to $`<m|[\chi ,H]|n>`$ is that $`H`$ is not hermitian in the full space of the eigenstates of $`\phi `$. The physical boundary condition on wave functions is that they be periodic functions of $`\phi `$. In a space of periodic functions of $`\phi `$, the Hamiltonian $`H`$ is hermitian. However in a space that includes aperiodic functions of $`\phi `$, like variable $`\phi `$ itself , $`H`$ is not hermitian. We evaluate $`<m|[\chi ,H]|n>`$ as follows

$$\begin{array}{c}<m|[\chi ,H]|n>=<m|\chi H|n><m|H\chi |n>\\ =(n+\frac{1}{2})\mathrm{}\omega <m|\chi |n><m|H\chi |n>,\end{array}$$
(6)
where $`<m|[\chi ,H]|n>`$ has been evaluated using (5). If $`H`$ were hermitian then also the last term in (6) could be evaluated directly using (5):

$$\begin{array}{c}<m|H\chi |n>=<n|\chi ^{}H^{}|n>^{}=(m+\frac{1}{2})\mathrm{}\omega <n|\chi ^{}|m>^{}=\\ (m+\frac{1}{2})\mathrm{}\omega <m|\chi |n>\end{array}$$
(7)
This would give an evaluation of $`<m|[\chi ,H]|n>`$ in agreement with (4). However, $`H`$ cannot be considered hermitian in the evaluation of $`<m|[\chi ,H]|n>`$ because the presence of $`\chi `$ in $`<m|[\chi ,H]|n>`$ leads to aperiodic functions of $`\phi `$. Therefore relation (4) is invalid in this case. The is consistent with the discussion in the introductory section of a paper by Cannata and Ferrari , where it is demonstrated that relation (4) holds for any operator $`B`$ whose action is confined to the space of eigenstates of $`H`$, i.e., for any $`B`$ that is representable in the form $`B=_{k,l}\beta _{kl}|k><l|`$. Because the operators $`\phi `$ and $`\chi `$ introduce aperiodic functions of $`\phi `$ ,they are not representable in this form. A correct evaluation of $`<m|[\chi ,H]|n>`$ yields the result that corresponds to the canonical commutation relation (2). This is subject of Sec.2.

## 2 Representation of the commutation relation

We now proceed to evaluate the matrix elements of $`<m|[\chi ,H]|n>`$ directly in terms of the phase variable representation. In terms of $`\phi `$, the commutation relation (2) is

$$[\phi ,H]=i\mathrm{}\omega $$
(8)
and the phase variable representation of $`H`$ implied by (8) is

$$H=i\mathrm{}\omega \frac{d}{d\phi }+f(\phi )$$
(9)
where $`f(\phi )`$ is some function of $`\phi `$. The normalized periodic eigenfunctions of $`H`$ are

$$<\phi |n>=\frac{1}{\sqrt{2\pi }}e^{in\phi }$$
(10)
whose eigenvalues are $`n\mathrm{}\omega +f(\phi )`$. In order that the eigenvalues of $`H`$ be correctly given as $`(n+\frac{1}{2}\mathrm{}\omega )`$, we must take $`f(\phi )=\mathrm{}\omega /2`$. Thus $`H`$ is given

$$H=i\mathrm{}\omega \frac{d}{d\phi }+\frac{\mathrm{}\omega }{2}$$
(11)
We now express the matrix elements of the commutator as

$$\begin{array}{c}<m|[\chi ,H]|n>=\frac{1}{\omega }<m|[\varphi ,H]|n>\\ =(n+\frac{1}{2})\mathrm{}<m|\varphi |n>+\frac{1}{\omega }<m|H\varphi |n>\end{array}$$
(12)
and evaluate the last term as follows

$$\begin{array}{c}\frac{1}{\omega }<m|H\varphi |n>=_0^{2\pi }<\phi |m>^{}H\phi <\phi |n>𝑑\phi \hfill \\ =\frac{\mathrm{}}{2}<m|\varphi |n>+i\mathrm{}_0^{2\pi }<\phi |m>^{}\frac{d}{d\phi }<\phi |n>𝑑\phi \hfill \\ =\frac{\mathrm{}}{2}<m|\varphi |n>+i\mathrm{}_0^{2\pi }\frac{d}{d\phi }<\phi |m>^{}\phi <\phi |n>𝑑\phi \hfill \\ i\mathrm{}_0^{2\pi }\phi <\phi |n>\frac{d}{d\phi }<\phi |m>^{}𝑑\phi =i\mathrm{}<\phi |m>^{}\phi <\phi |n>|_0^{2\pi }\hfill \\ +(m+\frac{1}{2})\mathrm{}<m|\varphi |n>=i\mathrm{}+(m+\frac{1}{2})\mathrm{}<m|\varphi |n>\hfill \end{array}$$
(13)
Using

$$\begin{array}{c}<m|\varphi |n>=_0^{2\pi }<\phi |m>^{}\phi <\phi |n>𝑑\phi =\frac{1}{mn},ifmn\\ =0,ifm=n\end{array}$$
(14)
we finally obtain for the matrix elements of the commutator

$$<m|[\chi ,H]|n>=i\mathrm{}+(mn)\mathrm{}<m|\varphi |n>=i\mathrm{}\delta _{mn}$$
(15)
Thus, the commutation relation (2) has a correct representation in the space of eigenstates of $`H`$, despite of the fact that $`\chi `$ cannot be represented in that space.


text:
# Black holes and a scalar field in an expanding universe

## 1 Introduction

Recently strongly motivated by the supernovae observations , general agreement has been reached that our universe is going to turn to accelerated expansion in the context of Friedmann universe model. As one possible model, a scalar field can be introduced in the Friedmann universe . Such a scalar field plays the role as the quintessence of accelerated expansion. The investigation of scalar field for an expanding universe is an important current issue.

On the other hand, our universe consists of celestial objects like galaxies, stars and possibly black holes. In order to draw a more precise picture of our universe than that of the Friedmann model, it is useful to investigate the effects of the inhomogeneities on the expansion of the universe 333 This point of view may lead us to study black holes in a non-flat background. So far only a few works related to this issue have been done numerically or in a toy model  . In ref., the expansion law of inhomogeneous universe including many black holes has already been investigated. We call the model of an inhomogeneous universe used in ref. the cell lattice universe. It expresses an averaged space time of inhomogeneous universe consisting of many identical black holes, which is constructed as follows: considering a regularly tessellated homogeneous and isotropic universe, each cell of the tessellation is replaced by a spherically symmetric black hole. Here it is assumed that the averaged values of any quantities are defined by the values on the junction surface where a homogeneous universe is connected with a spherical black hole. The behavior of the averaged scale factor is determined by the junction condition of the cell lattice universe, and reproduces the same expansion law that the dust-dominated Friedmann model obeys, as has already been shown in ref.. We believe that the cell lattice universe model gives us the well-defined averaged quantities of an inhomogeneous universe, and that this model coincides with the matter-dominated universe.

In this paper, we extend the cell lattice universe model in order to include a scalar field. Motivated by ref., we assume the averaged scalar field and scale factor are given by the values on the junction surface. Furthermore, we require the averaged quantities to agree with those of the dust-dominated Friedmann model including a scalar field. Here the undetermined quantity is the form of the metric of a spherical black hole. With these assumptions we are interested in the question: how should the black hole be affected by the expansion of the universe and the existence of a scalar field in order to retain the consistency of the cell lattice universe with the dust-dominated Friedmann model. That is, in contrast to ref. where the expansion law of the cell lattice universe has been investigated, we turn our interest to the black holes in the cell lattice univese. It is well known that the system of Einstein gravity with scalar fields can be transformed to the framework of scalar-tensor gravity. For simplicity, we consider the scalar field to be a massless one and carry out all of the calculations within the framework of the scalar-tensor gravity, especially of the Brans-Dicke gravity. In constructing the cell lattice universe in Brans-Dicke gravity in the way mentioned at the begining of this paragraph, we define the mass of the black hole by making use of the undetermined form of the metric. If the time dependence of our mass is adiabatic and the universe is expanding, the mass is equivalent to the ADM mass of a Schwarzschild black hole without regard to the details of the expansion law of the universe. The junction condition in replacing a cell by a black hole gives us an equation relating the mass to the averaged scalar field and scale factor. Through this equation we can investigate whether and how the mass evolves with time, which is an effect of the expansion of the universe and the existence of the scalar field on black holes.

In sections 2 and 3, Friedmann universe in Brans-Dicke gravity and the method of constructing the cell lattice universe are reviewed, respectively. Section 4 is devoted to the construction of the cell lattice universe in Brans-Dicke gravity and to analyses of the mass which is defined in the section 4. Finally we give a summary and discussion in section 5. Throughout this paper, we set $`c=1`$.

## 2 Friedmann Universe in Brans-Dicke Gravity

The action of Brans-Dicke gravity is

$`S={\displaystyle d^4x\sqrt{g}\left[\phi R\frac{\omega }{\phi }(\phi )^2+_m\right]},`$ (1)
where $`\omega `$ is the constant parameter of this theory, $`\phi `$ is the scalar field coupling to gravity and $`_m`$ is the matter Lagrangian which does not include $`\phi `$. The effective Newton “constant”, $`G_{eff}`$, is related to the scalar field as $`\phi =(16\pi G_{eff})^1`$. The field equations derived from this action are

$`R_{\mu \nu }{\displaystyle \frac{1}{2}}g_{\mu \nu }R`$ $`=`$ $`{\displaystyle \frac{\omega }{\phi ^2}}\left[(_\mu \phi )(_\nu \phi ){\displaystyle \frac{1}{2}}g_{\mu \nu }(\phi )^2\right]{\displaystyle \frac{1}{\phi }}\left[g_{\mu \nu }\mathrm{}_\mu _\nu \right]\phi +{\displaystyle \frac{1}{2\phi }}T_{\mu \nu },`$
$`\mathrm{}\phi `$ $`=`$ $`{\displaystyle \frac{1}{4\omega +6}}T_\mu ^\mu ,`$ (2)
where $`T_{\mu \nu }`$ is the energy-momentum tensor derived from $`_m`$, which is automatically divergenceless: $`_\mu T^{\mu \nu }=0`$.

We consider the Friedmann universe in Brans-Dicke gravity. The metric is spatially homogeneous and isotropic,

$`ds^2=dt^2+a(t)^2\left[{\displaystyle \frac{dr^2}{1kr^2}}+r^2(d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2)\right],`$ (3)
where $`k=1`$,$`0`$,$`1`$ and $`a(t)`$ is the scale factor. The matter in the universe is of perfect fluid type,

$`T_{\mu \nu }=ϵu_\mu u_\nu +p(g_{\mu \nu }+u_\mu u_\nu ),`$ (4)
where $`u_\mu `$ is the 4-velocity of comoving observer, $`ϵ`$ is the energy density and $`p`$ is the pressure. Then the field eqs.(2) give three independent equations,

$`H^2+{\displaystyle \frac{k}{a^2}}`$ $`=`$ $`{\displaystyle \frac{ϵ}{6\phi }}H{\displaystyle \frac{\dot{\phi }}{\phi }}+{\displaystyle \frac{\omega }{6}}\left({\displaystyle \frac{\dot{\phi }}{\phi }}\right)^2,`$ (5)
$`\ddot{\phi }+3H\dot{\phi }`$ $`=`$ $`{\displaystyle \frac{1}{4\omega +6}}\left(ϵ3p\right),`$ (6)
$`\dot{ϵ}`$ $`=`$ $`3H\left(ϵ+p\right).`$ (7)
where $`H=\dot{a}/a`$ is the Hubble parameter. This system can be solved provided an equation of state is specified.

For a dust-dominated universe the equation of state is $`p=0`$, with which eq.(7) gives $`ϵ=ϵ_0/a^3`$, where $`ϵ_0`$ is an integration constant. Then we obtain from eq.(6)

$`\phi ={\displaystyle \frac{ϵ_0}{4\omega +6}}{\displaystyle ^t}𝑑t{\displaystyle \frac{t+t_0}{a^3}},`$ (8)
where $`t_0`$ is an integration constant. Eqs.(5) and (8) determine the time evolution of the scale factor and the scalar field of the dust-dominated Friedmann universe in Brans-Dicke gravity.

For a radiation-dominated universe the equation of state is $`p=ϵ/3`$, then eq.(7) gives $`ϵ=ϵ_0/a^4`$, where $`ϵ_0`$ is an integration constant. Using eq.(6) we obtain

$`\phi ={\displaystyle ^t}𝑑t{\displaystyle \frac{q}{a^3}},`$ (9)
where $`q`$ is an integration constant. Eqs.(5) and (9) determine the time evolution of the scale factor and the scalar field of the radiation-dominated Friedmann universe in Brans-Dicke gravity.

For the use of section 4, let us review the frame transformation. By a suitable conformal transformation of the metric $`g_{\mu \nu }`$ to the other one $`\stackrel{~}{g}_{\mu \nu }`$, the action $`S`$ of eq.(1) can be treated as Einstein gravity with a massless scalar field which does not couple to gravity. Such a transformation is given by

$`g_{\mu \nu }=\mathrm{exp}\left[{\displaystyle \frac{\sigma }{2\omega +3}}\right]\stackrel{~}{g}_{\mu \nu },`$ (10)
where $`\sigma `$ is a scalar field defined by $`\phi =(16\pi G_0)^1\mathrm{exp}[\sigma /(2\omega +3)]`$. Here $`G_0`$ is the ordinary Newton constant in the Einstein gravity. Then the action $`S`$ becomes

$`S={\displaystyle d^4x\sqrt{\stackrel{~}{g}}\left[\frac{1}{16\pi G_0}\left(\stackrel{~}{R}(\stackrel{~}{}\sigma )^2\right)+e^{2\sigma /(2\omega +3)}_m\right]}.`$ (11)
The action expressed by $`\stackrel{~}{g}_{\mu \nu }`$ is said to be in the Einstein frame, while the action of eq.(1) is in the original frame. We add tilde to quantities in the Einstein frame. It is obvious by this frame transformation that the Brans-Dicke gravity in the original frame becomes the Einstein gravity in the limit $`\omega \mathrm{}`$, that is, $`g_{\mu \nu }=\stackrel{~}{g}_{\mu \nu }`$.

It has already known analytically that, in the dust-dominated and flat universe case, $`p=0`$ and $`k=0`$, only the choice that $`t_0=0`$ lets the scale factor $`a(t)`$ in Brans-Dicke gravity go over smoothly to that of Einstein gravity in the limit $`\omega \mathrm{}`$ . With $`t_0=0`$ we can find the relation

$`{\displaystyle \frac{8\omega +12}{3\omega +4}}={\displaystyle \frac{ϵ_0t^2}{\phi a^3}}.`$ (12)
The scale factor and scalar field in this case can be easily obtained to be

$`a(t)`$ $`=`$ $`a_0t^{(2\omega +2)/(3\omega +4)},`$
$`\phi (t)`$ $`=`$ $`\phi _0t^{2/(3\omega +4)},`$ (13)
where $`a_0`$ and $`\phi _0`$ are constants. In the other curvature cases $`k=\pm 1`$ with a restriction $`t_0=0`$, eq.(5) indicates that the universe is dominated at early times by terms other than the curvature one $`k/a^2`$. That is, the curved universe is effectively flat at early times.

## 3 Cell Lattice Universe

### 3.1 Strategy to construct the cell lattice universe and Averaged quantities

It is a purely geometric problem to construct the cell lattice universe which expresses an averaged inhomogeneous expanding universe including many identical spherically symmetric objects . The procedure of the construction consists of two steps as follows: 1st step, regularly tessellate the spherical, flat or hyperbolic spatial section, $`\mathrm{\Sigma }_t`$, of a homogeneous and isotropic universe by a regular polyhedron  . 2nd step, replace each cell made of a polyhedron by a metric of the spherically symmetric object. For the junction condition in the 2nd step we consider the radially directional differential of any great circle’s circumference on the two dimensional spherical junction surface in $`\mathrm{\Sigma }_t`$, which can be expressed as:

$`D={\displaystyle \frac{d(\text{Circumference of the junction surface in }\mathrm{\Sigma }_t)}{d(\text{Physical radial distance along }\mathrm{\Sigma }_t)}}.`$ (14)
The junction condition is to connect such a differential measured in the homogeneous and isotropic universe, $`D_{univ}`$, with that measured in the metric of a spherically symmetric object, $`D_{BH}`$. That is the junction condition is given by

$`D_{univ}=D_{BH}.`$ (15)
Furthermore, in the 2nd step, there are some regions where spherically symmetric metrics overlap, and other regions where the metrics do not cover. Here we assume that the averaged values of any quantities are defined by the values on the spherical junction surface. For example, the averaged scale factor of the cell lattice universe is the scale factor on the spherical junction surface. The averaged cosmological time is set by the proper time of the observer staying on the spherical junction surface. We consider the cell lattice universe to be an averaged spacetime in such a sense.

### 3.2 Construction of the cell lattice universe

The metric of a homogeneous and isotropic universe can be given in the same form as eq.(3). The 1st step in constructing the cell lattice universe is to determine the radius of each spherical cell measured in this metric   by requiring that the volume of a spherical cell coincides with that of a regular polyhedron. There are some regular polyhedra which are candidates for the regular tessellation for each curvature case. The types of polyhedra are listed later in this section after the explanation of the junction condition. We denote the comoving radius of the spherical cell by $`r_c`$. Since the comoving number density of cells is constant, the cell expands comovingly. Therefore $`r_c=const.`$, and the averaged cosmological time is identical to the time coordinate of the metric (3), $`t_c`$, on the junction surface $`r=r_c`$. Hereafter we denote the physical radius of a cell by $`l_c(t_c)=a(t_c)r_c`$.

In proceeding to the 2nd step, we define precisely the “radial direction” in the differential of the great circle on the junction surface. In choosing the parameter along the radial direction, $`\chi `$, into a coordinate system such as $`(t,\chi ,\theta ,\varphi )`$, we define the coordinate parameter $`\chi `$ as follows: with the conditions that $`t`$, $`\theta `$ and $`\varphi `$ are constant, the line element along the radial direction is given by $`ds^2=a^2d\chi ^2`$, that is to say, $`\chi `$ is the “comoving radial distance”. For convenience we transform the coordinate $`(t,r,\theta ,\varphi )`$ to $`(t,\chi ,\theta ,\varphi )`$ by $`r=f_k(\chi )`$, where $`f_k(\chi )=\mathrm{sin}\chi `$, $`\chi `$ and $`\mathrm{sinh}\chi `$ for $`k=1`$, $`0`$ and $`1`$, respectively:

$`ds^2=dt^2+a(t)^2\left[d\chi ^2+f_k(\chi )(d\theta ^2+\mathrm{sin}\theta d\varphi ^2)\right].`$ (16)
In denoting $`r_c=f_k(\chi _c)`$, this $`\chi _c`$ should be constant because $`\chi _c`$ is the comoving coordinate on the junction surface. The circumference of great circle, $`C_{gc}`$, is calculated in this metric to be

$`C_{gc}={\displaystyle _{\mathrm{\Sigma }_t,\chi =\chi _c,\theta =\pi /2}}\sqrt{g_{\varphi \varphi }}𝑑\varphi =2\pi af_k(\chi _c)(=2\pi l_c).`$ (17)
Then we obtain the radially directional differential of the great circle in the metric (16), $`D_{univ}`$, to be

$`D_{univ}={\displaystyle \frac{dC_{gc}}{d(a\chi )}}|_{\mathrm{\Sigma }_t,\chi _c}=2\pi {\displaystyle \frac{df_k(\chi _c)}{d\chi _c}}.`$ (18)
The metric of a spherically symmetric object can be given by

$`ds^2=A(T,R)dT^2+B(T,R)dR^2+R^2(d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2),`$ (19)
where $`A(T,R)`$ and $`B(T,R)`$ are the arbitrary independent functions of $`T`$ and $`R`$. The radial coordinate $`R`$ is the physical radial length. Hereafter we assume that the radius of a spherical object is smaller than that of a cell, therefore the metric (19) is the exterior spacetime of the object. In the 2nd step of constructing the cell lattice universe, we should pay attention to the fact that the coordinates $`(T,R)`$ of spherical object are not necessarily identical to $`(t,r)`$ or $`(t,\chi )`$ of a homogeneous and isotropic universe. When the origins of the coordinates of both metrics coincide, the radii of the junction surface in both coordinate systems coincide. That is, we have $`R(t_c)=l_c(t_c)`$ on the junction surface. The time $`T`$ of the metric (19) on the junction surface can be expressed by a function of $`t_c`$ and $`r_c`$ as $`T_c=T(t_c,r_c)`$, which should reflect the junction condition and the time evolution of the universe. The radially directional differential of the great circle in the metric (19) should be calculated on $`\mathrm{\Sigma }_t`$. The circumference of the great circle, $`C_{gc}`$, is given by eq.(17). At the junction surface, the vector which is orthogonal to $`\mathrm{\Sigma }_t`$ can be given by $`e^\mu =(dT_c,dl_c,0,0)`$ in the coordinate of metric (19), then the vector which is orthogonal to the junction surface and parallel to $`\mathrm{\Sigma }_t`$ can be obtained to be $`v^\mu =n(g^{TT}dl_c,g^{RR}dT_c,0,0)`$, where $`n`$ is the normalization constant. This gives us the difference along the radial direction as $`d(a\chi )^2=g_{\mu \nu }v^\mu v^\nu `$. Here we determine the normalization, $`n`$, by requiring that the difference of the circumference be related to the vector $`v^\mu `$ as $`dC_{gc}=2\pi dl_c=2\pi v^R`$. For the averaged cosmological time, we have the relation $`dt_c^{\mathrm{\hspace{0.17em}\hspace{0.17em}2}}=g_{\mu \nu }e^\mu e^\nu `$. Then we can calculate the radially directional differential of the great circle in the metric (19), which we denote as $`D_{BH}`$, as follows:

$`D_{BH}={\displaystyle \frac{2\pi dl_c}{d(a\chi _c)}}=\mathrm{\hspace{0.17em}2}\pi {\displaystyle \frac{g^{RR}dT_c}{\sqrt{g_{\mu \nu }v^\mu v^\nu }}}=\mathrm{\hspace{0.17em}2}\pi \sqrt{{\displaystyle \frac{g_{TT}}{g_{RR}}}}{\displaystyle \frac{dT_c}{dt_c}}`$ (20)
With the preparations done above, we proceed to calculating the junction condition, $`D_{univ}=D_{BH}`$. Substituting eq.(18) and the second expression of $`D_{BH}`$ of eq.(20) for the junction condition, we find the relation: $`(dT_c/dl_c)^2=\alpha _k^{\mathrm{\hspace{0.17em}\hspace{0.17em}2}}B/[A(\alpha _k^{\mathrm{\hspace{0.17em}\hspace{0.17em}2}}1/B)]`$, where $`\alpha _k=\mathrm{cos}\chi _c`$, $`1`$ and $`\mathrm{cosh}\chi _c`$, for $`k=+1`$, $`0`$ and $`1`$, respectively. Wiht this relation and $`D_{univ}=D_{BH}`$ with the third expression of $`D_{BH}`$ in eq.(20), we obtain the equation :

$`\left(\beta _k{\displaystyle \frac{da(t_c)}{dt_c}}\right)^2=\alpha _k^{\mathrm{\hspace{0.17em}\hspace{0.17em}2}}{\displaystyle \frac{1}{B(t_c,l_c)}},`$ (21)
where $`\beta _k=f_k(\chi _c)`$. The parameters $`\alpha _k`$ and $`\beta _k`$ are determined by the cell radius $`\chi _c`$. Here we list the types of regular polyhedra, $`\chi _c`$ and $`N`$ in the Table 1, where $`N`$ is the number of the cells included in the spatial surface $`\mathrm{\Sigma }_t`$  .

For the cell lattice universe made using Schwarzschild black holes in Einstein gravity, the averaged scale factor is given by eq.(21) with $`B^1=12G_0M/R`$, where $`G_0`$ is the ordinary Newton constant.

It can be easily checked that, without regard to the polyhedron used for the tessellation, this equation reproduces the same expansion law as for a dust-dominated Friedmann universe in Einstein gravity . In this paper we assume that, even in any of the generalized theories of gravity, the avereged values of any quantities coincide with those in the dust-dominated Friedmann model.

## 4 Time evolution of Black Holes in an Expanding Universe

In this section we consider the spherical object in the cell lattice universe to be a black hole.

### 4.1 Black hole mass and cell lattice universe

As mentioned in section 1, our interest in this paper is the effect of the expansion of the universe and a massless scalar field on black holes. We consider this issue within the framework of Brans-Dicke gravity and try to extract such effects from the cell lattice universe which is assumed to be meaningful as an approximation of multi-black hole spacetime even in Brans-Dicke gravity. This issue can be discussed by eq.(21) once the function $`B(T,R)`$ is specified.

To specify the form of $`B(T,R)`$, we should note here that the uniqueness theorem  of the black holes in Brans-Dicke gravity has already been established with the asymptotically flat condition, not in the expanding universe . That is, the spherically symmetric non-rotating black hole in the asymptotically flat spacetime is of the Schwarzschild black hole even in Brans-Dicke gravity. Furthermore, it is important for specifying $`B(T,R)`$ to note the issue peculiar to the scalar-tensor gravity: gravitational memory effect. This issue was originally recognized in ref. and expressed as follows: when a black hole is formed in the expanding universe, does the scalar field on the event horizon keep the value at the black hole formation time during cosmological evolution? If the scalar field keeps its original value, it is said that the black hole has a gravitational memory effect, if not, the black hole does not have the memory 444 In ref., two possible scenarios of cosmology are discussed: the scenario that primordial black holes have a gravitational memory effect, and the other scenario where primordial black holes don’t have such effects.. Concerning this issue, ref. treats the perturbation of a stationary black hole with the boundary condition that the scalar field is proportional to the time in the region spatially far from the black hole. With this model, it is suggested in ref. that the gravitational memory effect does not occur in scalar-tensor gravities. That is, the scalar field on the event horizon seems to have a time dependence along with cosmological evolution. Then we adopt a naively acceptable ansatz as follows: the spherically symmetric black hole composing the cell lattice universe in Brans-Dicke gravity is of the same form as the Schwarzschild black hole except for the one point that its radius depends on the cosmological time, $`R_g(t_c)`$. With this ansatz we define the mass of the black hole in the cell lattice universe, $`M`$, in Brans-Dicke gravity 555 There have already been some discussions concerning the spatial profile of the scalar field and the mass of a black hole and boson star in the scalar-tensor gravity, see  for examples.  as,

$`M(t_c)8\pi \phi (t_c)R_g(t_c).`$ (22)
That is, we assume the black hole to be of “Schwarzschild-type” given by eq.(19) on the junction surface with

$`A(t_c,l_c)=B(t_c,l_c)^1=1{\displaystyle \frac{M(t_c)}{8\pi \phi (t_c)l_c(t_c)}}.`$ (23)
This mass, $`M`$, certainly has the dimensions of mass, but it is not obvious whether this definition can be consistent with the ADM mass defined by the generator of Killing time translation in a stationary space time. Without regard to the details of the expansion law of the averaged scale factor, provided that $`M`$ changes adiabatically over a cosmological time scale and that the black hole horizon size is much smaller than the cosmological one, we can consider that the black hole exists in a local asymptotically flat region as for an ordinary Schwarzschild black hole and that the cosmological time $`t_c`$ is equivalent to the time $`T`$ of the metric (19) within the scale of the black hole. In such a case our mass $`M`$ can be effectively treated as the ADM mass defined by the generator of the local asymptotically timelike Killing vector, where the Killing time is $`T`$ of the matric (19).

With specifying $`B(T,R)`$ as above, the junction condition (21) gives

$`M(t_c)=8\pi \beta _k^{\mathrm{\hspace{0.17em}\hspace{0.17em}3}}\phi (t_c)a(t_c)\left[\dot{a}(t_c)^2+k\right].`$ (24)
Clearly, since $`\beta _k`$ is only related to the normalization of the mass, the behavior of $`M(t_c)`$ is not affected by the value of $`\chi _c`$, that is, it is not affeced by the choice of the polyhedron for the tessellation in the 1st step of the construction of the cell lattice universe. Here, motivated by the comment in the final paragraph in the previous section, we assume the cell lattice universe to be a well-defined averaged model of an inhomogeneous universe, reproducing the expansion law of a dust-dominated Friedmann model even in Brans-Dicke gravity. Because of this assumption, the averaged scale factor, $`a(t_c)`$, and the scalar field, $`\phi (t_c)`$, are given by the dust-dominated Friedmann universe in Brans-Dicke gravity. The time evolution of $`M(t_c)`$ can be investigated using eqs.(5), (8) and (24). It is not a priori obvious whether or not the mass $`M`$ depends adiabatically on the cosmological time.

As has already been pointed out in ref., provided our mass $`M`$ can be considered as the ADM mass in a local asymptotically flat region, we can discuss whether the mass in the Einstein frame changes with time. Attaching a tilde to the quantity in the Einstein frame, we find the relation $`M/\stackrel{~}{M}=d\stackrel{~}{T}/dT`$, because the ADM mass is given by the generator of the local asymptotic Killing vector, which is normalized at the boundary of the local asymptotically flat region and scales inversely with the time $`T`$. Since the Killing time translation is given by the line element along integral curves of the timelike Killing vector, we obtain the relation at the boundary of the local asymptotically flat region: $`dT^2=\mathrm{exp}[\sigma /(2\omega +3)]d\stackrel{~}{T}^2`$, where $`dT`$ and $`d\stackrel{~}{T}`$ express the Killing time translations in the Brans-Dicke and the Einstein frames, respectively, and the frame transformation given by eq.(10) is used. Therefore we obtain the relation of the masses

$`M=\mathrm{exp}\left[{\displaystyle \frac{\sigma }{4\omega +6}}\right]\stackrel{~}{M}\sqrt{\phi }\stackrel{~}{M}.`$ (25)
If $`\stackrel{~}{M}`$ is constant, $`M`$ is proportional to $`\sqrt{\phi }`$. In the case that the $`M`$ defined by ep.(24) can be treated as the ADM mass, by comparing the $`M`$ with $`\sqrt{\phi }`$, we can obtain a suggestion concerning the time dependence of $`\stackrel{~}{M}`$ which is the black hole mass in an expanding universe in Einstein gravity with the presence of a scalar field.

### 4.2 Analysis of the black hole mass

In the numerical calculations shown here, our interests are the time dependence of $`M`$, not the absolute value of it. That is, we investigate whether or not $`M`$ changes adiabatically on a cosmological time scale and whether $`M`$ is proportional to $`\sqrt{\phi }`$ or not, since these issues reveal that our black hole mass in the Brans-Dicke frame is equivalent to the ADM mass and that the black hole mass in the Einstein frame has a time dependence. Hereafter we denote the averaged cosmological time $`t_c`$ simply by $`t`$.

In flat case $`k=0`$, it is easily found from eqs.(13) and (24) that the mass $`M`$ is constant, and that the radius $`R_g(t)1/\phi (t)`$ decreases in a decelerating fashion, while the averaged scale factor continues to expand. Therefore, the mass can be effectively treated as ADM mass after the universe expands enough to be much larger than the black hole. Clearly $`M\propto ̸\sqrt{\phi }`$. This means that the mass in the Einstein frame has a time dependence.

In two cases $`k=\pm 1`$, we calculate the mass $`M(t)`$ numerically using Mathematica from $`t=0.5`$ to $`t=50000`$ for the open case $`k=1`$, and to $`t=1000`$ for the closed case $`k=1`$. We attach the subscript “$`i`$” to the quantities evaluated at $`t=0.5`$, and the subscript “$`f`$” at $`t=50000`$ for $`k=1`$ and at $`t=1000`$ for $`k=1`$. These calculation time intervals of calculation should correspond to the dust-dominated era of our universe, since the cell lattice universe is assumed to be the averaged spacetime of a dust-dominated inhomogeneous universe as discussed in the previous subsection. We regard the end time of the calculation as our present time, because, as shown in the following, the rate of change of the mass $`M`$ in time evolution is sufficiently small in comparison with the averaged Hubble parameter at the end time for both cases $`k=\pm 1`$, that is $`M`$ evolves adiabatically. Here the averaged Hubble parameter $`H=\dot{a}/a`$ is defined by using the averaged scale factor. As mentioned at the end of section 2, the curved Friedmann universe can be effectively treated as being flat in the early stage of the time evolution. Therefore, we set the integration constant $`t_0`$ of eq.(8) to zero ($`t_0=0`$) in order to reproduce the expansion law of the Friedmann model in Einstein gravity as $`\omega \mathrm{}`$. The initial values of the averaged scale factor $`a_i`$ and scalar field $`\phi _i`$ are related through eq.(12). Our choices of the other parameters are as follows: the initial value of scale factor, $`a_i=10`$. The Brans-Dicke parameter $`\omega =500`$, which is the experimental lower bound of $`\omega `$. Another integration constant $`ϵ_0=100`$. The size of one cell in the cell lattice universe $`\chi _c=0.9747`$ (Icosahedron) and $`0.1993`$ (Tetrahedron) for $`k=1`$ and $`1`$, respectively  . The value of $`\chi _c`$ does not affect the behavior of $`M`$ as shown in eq.(24).

For the open case $`k=1`$, the averaged Hubble time at present $`t=50000`$ is $`1/H_f=53931`$. The table 2 shows some values of the averaged cosmological red shift which is defined by $`z(t)=a_f/a(t)1`$ with the averaged scale factor $`a`$. Fig.2 is for the rate of change of $`M(t)`$ normalized by the averaged Hubble parameter: $`H_m(t)=(\dot{M}/M)/H`$. The $`H_m`$ is of negative and asymptotes to zero for sufficiently late stage of cosmological evolution $`t>10000`$ where $`z<3.170`$. The rate of change, for example at $`t=2`$, takes the value $`H_m(2)=0.0001469`$. This indicates that the time evolution of $`M`$ is less significant than that of the averaged scale factor at least for $`t>2(z<1840)`$, where $`H_m<O(10^4)`$. That is, $`M`$ is adiabatic in this epoch. The radius of black hole $`R_g`$ decreases in a decelerating fashion with time, from $`R_{gi}=5210`$ to $`R_{gf}=5142`$, while the averaged scale factor continues to expand for all time. Therefore, the mass $`M`$ can be effectively treated as being ADM mass after the universe expands enough to be much larger than the black hole. Fig.3 is the plot of $`M`$ with $`SqM`$ which is defined by $`SqM(t)=M_f\sqrt{\phi (t)/\phi _f}`$. The square root of the averaged scalar field $`SqM`$ at the initial time is $`SqM_i=609.4`$, and increases in a decelerating fashion for all time. The rate of change of mass $`H_m`$ becomes zero at $`t=14.92`$. The mass $`M`$ increases from $`M_i=613.6`$ to $`M(14.29)=613.9`$, and decreases in a decelerating fashion after the time $`t=14.92`$ where $`z(14.92)=17.25`$. Fig.3 indicates that mass $`M`$ does not coincides with $`SqM`$, and that the mass in the Einstein frame has a time dependence.

For the closed case $`k=1`$, the averaged Hubble time at present $`t=1000`$ is $`1/H_f=1743`$. Table 3 includes some values of the averaged cosmological red shift $`z(t)=a_f/a(t)1`$. Fig.4 shows $`H_m`$, the rate of change of $`M`$ normalized by the averaged Hubble parameter. The rate of change $`H_m`$ decreases rapidly at early times, and turns to increase around $`t=20`$ where $`z=11.09`$. At present, the rate of change is $`H_{mf}=0.0006698`$. The increase of $`H_m`$ continues up to the turning time of the averaged scale factor from expansion to contraction, where the turning time is $`t=5483`$. The rate of change, for example at $`t=2`$, is $`H_m(2)=0.0001565`$. This indicates that the time evolution of the mass $`M`$ is adiabatic within the epoch from $`t=2`$ to the present $`t=1000`$, where $`H_m<O(10^4)`$. The black hole radius $`R_g`$ decreases in a decelerating fashion with time, from $`R_{gi}=27.54`$ to $`R_{gf}=27.29`$, while the averaged scale factor lasts to expand untill the turning time $`t=5483`$. Therefore, the mass $`M`$ can be effectively treated as ADM mass after the universe expands to a size much larger than the black hole. Fig.5 shows a plot of $`M`$ and $`SqM`$. The square root of the averaged scalar field $`SqM`$ increases in a decelerating fashion for all time from the initial value $`SqM_i=3.227`$. The mass $`M`$ takes the initial value $`M_i=3.224`$, and also increases in a decelerating fashion for all time. Fig.5 indicates that $`M`$ and $`SqM`$ do not coincide with each other, and that the mass in the Einstein frame has a time dependence.

## 5 Summary and Discussion

In order to investigate the effects of expansion of the universe and a scalar field on the celestial objects which compose the inhomogeneity of the universe, we constructed the cell lattice universe in Brans-Dicke gravity. The assumptions in constructing the universe were that the black holes included in the universe are of “Schwarzschild-type”, eqs.(19) and (23), and that the averaged scale factor and scalar field are given by those of the Friedmann universe. Furthermore, we defined the mass $`M`$ by eq.(22) which is equivalent to the ordinary ADM mass in the case of adiabatic time evolution of $`M`$. The junction condition of the cell lattice universe gives eq.(24) to calculate the time evolution of the black hole mass $`M`$. It is not a priori known whether and how the mass depends on time.

As a result, it turns out that $`M`$ behaves in a qualitatively different way with respect to the value of $`k`$. The mass $`M`$ decreases in a decelerating fashion for the open case $`k=1`$, stays completely constant for the flat case $`k=0`$ and increases in a decelerating fashion for the closed case $`k=1`$. The rate of change of the mass $`H_m`$ is very small for both cases of $`k=\pm 1`$, $`H_m<O(10^4)`$ with our numerical results. This means that $`M`$ evolves adiabatically in time. Because here we assume that the averaged scale factor in eq.(24) is given by that of the Friedmann model in Brans-Dicke gravity, the adiabaticity of $`M`$ is consistent with the fact that the cell lattice universe constructed with the Schwarzschild black hole, whose mass is completely constant, reproduces the expansion law of the Friedmann universe in Einstein gravity . Furthermore, the radius $`R_g`$ decreases in a decelerating fashion for every case of $`k=\pm 1`$, $`0`$, while the averaged scale factor continues to increase. This behavior of $`R_g`$ can be easily understood for the case $`k=0`$, which is $`R_gt^{2/(3\omega +4)}`$. Therefore, it is reasonable to consider the black hole is in a local asymptotically flat region after the universe expands enough to be much larger than the black hole, consequently mass $`M`$ is equivalent to the ADM mass defined by the local asymptotic Killing time translation.

Furthermore, we can recognize using eq.(25) that the black hole mass in the Einstein frame has a time dependence for every case of $`k=\pm 1`$, $`0`$. According to the uniqueness theorem  , non-rotating and non-charged black hole in asymptotically flat spacetime is specified only by the ADM mass, which is completely constant in this theorem. The time dependence of the mass may indicate that the uniqueness theorem in an expanding universe will be broken on a cosmological time scale. It is a very interesting problem which remains to be solved.

In the numerical calculations shown in the section 4, we set $`\omega =500`$ which is the lowest experimental value. If the parameter is set as $`\omega <500`$, the absolute value of the rate of change of the mass $`M`$ tends to increase. However the total behavior of the mass $`M`$ and the radius $`R_g`$ for the case $`\omega <500`$ are the same as for the case $`\omega =500`$. Consequently, for arbitrary value of $`\omega `$, the mass can be equivalent to the ADM mass once the universe expands enough to be much larger than the black hole size.

So far we have only treated the Brans-Dicke gravity as a representative case of scalar-tensor gravities in which the parameter $`\omega `$ depends on the scalar field, and have not introduced any potential of the scalar field. When, instead of Brans-Dicke gravity, we consider a general scalar-tensor gravity with a potential of the scalar field in constructing the cell lattice universe, the averaged scale factor and scalar field may behave quite differently from those in Brans-Dicke gravity. However, the same arguments and results obtained until the previous paragraph are true of this case if the averaged quantities asymptote to those in the Brans-Dicke gravity at least in a sufficiently late stage of the expansion of the universe. Furthermore, in the case that the averaged scale factor or scalar field in a general theory of gravity always behaves differently from that in the Brans-Dicke gravity, though we need to reanalyze eq.(24) in order to know the details of the time dependence of $`M`$, however, it is natural to propose that the mass $`M`$ depends on time in this case too. That is, it seems to be quite a general point that the black hole mass in expanding universe with the presence of scalar field has some time dependence.

Let us comment on our assumption, i.e. the “Schwarzschild-type” ansatz of the black hole composing the cell lattice universe given by eqs.(19), (22) and (23). The construction of the cell lattice universe where a cell is replaced by a black hole, means that the dust-matter on homogeneous and isotropic universe in a cell is concentrated at a center point of the cell. Therefore, it is appropriate to consider the black hole spacetime replacing the cell includes only the scalar field coupling with gravity in the Brans-Dicke frame. The field equation of such a scalar field is obtained by substituting $`T_{\mu \nu }=0`$ into eqs.(2),

$`G_{\mu \nu }={\displaystyle \frac{\omega }{\phi ^2}}\left[(_\mu \phi )(_\nu \phi ){\displaystyle \frac{1}{2}}g_{\mu \nu }(\phi )^2\right]+{\displaystyle \frac{1}{\phi }}_\mu _\nu \phi ,`$ (26)
where the non-zero components of the Einstein tensor are calculated from eqs.(19), (22) and (23),

$`G_{01}`$ $`=`$ $`{\displaystyle \frac{\dot{R_g}}{R^2RR_g}},`$
$`G_{22}`$ $`=`$ $`R^3\left[{\displaystyle \frac{\ddot{R_g}}{2(RR_g)^2}}+{\displaystyle \frac{\dot{R_g}^2}{(RR_g)^3}}\right],`$ (27)
$`G_{33}`$ $`=`$ $`\mathrm{sin}^2\theta G_{22},`$
where $`\dot{R_g}=dR_g/dT`$ with the time $`T`$ of metric (19). The scalar field on the Schwarzschild-type spacetime should satisfy the above field eqs.(26) and (27).

Our ansatz of the Schwarzschild-type black hole is motivated by the uniqueness theorem of the black hole in asymptotically flat spacetime  and the indication in ref. that the gravitational memory effect seems not to occur. Conversely, in paying attention to the latter motivation, the gravitational memory effect can be discussed using eqs.(26) and (27) with the boundary condition $`\phi (T)`$ at $`R=l_c(T)`$ given by the scalar field of Friedmann model, where $`l_c(T)`$ is determined by the junction of our cell lattice universe. By investigating this system we can know whether or not the scalar field $`\phi `$ has a spatial dependence. If not, it means that the graviational memory effect does not occur in Brans-Dicke gravity, and that the indication in ref. and our discussion in this paper are supported. The model of spacetime in such an approach to the gravitational memory effect, can be considered as a modified swiss cheese universe . In the ordinary swiss cheese model, a spherically symmetric region in the Friedmann universe is replaced by an ordinary Shcwarzschild black hole of constant radius and mass, but here we replace the spherical region by a “Schwarzschild-type” black hole. This is an interesting and solvable problem.

## Acknowledgments

We would like to Thank M.Sakagami and A.Ohashi for their useful comments and discussions.