Gravitational Interaction in the Chimney Lattice Universe

We investigate the influence of the chimney topology $T\times T\times R$ of the Universe on the gravitational potential and force that are generated by point-like massive bodies. We obtain three distinct expressions for the solutions. One follows from Fourier expansion of delta functions into series using periodicity in two toroidal dimensions. The second one is the summation of solutions of the Helmholtz equation, for a source mass and its infinitely many images, which are in the form of Yukawa potentials. The third alternative solution for the potential is formulated via the Ewald sums method applied to Yukawa-type potentials. We show that, for the present Universe, the formulas involving plain summation of Yukawa potentials are preferable for computational purposes, as they require a smaller number of terms in the series to reach adequate precision.


Introduction
The shape of the space, whether it is positively curved, negatively curved, or flat, and whether there is a limit to the size of the Universe are all among essential topics of contemporary debate in theoretical physics and cosmology. Spatial topology of the Universe, its function at the very early stages of evolution (in the quantum gravity regime), and in the later process of large scale structure formation are quite interesting questions yet to be answered. General Relativity does not favor any particular topology; hence, on theoretical grounds, the space might be simply connected, in agreement with concordance cosmology, or equally as well, multiply connected. It is worth noting that some "exotic" non-simply connected spacetimes are timelike geodesically incomplete, since they have singularities [1]. Hence, such topologies are not viable.
If the Universe is multiply connected, it may have a finite volume and yet be negatively curved or flat [2]. The current available data cannot reveal the finiteness of its volume if the Universe covers a much wider region than the observable sector. However, a rather smaller volume points at the possibility of finding observational indications of its topological features [3]. For instance, a photon can travel plenty of times across the volume of multiply connected space and, thus, generate multiple images of the emitting source as a signature [4,5]. Spaces with toroidal topology in one to three dimensions may be presented as common examples of multiply connected spaces. To this class belongs the three-torus T × T × T, chimney T × T × R, and slab T × R × R topologies.
There are various comprehensive studies in the literature on potential indicators of the shape of the space [6][7][8][9][10], and the majority of research is focused on their relation to the Cosmic Microwave Background (CMB) data. Indeed, there exists a very appealing conjectural relation which suggests that CMB anomalies in large angular scale observations, e.g., the suppression of the quadrupole moment and the quadrupole and octopole alignment, are imprints of spatial topology [11,12]. The weak wide-angle temperature correlations in the CMB can be also explained, e.g., with the help of dodecahedral topology of the Universe [13].
In the present work we study the chimney topology, which admits a single infinite axis subject to interpretations such as the preferred direction of the quadrupole and octopole alignment and the commonly named "axis of evil" [14] (see [15] for additional observable signatures of a preferred axis).
In connection with the investigation of possible topological imprints in CMB observations, Planck 2013 data [2] place the constraint R i > 0.71χ rec on the radius of the largest sphere that may be inscribed in the topological domain for a flat Universe with the equal-sided chimney topology. The parameter χ rec specifies the distance from the recombination surface and it is of the same order with the particle horizon, that is, χ rec ∼ 14 Gpc. For the toroidal topologies, the former restrictions on the size of the Universe from the seven-and nine-year WMAP temperature map analyses are presented in [14,16]. The smallest possible size of the fundamental topological domain for flat space, according to the seven-year WMAP results, is d = 2R LSS cos(α min ) 27.9 Gpc [12], where R LSS stands for the distance from the last scattering surface.
In this paper, we study the chimney topology T × T × R in terms of the gravitational characteristics of the Universe, manifested in the shape of the gravitational potential and force. In the cosmological setting, the inhomogeneous gravitational field is sourced by fluctuations in the matter density [17] and, as expected, in the Newtonian limit, the potential satisfies the Poisson equation. The form of the gravitational potential in the case of toroidal topologies was previously studied in [18]. Particularly, the authors have shown that there exists no physically justified nontrivial solution of the Poisson equation for the T × T × R model. On the other hand, by employing the perturbed Einstein equations from the very beginning, one automatically includes the essential relativistic effects in the formulation and, for the gravitational potential, obtains a Helmholtz-type equation instead of the Poisson one [19][20][21]. Quite remarkably, as we show in the present work for the chimney topology, it then becomes possible to obtain exact solutions of this equation that are nontrivial and physically meaningful. Herein, we derive distinct expressions for the gravitational potential and force through alternative methods and point out the particular solutions appearing in the form of summed Yukawa potentials as ready-to-use notable sources for numerical computations. It is worth mentioning that, in the above outlined approach, we make no presumptions regarding the spatial distribution of gravitating bodies.
The outline of the paper is as follows. In Section 2, following [22], we introduce the main equations and derive alternative expressions for the gravitational potential using distinct methods, now including the Ewald technique. Subsequently, in Section 3.1, extending the results of [22], we compare these expressions in view of their usefulness for numerical computations. In Section 3.2, we obtain the gravitational force expressions for each form of the potential. We briefly review the results of our work in the concluding Section 4.

The Model and Basic Equations
It is well known that the gravitational potential Φ, created by fluctuations in the matter density, is defined by scalar perturbations of the metric coefficients [23] and that in the framework of General Relativity, it satisfies the linearized Einstein equation (see, e.g., [24,25]). Ignoring peculiar velocities, in the case of the ΛCDM cosmological model, this equation reads [19][20][21] where κ ≡ 8πG N /c 4 (with G N and c being the Newtonian gravitational constant and the speed of light, respectively), a is the scale factor, while ∆ represents the Laplace operator in comoving coordinates. Here, ρ andρ = const are the comoving mass density and its averaged value, respectively. As we operate within the ΛCDM model, matter is pressureless, and we consider it in the form of discrete point-like gravitating bodies with masses m n to represent, e.g., galaxies. Therefore, the comoving mass density The 0 subscript of Φ in Equation (1) refers to the fact that peculiar velocities have been disregarded (see also [26]).
The shifted gravitational potential fulfils the equation where the screening length [19] λ ≡ 3κρc 2 2a 3 The presence of the term ∝ Φ 0 in Equation (1) (consequently, ∝ Φ 0 in Equation (4)) results in the Yukawa-type cutoff of the potential with the characteristic length λ. The term ∝ Φ enters as a summand into energy-momentum fluctuations generating metric perturbations [19].
In what follows, the overhat indicates that the gravitational potential is shifted. A significant bonus of working with the shifted potential is that it is now possible to employ the superposition principle in solving Equation (4): once we find a solution for a single particle that is located at the center of Cartesian coordinates, we may immediately generalize it for a collection of particles at random positions.
We consider the space with chimney topology T 1 × T 2 × R, where the tori T 1 and T 2 have periods l 1 and l 2 along, e.g., the xand y-axes, respectively. Hence, each gravitating body has its images positioned away from the original point in multiples of periods l 1 and l 2 along the corresponding axes. Now, let us place a particle with mass m at the center of Cartesian coordinates. For the above indicated topology, the delta functions δ(x) and δ(y) may be presented as which implicitly include the contribution from the images of the particle. Consequently, Equation (4) for this particle reads so, we are motivated to consider the solution where the coefficients C k 1 k 2 (z) satisfy the equation Using the condition d 2 |z|/dz 2 = 2δ(z), we can easily obtain the explicit expressions for the coefficients C k 1 k 2 (z), so that the shifted gravitational potential for the selected particle and all its images eventually reads The above expression has the correct behavior in the Newtonian limit in the neighborhood of the considered particle, where it is no longer possible to distinguish between different (infinite and periodic) axes. For such regions, the summations in (10) may be replaced by the integrals: for k x ≡ k 1 /l 1 and k y ≡ k 2 /l 2 . Introducing the vectors k = (k x , k y ), r = (x, y) with the absolute values k = k 2 x + k 2 y and r = x 2 + y 2 , and assuming an angle ϕ between them, we obtain where Z = az and R = ar represent the physical distances, and the last integration is performed by using the formula 2.12.10(10) of [27]. Evidently, for a system of randomly positioned gravitating bodies, we have For linear fluctuations, the averaged value of this expression is equal to zero, as it should be (also see [28]). Indeed, and, hence, the spatial average of the total gravitational potential is For the sake of simplicity, here we have considered the particular configuration in which all N bodies in the volume V = l 1 l 2 L z are assigned identical masses m.
Equation (4) is of Helmholtz type and we can likewise solve it by considering the contribution of periodic images. In this case, the resulting expression consists of summed Yukawa potentials attributed to each one of them: As we have noted previously, the peculiar motion of gravitating bodies is disregarded in Equation (1) and, consequently, in (4). Nevertheless, the significance of such a contribution has recently been pointed out in [29], where the authors have also shown that peculiar velocities may be effectively restored by employing the effective screening length λ eff (given by the Formula (41) of [29]) instead of the screening length λ in Equations (1) and (4). Specifically, in the matter-dominated epoch, the two quantities λ eff and λ are related to one another as λ eff = √ 3/5λ. Returning to our formulation, the effect of peculiar motion is included by replacing λ with λ eff in the Formulas (10) and (16), which yields For simpler demonstration, we have assumed l 1 = l 2 = l and introduced the rescaled quantities Two alternative solutions are labeled with the subscripts "cos" and "exp" in Equations (17) and (18). Now that the peculiar velocities are included in the calculations, the 0 subscript is omitted in the new formulas. There is also a third way to express the gravitational potential for the given topology. Indeed, Yukawa-type interactions that are subject to periodic boundary conditions can be formulated via Ewald sums, so that the expression for the potential consists of two rapidly converging series, one in each of the real and Fourier spaces. The technique is commonly employed while modeling particle interactions in plasma and colloids, and, in such a context, the corresponding potential for quasi two-dimensional systems, i.e., three-dimensional systems with two-dimensional periodicity, has previously been derived in [30,31]. Being implemented in the cosmological setting considered in our paper, the discussed expression for the "Yukawa-Ewald" potential reads where and In these formulas, erfc is the complementary error function and the free parameter α, as indicated in [31], is to be chosen in such a way that a balanced interplay of computational cost and satisfactory precision is achieved. For definiteness, we set α equal toλ eff .
In the forthcoming section, we will compare three expressions and present the optimum formula in view of its efficiency in use for numerical analysis.

Gravitational Potentials
Formulas (17), (18) and (20) describe the gravitational potential due to a point-like body, with mass m, placed at (x, y, z) = (0, 0, 0), and by the accompanying images placed at (x, y, z) = (k 1 l, k 2 l, 0), where k 1,2 = 0, ±1, ±2, . . . All three forms of the rescaled potential are composed of infinite series. Hence, for any desired precision, one needs to determine the minimum number of terms that are needed to numerically calculate the potential. The criterion that we use to specify this number n is the following: the ratio |exactΦ − approximateΦ|/|exactΦ| should be less than 0.001. This defines the order of accuracy in our analysis. Evidently, for each form of the potential the number n can be different, so we denote these as n exp , n cos and n mix , correspondingly. The formula requiring the smallest number of terms to defineΦ (i.e., to calculate "approximateΦ") up to the adopted accuracy is clearly the best alternative for numerical computation purposes. In this connection, we are interested in comparing (17), (18) and (20) here. Because the formulas include double series, the accompanying numbers n are to be ascribed the smallest possible number of combinations (k 1 , k 2 ) that can provide the necessary precision. We find these by listing the summands in increasing order for k 2 1 + k 2 2 and assigning to n the total number of terms that are included in the list eventually. This procedure of generating a sequence of combinations (k 1 , k 2 ) and finding n is performed using Mathematica [32]. Tables 1 and 2 show the outputs for eight selected points. As to the adopted accuracy, for all n n exp , the approximate value ofΦ exp (calculated by (18)) differs from the exact value by less than one tenth of a percent. In both tables, the exact valueΦ is calculated from (18) for n n exp . The quantities n cos and n mix indicate the numbers of terms in formulas (17) and (20), respectively, which one needs to keep in order to obtain the same values of the potential at the selected points with the same precision as attained by using (18). In the n cos column, the dash reflects either incorrect outputs that are produced because of complications in the computational process, or the fact that unreasonably large number of summands is necessary. Because our results depend on the ratio of λ eff to the physical size al of the periods of tori, i.e., λ eff = λ eff /(al), we present the results that were obtained for small and large values ofλ eff separately in Tables 1 and 2, which includeλ eff = 0.01, 0.1 andλ eff = 1, 3, respectively. According to these tables, both expressions (18) and (20) seem preferable for numerical calculations in the caseλ eff < 1 since n exp , n mix n cos , although Equation (20) is, of course, much more complicated than Equation (18), and the computation of its every single summand takes longer. However, forλ eff 1, the Yukawa-Ewald formula (20) alone becomes superior to the remaining two and the distinction grows asλ eff becomes larger. According to Planck 2013 data [2], the lower limit on the periods of tori (in the case of chimney topology) is ∼20 Gpc. Meanwhile, the current value of the effective cosmological screening length, as indicated in [29], is 2.6 Gpc. Thus, the region that is defined byλ eff < 1 depicts the observable Universe and, here, as we have just discussed, Equation (18) is more convenient for numerical analysis. Concluding this section, we also present Figures 1-4, demonstrating the shape of the rescaled potentialΦ for the same values ofλ eff as those picked for Tables 1 and 2. To plot these figures (using Mathematica [32]), we use (18) for n n exp .

Gravitational Forces
It is also interesting to study the forces (per unit mass) associated with the alternative forms of the gravitational potential derived in the previous section. We intend to consider the projections of these forces on the xand z-axes. Owing to the symmetry of the model, the x and y projections are similar. We calculate the gravitational forces for the same points as for the potentials and, among these, naturally, we only investigate the points at which projections on the axis of interest are nonzero. In this connection, the points A 1 , A 2 , A 3 , C 1 , C 2 and the points A 3 , B 3 are omitted for the xand z-components, respectively. The accuracy of force calculations is of the same level as that of the potentials. Here, once again, we compare the number of terms needed to achieve this accuracy, but now for three different forms of the gravitational force presentation.

x-Component of the Gravitational Force
From (17), (18) and (20), we derive three alternative expressions for the x-component of the rescaled force: where We present the results of the calculations that were performed in Mathematica [32] in Tables 3 and 4 forλ eff = 0.01, 0.1 andλ eff = 1, 3, respectively. As in the case for the gravitational potential, a straightforward analysis shows that Formulas (24) and (25) that are related to the Yukawa and Yukawa-Ewald potentials, respectively, are preferable over (23) for the physically significant caseλ eff < 1 (although the structure of the expression (25) is again much more complicated compared to (24)). Meanwhile, whenλ eff 1, the Yukawa-Ewald force becomes superior. In both tables, Equation (24) was employed for n n exp while computing the values of the rescaled x-componentΦ x . Table 3. Numerical values of the x-component of the rescaled forceΦ x ≡ ∂Φ/∂x as well as the numbers n exp , n cos and n mix for points B 1 , B 2 and B 3 forλ eff = 0.01 andλ eff = 0.1 in the left and right tables, respectively.

z-Component of the Gravitational Force
For the z-component of the rescaled force, three alternative formulas are: where, for simplicity,z > 0, and where and C ∓ are given by (26). Now, we employ these formulas to calculate the nonzero z-components of the gravitational force at the previously selected set of points and, again, for the desired precision. The results that were obtained in Mathematica [32] are presented in Tables 5 and 6, which show that whileλ eff < 1, depicting well the observational restrictions, two Formulas (28) and (29) are favorable (as before, the latter is much more cumbersome). On the other hand, forλ eff 1, the Yukawa-Ewald Formula (29) gives the best results. Herein, the quantityΦ z stands for the z-component of the rescaled force, calculated from Equation (28) for n n exp . We depict the behavior of this component in Figures 9 and 10 for the section y = 0. Obviously, the projection of the gravitational force on the z-axis is absent for the section z = 0 due to the symmetry of the model. Table 5. Numerical values of the z-component of the rescaled forceΦ z ≡ ∂Φ/∂z as well as the numbers n exp , n cos and n mix for points A 1 , A 2 , B 1 , B 2 , C 1 and C 2 forλ eff = 0.01 andλ eff = 0.1 in the left and right tables, respectively.

Conclusions
In this work, we have studied how the chimney topology T × T × R of the Universe affects the form of the gravitational potential and, consequently, that of the gravitational force. In this connection, we have proposed three alternative forms for each of the solutions. One of them (see Equation (17)) relies on the Fourier expansion of the delta functions into series while using periodicity in two toroidal dimensions in the model. The second one (see Equation (18)) follows from the summation of solutions of the Helmholtz equation, each in the form of the Yukawa potential, for a source mass and all of its periodic images. Finally, the third form of the potential (see Equation (20)) is formulated via the Ewald sums for Yukawa potentials. Subsequently, we have presented three alternative forms of the gravitational force (see Equations (23)- (25) and (27)-(29) for the xand z-components, respectively) derived from the potential expressions.
In all three alternative forms, the screening lengthλ eff serves as a crucial parameter, as it specifies the distance (from the source or the periodic images) where the gravitational potential undergoes exponential cutoff. This fact is most clearly demonstrated in the case where the solution takes on the form of summed Yukawa potentials (18). The observational data show that this parameter should be less than 1 (λ eff < 1) in today's Universe.
One of the main goals of this work was to reveal which of the obtained alternative formulas would serve better as a tool to be employed in numerical calculations. Namely, to show which formula would require less terms in the series to reach the desired precision. Our calculations have demonstrated that, for both the gravitational potentials and forces, two expressions involving plain summations of Yukawa potentials are preferable for the physically significant caseλ eff < 1. However, in the caseλ eff 1, the Yukawa-Ewald presentation stands out as the best alternative.

Conflicts of Interest:
The authors have no conflicts of interest to declare that are relevant to the content of this article.