The non-Gaussian distribution of galaxy gravitational fields

In models of galaxies and their structure formation, the distribution of gravitational fields of their constituents plays the decisive role. Many scenarios of such formation have been around for some time (Peebles 1969; Sunyaev & Zeldovich 1972; Zel'dovich 1970; Doroshkevich 1973; Shandarin 1974; Dekel 1985; Efstathiou & Silk 1983). Under the influence of new observational data, these scenarios are constantly being revised and improved, see Shandarin et al. (2012); Giahi-Saravani & Schäfer (2014) and references therein for the latest discussion. The main controversy here is how galaxies acquire their angular momentum, which subsequently influences those of galaxy clusters and larger structures. On the other hand, this angular momentum acquisition is intimately related to the above distribution of gravitational fields.

The commonly accepted model of the Universe is the spatially flat homogeneous and isotropic lambda cold dark matter (ΛCDM) model. Galaxy clusters in this model are formed as a result of adiabatic and almost scale invariant Gaussian fluctuations (Silk 1968; Peebles & Yu 1970; Sunyaev & Zeldovich 1970). This assumption is the basis of the so-called hierarchical clustering model (Doroshkevich 1970; Dekel 1985; Peebles 1969), the most popular scenario of galaxy formation. Note, however, the presence of models with non-Gaussian initial fluctuations, see Bartolo et al. (2004) and references therein. This non-Gaussianity, however, has been postulated to be in a certain form rather than being calculated. At the same time, the non-Gaussian distributions can be obtained from initial Gaussian ones as a result of time evolution in generalized stochastic models, where probability density functions (pdfs) are obtained from the solutions of Fokker-Planck type differential equations with so-called fractional derivatives (Garbaczewski & Stephanovich 2009; Garbaczewski et al. 2011). In other words, the initial Gaussian fluctuations (if any) may become non-Gaussian as a result of primordial, fast time evolution. After that, slower evolution, dictated by the ΛCDM scenario, occurs. Although here we do not present the details of this primordial time evolution, one of the aims of the present paper is to draw attention to the method, which permits calculating the non-Gaussian distribution function, based solely on the form of interaction between astronomical objects. This distribution function is a terminal function for the above initial fast time evolution process.

In hierarchical clustering type scenarios, the large scale structure forms from bottom to top as a consequence of gravitational interactions between the constituents. This means that smaller structures like galaxies are formed first with their later merger into larger ones. Consequently, the spin angular momenta of galaxies arise as a result of tidal interaction with their neighbors (Schäfer 2009). In the original hierarchical clustering scenarios, the completely random distribution of galaxy angular momenta was obtained. Note however, that it was shown later that the local tidal shear tensor can generate a local alignment of galaxy angular momenta (Catelan & Theuns 1996a,b; Lee & Pen 2002; Navarro et al. 2004). The mechanisms of this type are known as tidal torque mechanisms, which had first been introduced by Hoyle (1951) and later by White (1984), see also the review paper Kiessling et al. (2015) that discusses galaxy alignments in the context of gravitational lensing. Note that in our model, angular momentum is the result of tidal interaction with the entire environment, which occurs via interaction transfer from close to distant galaxies, see below. In this sense our approach is an improvement on those considered by Schäfer & Merkel (2012); Catelan & Theuns (1996a,b); Lee & Pen (2002), where the "mean" tidal interaction with the entire environment has been considered.

The above tidal torque mechanism has an opposite idea, which is constructed on the basis of the Zeldovich pancake model (Sunyaev & Zeldovich 1972; Doroshkevich 1973; Shandarin 1974). In this model, structures in the Universe arise from top to bottom. The crucial role here is played by a magneto-hydrodynamic shock wave, which makes a large structure fragment. This shock wave arises as the result of asymmetrical collapse of an initial large structure and also imparts galaxies with spin angular momentum. The model predicts a coherent, non-random spatial orientation of galaxy planes with rotational axes of galaxies being parallel to the main plane of a large structure.

In the model of primordial turbulence, the spin angular momentum is a remnant of the primordial whirl (von Weizsäcker 1951; Gamow 1952; Ozernoi 1978; Efstathiou & Silk 1983). As a result it is obtained that the rotational axes of galaxies are not oriented randomly. The preferred direction of angular momentum of galaxies is perpendicular to the initial large structure's main plane.

It had been pointed out in Gamow (1946); Gödel (1949) and later in Collins & Hawking (1973) that, if the Universe is rotating, the emerging galaxy angular momentum is a consequence of its conservation in a rotating Universe. At that time, the argument against was that this model predicts the alignment of galaxy rotational axes, which has not been confirmed observationally (see Godłowski 2011 for details). Based on this idea, Li (1998) introduced a model in which galaxies form in a rotating Universe.

We emphasize that simple picture, in which each of the above approaches (primordial turbulences, hierarchical clustering and Zeldovich pancakes) predicts different ways of ordering galaxy rotational axes that are not completely true. The point is that in each of the above models, including hierarchical clustering, a phase with a shock wave can appear. Later, the phase is usually accompanied by the collapse of structures or substructures (Melott & Shandarin 1989; Sahni et al. 1994; Pauls & Melott 1995; Mo et al. 2005; Shandarin et al. 2012), which may generate ordering of the rotational axes. Apparently, the scale of such orientation is different in different models. For instance, in Bower's scenario (Bower et al. 2006), we do not have hierarchical clustering for all scales of masses. Instead, we have anti-hierarchical clustering on small scales as tidal interaction effects yield the Zeldovich pancake-like emergence of objects (Zel'dovich 1970) rather than spherically collapsing halos. There is, however, a fundamental difference with the above classical pancake scenario. Namely, anti-hierarchical clustering is local as it occurs on a small scale.

The model of hierarchical clustering is the only model explicitly taking into account the existence of dark matter. The Li model originally considered the Universe as a dust fluid, however, nothing prevents introducing dark matter as a background. As a result, in this model, dark matter is not interacting with observable matter in any other way than gravitational forces. In the remaining models, namely primordial turbulences and the Zeldovich pancake model, only the dust component has been considered so that there are no clear and successful attempts to introduce dark matter there. Therefore, we exclude both models from present consideration.

Theoretical models of galaxy formation have problems with explaining the observational dependence between structure angular momentum and its mass. This dependence can only be seen in two classes of models. There are the tidal torque scenario (Heavens & Peacock 1988; Catelan & Theuns 1996a; Hwang & Lee 2007; Noh & Lee 2006a,b) and the Li model (Li 1998; Godłowski et al. 2005). The remaining models do not anticipate such dependence.

Comparing the two models, we should note that the Li model needs a global or at least large scale rotation of the Universe. Li (1998) studied the dependence between angular momentum and the mass of spiral galaxies and he estimated the rotation of the Universe to be close to the value obtained by Birch (1982). However, the obtained value is too large compared to observed anisotropy in cosmic microwave background radiation. Hence, in the present paper, we only consider the tidal torque scenario.

In the present work we perform a comprehensive theoretical analysis of the influence of tidal interaction between astronomical objects on formation of the larger (than initial constituents) structures. The idea of our approach is to use the statistical method originally proposed by Chandrasekhar (1943), where we also account for dark matter halos. The statistical method of Chandrasekhar (1943) permits deriving the distribution functions of gravitational fields and angular momenta of stellar components. Our main result is that in stellar systems with multipole (tidal) gravitational interaction, the derived distribution function cannot be Gaussian. Instead we obtain the pdf which rather belongs to the family of so-called "heavy-tailed distributions" (Garbaczewski & Stephanovich 2009; Garbaczewski et al. 2011; Van Kampen 1992). As we have mentioned above, the obtained non-Gaussian pdf is a result of fractional time evolution for initial Gaussian fluctuations. This function allows us to calculate the distribution of virtually any observable (like angular momentum) of astronomical structures (not only galaxy clusters but also smooth components like halos, whose mass dominates the total mass of the cluster, see Kravtsov & Borgani (2012) in any (linear or nonlinear) Eulerian approach.

The paper is organized as follows. To make the paper self-contained, in Section 2 we shortly recollect our method (Stephanovich & Godłowski 2015), emphasizing its points, important for present consideration. Some technical details are described in the Appendix. In Section 3 we discuss the problem of angular momentum pdfs. We show that different (physically reasonable) assumptions about the structure of galaxy clusters generate different relations between their mass M and average angular momentum L. We demonstrate that it is possible to derive not only the relation $L\sim {M}^{4/3}$ (like in Stephanovich & Godłowski 2015) but also recover the well-known empirical relation $L\sim {M}^{5/3}$. We also show that while it is possible to discriminate between the above model assumptions theoretically, present observational data are not sufficient to come to an unambiguous conclusion. In addition, we discuss the possibility of observationally testing our theoretical results related to the time evolution of the angular momentum pdf and its mean value L. We conclude our article with Section 4.

We consider tidal interaction in the ensemble of galaxies and their clusters in a Friedmann-Lemaître-Robertson-Walker Universe with Newtonian self-gravitating dust fluid (p = 0) containing both luminous and dark matter. The tidal (shape distorting) interaction between astronomical objects can be derived by the multipole expansion of Newtonian interaction potential between fluid elements (Poisson 1998). Limiting ourselves to the quadrupolar term, we write the Hamiltonian function of interaction between above elements in the form

$$\begin{eqnarray}{\mathcal{H}} & = & -G\displaystyle \mathop{\sum }\limits_{{ij}}{Q}_{i}{m}_{i}V({{\boldsymbol{r}}}_{{ij}}),\\ V({\boldsymbol{r}}) & = & \displaystyle \frac{1}{2}\displaystyle \frac{3\,{\cos }^{2}\theta -1}{{r}^{3}},\end{eqnarray} \tag{ 1 }$$

where G is the gravitational constant, Q_i and m_i are, respectively, the quadrupole moment and mass of the i-th object, ${r}_{{ij}}\equiv$ $|{{\boldsymbol{r}}}_{{ij}}|$, ${{\boldsymbol{r}}}_{{ij}}={{\boldsymbol{r}}}_{j}-{{\boldsymbol{r}}}_{i}$ is a relative distance between objects while θ is the apex angle. Our Hamiltonian function, Equation (1), is obtained for the ensemble of N objects, thus generalizing the result of Poisson (1998) for two particles.

Note that the Hamiltonian function Equation (1) describes the interaction of quadrupoles, formed both from luminous and dark matter. This is important as real galaxies, formed from luminous matter, reside inside dark matter halos that are much more extended and massive. In other words, the Hamiltonian function, Equation (1) (and subsequent results), already contains information about dark matter halos. We have discussed this question in our previous work (Stephanovich & Godłowski 2015). The main point was that properties of luminous matter (like galaxies and their clusters) give us information about dark matter (sub)structures. This point is corroborated by observations (see, e.g. Paz et al. 2008; Bett et al. 2010; Kimm et al. 2011; Varela et al. 2012) that angular momentum of luminous matter is correlated with that of corresponding dark matter halos. Below we will calculate the angular momentum of luminous astronomical structures. Our formalism can be generalized to describe not only this situation, but also the structures with a larger smooth component. Namely, in general, luminous galaxy matter is not only surrounded by dark matter halos, but also (along with latter halos) submerged in the "mud," which is hypothetical intergalaxy dark matter. We plan to fulfill this interesting generalization in our subsequent publications.

In the function defined in Equation (1), we split the interaction between many stellar objects (particles) in pairs, see Appendix for details. Such splitting is usual, for instance in the theory of magnetism, where the spins of the interacting ensemble are represented by the sum of all possible couplings between particle pairs i and j. For instance, three particle interaction may be decomposed as 123 = 12 + 13 + 23, see, e.g. Mattis (2006).

The Hamiltonian function in Equation (1) describes the pairwise, shape-distorting interaction between the structures. Namely, this interaction distorts the shape of a given i-th object, which alters its density field ${\rho }_{i}({\boldsymbol{r}})$. As the objects have random shapes, their masses m_i and quadrupole moments Q_i vary randomly, as does the gravity field ${{\boldsymbol{E}}}_{\mathrm{quad}}$ from these quadrupoles. One should note that the latter field is in fact a gradient of the potential, given by Equation (1). It has the form

$$\begin{eqnarray}&&{{\boldsymbol{E}}}_{\mathrm{quad}}({\boldsymbol{r}})={{\boldsymbol{i}}}_{r}{E}_{0}\displaystyle \frac{3\,{\cos }^{2}\theta -1}{{r}^{4}},\end{eqnarray} \tag{ 2 }$$

where ${E}_{0}={GQ}/2$ and ${{\boldsymbol{i}}}_{r}$ is the unit vector in the radial direction.

According to the statistical method of Chandrasekhar (1943), the distribution function of random quadrupolar fields is

$$\begin{eqnarray}&&f({\boldsymbol{E}})=\overline{\delta ({\boldsymbol{E}}-{{\boldsymbol{E}}}_{i})},\end{eqnarray} \tag{ 3 }$$

where $\delta (x)$ is a Dirac δ function, while ${{\boldsymbol{E}}}_{i}\equiv {{\boldsymbol{E}}}_{\mathrm{quad}}({{\boldsymbol{r}}}_{i})$ is given by Equation (2) where the bar means averaging over spatial (and any other possible) disorder. Moreover, if all objects in the ensemble are similar (no randomness), then the distribution function is represented by the simple delta function, centered at the field ${{\boldsymbol{E}}}_{i}$. The disorder broadens this delta function, giving rise to a "bell-shaped" continuous probability distribution, see Stephanovich (1997); Semenov & Stephanovich (2003, 2002) and references therein.

The explicit averaging in Equation (3) is performed with the help of the integral representation of a Dirac δ function, see Stephanovich & Godłowski (2015) for details. The idea is that the mass and quadrupole moment of the object in volume V obey a uniform distribution with probability density equal to 1/V. In such a case we introduce the number of objects N so that in the limit $N\to \infty$ and $V\to \infty$, their density n = N/V remains constant. The final expression for the distribution function Equation (3) reads (Stephanovich & Godłowski 2015)

$$\begin{eqnarray}&&f({\boldsymbol{E}})=\displaystyle \frac{1}{(2\pi {)}^{3}}\displaystyle {\int }_{-\infty }^{\infty }{e}^{i{\boldsymbol{E}}{\boldsymbol{\rho }}-F(\rho )}{d}^{3}\rho,\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&F(\rho )=n\displaystyle {\int }_{V}\left[1-\displaystyle \frac{\sin \rho E({\boldsymbol{r}})}{\rho E({\boldsymbol{r}})}\right]{d}^{3}r.\end{eqnarray} \tag{ 5 }$$

In this case $F(\rho )$ is in fact the characteristic function for the distribution of random gravitational fields. Note also that characteristic function $F(\rho )$ depends only on modulus ρ and not on its angles. This will result (see Eq. (6) below) in the only field modulus dependence on the pdf of random gravitational fields. The reason is that we take only the zz component of the quadrupolar field in Equation (2). If we need the complete (i.e. including its possible angular dependence) distribution function of vector ${\boldsymbol{E}}$, we should account for the complete tensor structure of the Hamiltonian (1), ${\mathcal{H}}=-G{\sum }_{{ij}\alpha \beta }{Q}_{i\alpha \beta }{m}_{j}{V}_{\alpha \beta }({{\boldsymbol{r}}}_{{ij}})$, $\alpha,\beta =x,y,z$. Such account (Stephanovich 1997; Semenov & Stephanovich 2003, 2002), while not changing our conclusions qualitatively (and in many cases quantitatively, see below), will make the problem tractable only numerically. At the same time our present approach permits us to gain analytical insights into the problem (for example investigating the implication of non-Gaussian character in the distribution function for gravitational fields), which is a good starting point for future numerical simulations. One more justification of the radial distribution is the results of numerical simulations in the halo model (Schneider & Bridle 2010), where the axes of galaxies embedded in a dark matter halo were preferentially radially oriented.

Moreover, the spin angular momentum is usually known only for very few galaxies and other structures. For this reason, the spatial orientation of galaxies (see, for example, Flin & Godlowski 1986; Romanowsky & Fall 2012) is studied instead of their angular momenta. Alternatively, only the distribution of position angles of galaxy planes is analyzed in Hawley & Peebles (1975).

In more realistic models of galaxy clustering, we can assume that the density of stellar objects, like galaxy density, is not a constant but rather depends on their separation $n=n({\boldsymbol{r}})$. The other factor, which may improve the coincidence with observational results, is to consider galaxy clustering within a model of the inhomogeneous distribution of masses (and/or quadrupolar moments) in the large scale structure. The idea here is to introduce the distribution function of masses $\tau (m)$, which had been put forward by Chandrasekhar (1943).

It is important that distribution function $f({\boldsymbol{E}})$ (4) in a general case could be much more complicated than simply Gaussian. We had shown in Stephanovich & Godłowski (2015) that for the multipole interaction between astronomical objects, the function in Equation (4) does not admit a Gaussian limit. The calculation of $F(\rho )$ (5) generates the following explicit form of f(E) (Stephanovich & Godłowski 2015)

$$\begin{eqnarray}f(E) & = & \displaystyle \frac{1}{2{\pi }^{2}E}{\displaystyle \int }_{0}^{\infty }\rho {e}^{-\alpha {\rho }^{3/4}}\sin {\rho E}\,{d\rho },\\ \alpha & = & 2\pi n\cdot 0.41807255\cdot {E}_{0}^{3/4}.\end{eqnarray} \tag{ 6 }$$

The expression in Equation (6) is the chief theoretical result of our studies. The distribution function (6) depends parametrically on the objects (i.e. both luminous and dark matter) density n, and on average quadrupole moment Q.

The normalization condition for distribution function (6) reads

$$\begin{eqnarray}&&4\pi \displaystyle {\int }_{0}^{\infty }{E}^{2}f(E){dE}=1.\end{eqnarray} \tag{ 7 }$$

As we have shown previously (Stephanovich & Godłowski 2015), the distribution function of the gravitational fields cannot be Gaussian for multipole interaction between galaxies or any other astronomical objects including elements of dark matter halos. However, all previous theories postulated the distribution function to have a Gaussian form rather than calculated it. We mention here that a non-Gaussian distribution has also been postulated rather than calculated in Bartolo et al. (2004). In our opinion, the non-Gaussian, heavy-tailed nature of the above pdf captures the essential physics of systems with long-range gravitational multipole interaction. Namely, the long-range interaction in such systems makes the objects (galaxies, their clusters and even the dark matter halos) also interact with each other at very large separations. This, in turn, implies nonzero probabilities of such configurations, contrary to the case of the Gaussian distribution, generated by short-range interactions. Below we will see the important implications of this fact.

To plot the function f(E), we define the dimensionless variables $\rho E=x$ and $\beta =E/{\alpha }^{4/3}$. In these variables the integral (6) assumes the form

$$\begin{eqnarray}&&f(\beta )=\displaystyle \frac{H(\beta )}{4\pi {\beta }^{2}{\alpha }^{4}},\quad H(\beta )=\displaystyle \frac{2I(\beta )}{\pi \beta },\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&I(\beta )=\displaystyle {\int }_{0}^{\infty }x\sin x\exp \left[-{\left(\displaystyle \frac{x}{\beta }\right)}^{3/4}\right]{dx}.\end{eqnarray} \tag{ 9 }$$

The physical interpretation of function $H(\beta )$ is the following. This function gives the effective one-dimensional (1D) distribution function of random gravitational fields. This is because the normalization condition for $H(\beta )$ is effectively the 1D form ${\int }_{0}^{\infty }H(\beta )d\beta =1$, see (7). In this case, the average value $\bar{\beta }$ of the dimensionless random field β has the form $\bar{\beta }={\int }_{0}^{\infty }\beta H(\beta )d\beta$. The mean value $\bar{\beta }$ exists if the integral $H(\beta )$ is convergent. Asymptotic analysis of the function $f(\beta )$, which had been performed in our previous work (Stephanovich & Godłowski 2015), shows that $f(\beta )$ does not depend on β for small β and decays as ${\beta }^{-7/4}$ at large β. The character of decay at large β shows that although the normalization integral is convergent, already the first moment does not exist. Such behavior is a characteristic feature of so-called heavy-tailed distributions (Garbaczewski & Stephanovich 2009; Garbaczewski et al. 2011).

Our aim is to derive the distribution function of angular momenta. For that we need to calculate how the angular momentum ${\boldsymbol{L}}$ of a stellar object depends on its gravitational field ${{\boldsymbol{E}}}_{\mathrm{quad}}({\boldsymbol{r}})$ (2). The expression for angular momentum components L_α ($\alpha =x,y,z$) could be obtained perturbatively in small Lagrangian coordinate ${\boldsymbol{q}}$. One should note that the first order terms were obtained in equation (11) of Catelan & Theuns (1996a), while the second order ones in their next article, Catelan & Theuns (1996b) (eq. (28)). Note that both equations have identical structure, i.e.

$$\begin{eqnarray*}&&{L}_{\alpha }^{(i)}={f}_{i}(t){\varepsilon }_{\alpha \beta \gamma }{E}_{i\beta \sigma }{I}_{\sigma \gamma },\quad \alpha,\beta,\gamma,\sigma =x,y,z,\end{eqnarray*}$$

where index i = 1, 2 defines the order of perturbation theory, ${\varepsilon }_{\alpha \beta \gamma }$ is the Levi-Civita symbol, ${E}_{\beta \sigma }$ are components of quadrupole (tidal) field (2) while ${I}_{\sigma \gamma }$ represent components of the inertia tensor.

In order to obtain the distribution function of modulus of E (and subsequently L), it is sufficient to take the zz component in (2). If we need the complete distribution function of vector ${\boldsymbol{E}}$, we should account for the complete tensor structure of Hamiltonian (1)

$$\begin{eqnarray*}{\mathcal{H}} & = & -G\displaystyle \mathop{\sum }\limits_{{ij}\alpha \beta }{Q}_{i\alpha \beta }{m}_{j}{V}_{\alpha \beta }({{\boldsymbol{r}}}_{{ij}}),\\ \alpha,\beta & = & x,y,z.\end{eqnarray*}$$

Also, as ${\boldsymbol{L}}$ is a function of time t by means of the functions ${f}_{i}(t)$, the distribution function will be time dependent. With respect to symmetry relations ${I}_{{ab}}={I}_{{ba}}$ and ${E}_{{ab}}={E}_{{ba}}$ and leaving only E_zz, we obtain

$$\begin{eqnarray*}{L}_{x} & = & -b(t){E}_{{zz}}{I}_{{yz}},\\ {L}_{y} & = & b(t){E}_{{zz}}{I}_{{xz}},\\ {L}_{z} & = & 0,\\ L & = & \sqrt{{L}_{x}^{2}+{L}_{y}^{2}+{L}_{z}^{2}}={L}_{0}E,\\ {L}_{0} & = & {L}_{0}(t)={f}_{i}(t)\sqrt{{I}_{{xz}}^{2}+{I}_{{yz}}^{2}}.\end{eqnarray*}$$

The above equations constitute a linear relation between angular momentum and tidal field moduli. They are valid both in linear (i = 1) and nonlinear (i = 2) regimes. Because the above relation between gravitational field modulus and angular momentum is linear in both cases, it is easy to see that the shape of the distribution function of angular momenta f(L) repeats that of gravitational fields. In an explicit form, the expression for f(L) can be derived using a well known relation from probability theory $f(L)=f[E(L)]\left|\tfrac{{dE}}{{dL}}\right|$, which yields

$$\begin{eqnarray}&&f(L)=\displaystyle \frac{1}{2{\pi }^{2}L}\displaystyle {\int }_{0}^{\infty }\rho {e}^{-\alpha {\rho }^{3/4}}\sin \left(\rho \displaystyle \frac{L}{{L}_{0}(t)}\right)\ d\rho,\end{eqnarray} \tag{ 10 }$$

where ${L}_{0}(t)$ is defined above. Dimensionless variables $\rho (L/{L}_{0})=x$, $\lambda =L/({L}_{0}{\alpha }^{4/3})$ generate the pair of functions which are similar to those obtained for the distribution of gravitational fields. They read

$$\begin{eqnarray}&&f(\lambda )=\displaystyle \frac{H(\lambda )}{4\pi {\lambda }^{2}{\alpha }^{4}{L}_{0}},\quad H(\lambda )=\displaystyle \frac{2I(\lambda )}{\pi \lambda },\end{eqnarray} \tag{ 11 }$$

where $I(\lambda )$ is defined by the expression in Equation (9) and is usually referred to as the spin parameter.

The effective 1D distribution function for gravitational fields or angular momenta is presented in the left panel of Figure 1. It is seen that while initial 3D function $f(\lambda )$ decays monotonically (right panel), this function is strongly asymmetric and has a characteristic bell shape. Note that as the initial expression in Equation (1) allows for interaction between all astronomical objects in an ensemble, it also naturally considers the interaction with surrounding structures and dark matter halos. This fact renders the distribution functions of gravitational fields (8) and angular momenta (11) to account not only for isolated cluster regions, but for long-range interactions with surrounding structures as well. To be specific, the narrow peak of the distribution function in the left panel of Figure 1 stems from the closely situated cluster region, while its long tail stems from the long-range (quadrupole) interaction with surrounding structures. In other words, the interaction with surrounding structures is essential (and our distribution functions take this fact into account) as the interaction between objects in stellar ensembles have long-range multipole character.

Zoom In Zoom Out Reset image size

As we have shown in the previous article (Stephanovich & Godłowski 2015), the integral for the first moment of angular momentum pdf diverges. It is well known (see, for example, Van Kampen 1992) that for the distribution functions, which decay slowly at infinities, the corresponding mean value can be approximately estimated as the maximum of such a function. In this spirit we calculate ${\lambda }_{\max }$, corresponding to the maximum of distribution function $H(\lambda )$, as presented in Figure 1. The analysis of ${\lambda }_{\max }$ in dimensional units makes it possible to obtain some useful relations, which earlier had been guessed only empirically. To consider the characteristics of galaxies, i.e. luminous matter, here we use the idea of a halo model (Schneider & Bridle 2010), which states that galaxies (i.e. "pieces" of luminous matter) are embedded in dark matter halos so that their observable characteristics like angular momentum emerge from the mass and hence gravitational field of dark matter. Also, as the galaxies and their clusters reside in larger structures like voids and filaments, the gravitational field of the latter large objects also influences galaxies, see, e.g. Joachimi et al. (2015). As our distribution function (11) takes these effects into account by virtue of model (1), our subsequent calculations of mean angular momentum of the galaxies take the above effects into account.

Let us first consider the simplest possible cold dark matter (CDM) model in the first order of perturbation theory. In such model the evolution of scale factor is given by the equation $a{(t)=D(t)=(t/{t}_{0})}^{2/3}$ (Doroshkevich 1970) so that ${L}_{0}=\tfrac{2I}{3{t}_{0}}\tau$, $\tau =t/{t}_{0}$ and $I=\sqrt{{I}_{{xz}}^{2}+{I}_{{yz}}^{2}}$. The equation ${dH}/d\lambda =0$ has solution ${\lambda }_{\max }=0.602730263$, which gives in dimensional units

$$\begin{eqnarray}{L}_{\max } & = & 0.7281884{n}^{4/3}\displaystyle \frac{t}{{t}_{0}^{2}}{GIQ}\\ & \approx & \kappa {n}^{4/3}\displaystyle \frac{t}{{t}_{0}^{2}}{{GR}}^{4}{m}^{2},\end{eqnarray} \tag{ 12 }$$

where $n=N/V$ and $\kappa \sim 1$ is a constant. To derive Equation (12), we estimate galaxy quadrupole moment Q and its mean moment of inertia I as being proportional to ${{mR}}^{2}$, where m is mass of the galaxy while R is its mean radius. In our approach we represent volume V as $V={R}^{3}$, then R cancels in Equation (12) so that ${L}_{\max }\sim (t/{t}_{0}^{2}){m}^{2}{N}^{4/3}$. Then, we introduce the mass of a galaxy cluster $M={mN}$ and obtain

$$\begin{eqnarray}&&{L}_{\max }\sim \displaystyle \frac{t}{{t}_{0}^{2}}{M}^{5/3}{\left(\displaystyle \frac{m}{N}\right)}^{1/3}\equiv \displaystyle \frac{t}{{t}_{0}^{2}}{M}^{5/3}\displaystyle \frac{{\rho }^{1/3}}{{n}^{1/3}},\end{eqnarray} \tag{ 13 }$$

where $\rho =m/V$ is a mass density and $n=N/V$ is galaxy density. Following Catelan & Theuns (1996a), we assume that mass density ρ is a function of time, defined by the Friedmann equation in the CDM model $\dot{a}/a={H}_{0}=\sqrt{8\pi G\rho /3}$, where H₀ is the Hubble constant. This generates the dependence $\rho \propto$ t⁻², which, being substituted in Equation (13), yields

$$\begin{eqnarray}&&{L}_{\max }\sim \displaystyle \frac{{t}^{1/3}}{{t}_{0}^{2}}\displaystyle \frac{{M}^{5/3}}{{n}^{1/3}}\sim {t}^{1/3}{M}^{5/3}.\end{eqnarray} \tag{ 14 }$$

To derive Equation (14), we assume that n = const. We see that Equations (13) and (14) recover expression (27) of Catelan & Theuns (1996a), giving the theoretical derivation of the well-known empirical relation between the mean angular momentum of a galaxy ensemble (galaxy clusters) and their mass ${L}_{\max }\sim {M}^{5/3}$ (see Catelan & Theuns 1996a and references therein). Note that within the tidal torque model, the ${M}^{5/3}$-power law was first obtained by Heavens & Peacock (1988) while reasonable values for lambda in Equation (11) within the tidal torque approach were derived by Schäfer & Merkel (2012), who followed Heavens & Peacock (1988).

There is also another approach to interpreting the dependence of ${L}_{\max }$ on stellar parameters. Namely, suppose that volume $V={R}_{A}^{3}$, where R_A is a mean cluster radius, proportional to the autocorrelation radius (see Longair 2008 and references therein). In such an approach (see Stephanovich & Godłowski 2015 for details) n is still a constant for any particular cluster, but now it varies from cluster to cluster with increasing richness N. In this case we may rewrite $N=M/m$ to obtain the alternative (to Eq. (14)) form of expression for ${L}_{\max }$

$$\begin{eqnarray}&&{L}_{\max }\sim \displaystyle \frac{t}{{t}_{0}^{2}}\ {\left(\displaystyle \frac{R}{{R}_{A}}\right)}^{4}\ {m}^{2/3}{M}^{4/3},\end{eqnarray} \tag{ 15 }$$

which does not contain ρ.

It is instructive to comment on time dependence ${L}_{\max }(t)$ in Equation (15). At first sight, it follows from (15) that ${L}_{\max }\sim t$, but the problem is complicated a lot by the intricate time dependence of quantities R and R_A (Longair 2008). We plan to study this question in future works.

It is clear from Equation (12) that mean orbital momentum of a galaxy increases with the number of galaxies N and it is proportional to ${N}^{4/3}$. Moreover, even in the model with constant galaxy density n, number (richness) N varies from cluster to cluster so that the dependence ${L}_{\max }(N)={\kappa }_{2}{N}^{4/3}$ holds and shows that angular momenta increase with number of galaxies N in the analyzed structure.

The sample of 247 Abell clusters has been analyzed by Godłowski et al. (2010). Namely, the orientation of galaxies in particular clusters has been studied. The idea was to test hypotheses that the galaxy angular momenta increase with cluster richness. If galaxy clusters do not rotate (see Regos & Geller 1989; Hwang & Lee 2007), then increasing alignment of galaxies in clusters means an increase of angular momentum for the whole cluster. In the paper Godłowski et al. (2010), the orientation of galaxies was quantified by distribution of the angles. Specifically, the position angle of the galaxy plane p and two angles δ_d, giving spatial orientation of the normal to the galaxy plane, have been considered. The authors have also studied two additional angles. One is the angle between the normal to the galaxy plane and the main plane of the coordinate system. The second is the angle η between the projection of this normal onto the main plane and the direction toward the zero initial meridian (Flin & Godlowski 1986).

The entire range of all investigated angles was arranged into n bins. As we would like to detect a non-random effect in the orientation of galaxies, we first check whether the orientation is isotropic. To be specific, we check if the distribution of analyzed angles in the clusters under investigation is isotropic. The distribution of the above angles has been investigated using statistical tests. They were the χ² and Fourier tests. However, in the present paper we extend the analysis to first autocorrelation and Kolmogorov-Smirnov (K-S) tests (Hawley & Peebles 1975; Flin & Godlowski 1986; Godłowski et al. 2010; Godłowski 2012).

The χ² statistic is

$$\begin{eqnarray}&&{\chi }^{2}=\displaystyle \underset{k=1}{\overset{n}{\sum }}\displaystyle \frac{({N}_{k}-N\,{p}_{k}{)}^{2}}{N\,{p}_{k}}=\displaystyle \underset{k=1}{\overset{n}{\sum }}\displaystyle \frac{({N}_{k}-{N}_{0,k}{)}^{2}}{{N}_{0,k}},\end{eqnarray} \tag{ 16 }$$

where p_k is the probability that a chosen galaxy falls into the k-th bin, N is the total number of galaxies in a sample (in our case in a cluster), N_k is the number of galaxies within the k-th angular bin and ${N}_{0,k}$ is the expected number of galaxies in the k-th bin.

The first autocorrelation test quantifies the correlations between galaxy numbers in neighboring angle bins. The C statistic reads

$$\begin{eqnarray}&&C=\displaystyle \underset{k=1}{\overset{n}{\sum }}\displaystyle \frac{({N}_{k}-{N}_{0,k})({N}_{k+1}-{N}_{0,k+1})}{{\left[{N}_{0,k}{N}_{0,k+1}\right]}^{1/2}},\end{eqnarray} \tag{ 17 }$$

where ${N}_{n+1}={N}_{1}$.

If the deviation from isotropy is a slowly varying function of the analyzed angle θ, one can use the Fourier test

$$\begin{eqnarray}{N}_{k} & = & {N}_{0,k}(1+{{\rm{\Delta }}}_{11}\cos 2{\theta }_{k}+{{\rm{\Delta }}}_{21}\sin 2{\theta }_{k}\\ & & +{{\rm{\Delta }}}_{12}\cos 4{\theta }_{k}+{{\rm{\Delta }}}_{22}\sin 4{\theta }_{k}+\ldots ).\end{eqnarray} \tag{ 18 }$$

In this test, the crucial statistical quantities are amplitudes

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{1}={\left({{\rm{\Delta }}}_{11}^{2}+{{\rm{\Delta }}}_{21}^{2}\right)}^{1/2},\end{eqnarray} \tag{ 19 }$$

(only the first Fourier mode is taken into account) or

$$\begin{eqnarray}&&{\rm{\Delta }}={\left({{\rm{\Delta }}}_{11}^{2}+{{\rm{\Delta }}}_{21}^{2}+{{\rm{\Delta }}}_{12}^{2}+{{\rm{\Delta }}}_{22}^{2}\right)}^{1/2},\end{eqnarray} \tag{ 20 }$$

where the first and second Fourier modes are analyzed together. During our investigations we analyzed ${{\rm{\Delta }}}_{1}/\sigma ({{\rm{\Delta }}}_{1})$ and ${\rm{\Delta }}/\sigma ({\rm{\Delta }})$ statistics (see Godłowski et al. 2010 for details).

In the case of the K-S test we investigate the λ statistic

$$\begin{eqnarray}&&\lambda =\sqrt{N}\,{D}_{n}\end{eqnarray} \tag{ 21 }$$

which is given by the limit of the Kolmogorov distribution, where

$$\begin{eqnarray}&&{D}_{n}=\sup | F(x)-S(x)| \end{eqnarray} \tag{ 22 }$$

and F(x) and S(x) are theoretical and observational distributions of the investigated angle respectively.

The aim of the paper Godłowski et al. (2010) was to test the hypothesis that alignment of galaxies increases with cluster richness (Godłowski et al. 2005). The main result of Godłowski et al. (2010) was that the values of investigated statistics increase with increasing number of cluster galaxy members. This permits us to conclude that a relation exists between anisotropy and the number of galaxies in a cluster. Note that the above testing has been performed assuming a linear model $y={aN}+b$, where y is a value of the investigated statistic, N is the number of cluster members while a and b are linear regression coefficients. In the paper Godłowski et al. (2010), the null hypothesis H₀ (that the investigated statistic is random, i.e. neither increasing nor decreasing so that parameter $t=a/\sigma (a)=0$) has been confronted against hypothesis H₁ that a statistic increases with the cluster richness i.e. $t\gt 0$. In our previous paper (Stephanovich & Godłowski 2015), as well as in the present paper, we show that the dependence between alignment of galaxies in clusters and number of member galaxies is not necessarily linear, but could be according to the above assumptions of either ${N}^{4/3}$ or ${N}^{5/3}$.

For this reason, in Figure 2, we present statistics (χ², Fourier and first autocorrelation tests (Hawley & Peebles 1975; Flin & Godlowski 1986; Godłowski et al. 2010)) for the case of position angles obtained from the sample of 247 rich Abell clusters, analyzed by Godłowski et al. (2010). We present linear dependence $\sim N$ as well as the cases when analyzed statistics increase as ${N}^{4/3}$ and ${N}^{5/3}$. The error bars displayed in Figure 2 suggest that the data points seem to be not sufficient do discriminate between models. For this reason we analyze the dependence between the number of galaxies in a cluster and the value of analyzed statistics in more detail. We investigate the linear regression given by $y={aN}+b$ as indicated by various statistics. Namely, we have studied the linear regression between different statistics χ², ${{\rm{\Delta }}}_{1}/\sigma ({{\rm{\Delta }}}_{1})$, ${\rm{\Delta }}/\sigma ({\rm{\Delta }})$, C or λ and the number of analyzed galaxies in each particular cluster. This has been done for the case of linear dependence $\sim N$ or power laws $\sim {N}^{4/3}$ and $\sim {N}^{5/3}$ in the case of remaining models.

Zoom In Zoom Out Reset image size

Now we assume that the theoretical, uniform, random distribution contains the same number of clusters as the observed one. To be specific, we consider null hypothesis H₀ that the distribution is random and neither increases nor decreases. This means that the expected value of statistic $t=a/\sigma (a)=0$, while the t statistic has a Student's distribution with $u-2$ degrees of freedom, where u is the number of analyzed clusters. In other words, we test hypothesis H₀ that t = 0 against hypothesis H₁ that $t\gt 0$. Of course, in order to reject hypothesis H₀, the value of the observed t statistic should be greater than ${t}_{\mathrm{cr}}$ which we could obtain from statistical tables. Note that for our sample containing only 247 clusters, the critical value at the significance level α = 0.05 is equal to ${t}_{\mathrm{cr}}=1.651$.

The result of our statistical analysis is presented in Table 1. We analyzed two samples of data. In Sample A, all galaxies lying in the area regarded as a cluster were taken into account. In Sample B, to avoid "contamination" by background objects, we restrict ourselves by only considering galaxies brighter than ${m}_{3}+3$.

Table 1.  The $t=a/\sigma (a)$ statistic for 247 rich Abell clusters. Sample A means all galaxies while Sample B means galaxies brighter than ${m}_{3}+3$.

Test	N	${N}^{4/3}$	${N}^{5/3}$
Sample A
χ2	1.872	1.766	1.667
${{\rm{\Delta }}}_{1}/\sigma ({{\rm{\Delta }}}_{1})$	1.613	1.588	1.580
${\rm{\Delta }}/\sigma ({\rm{\Delta }})$	1.964	1.941	1.821
C	1.352	1.381	1.417
λ	2.366	2.500	2.400
Sample B
χ2	1,979	1.801	1.625
${{\rm{\Delta }}}_{1}/\sigma ({{\rm{\Delta }}}_{1})$	2.182	1.962	1.702
${\rm{\Delta }}/\sigma ({\rm{\Delta }})$	2.104	1.885	1.596
C	1.225	1.170	1.125
λ	2.421	2.000	1.765

Note that the cases of linear dependence for statistics χ², ${{\rm{\Delta }}}_{1}/\sigma ({{\rm{\Delta }}}_{1})$ and ${\rm{\Delta }}/\sigma ({\rm{\Delta }})$ have generally been analyzed in the paper Godłowski et al. (2010) (Table 1). Note, however, that our present results are somewhat different from those obtained by Godłowski et al. (2010). For example, in the case of χ² instead of t = 0.025/0.015 = 1.67 we obtain t = 1.87. The reason is that in the paper Godłowski et al. (2010) the error bars of data points (i.e. statistics for individual clusters) have been estimated from the sample, but now it is taken from the exact theoretical analysis (Godłowski 2012; Wang et al. 2003).

In the majority of cases, except for the first autocorrelation test, the values of obtained statistics are greater than the critical case of ${t}_{\mathrm{cr}}=1.651$. One could observe that for all three analyzed models (i.e. linear dependence $\sim N$ and the increased exponent like $\sim {N}^{4/3}$ or $\sim {N}^{5/3}$) we can eliminate hypothesis H₀ (that statistic $t=a/\sigma (a)=0$) in favor of hypothesis H₁ that $t\gt 0$. The effect increases if we analyze Sample B which means that we restrict the cluster membership to galaxies brighter than ${m}_{3}+3$. The significance of the effect decreases with increasing powers of m in models like N^m, but in the majority of cases the effect is significant. The above results allow us to conclude that the presented data are not sufficient to discriminate between the above three models, so we need future investigations based on larger cluster samples.

In our investigations, we have also studied time dependence of the pdf representing galaxy gravitational fields (Stephanovich & Godłowski 2015). The distribution function (11) evolves in time. It relies on explicit dependences ${f}_{1}(t)$ and ${f}_{2}(t)$. The functions ${f}_{1}(t)={a}^{2}(t)\dot{D}(t)$ while ${f}_{2}(t)=\dot{E}(t)$ (we use standard notations where dot means time derivative) could be obtained from the set of differential equations derived in the i-th order of perturbation theory by Bouchet et al. (1992)

$$\begin{eqnarray}&&{t}_{0}^{2}\ddot{D}(t)+a(t)D(t)=0,\end{eqnarray} \tag{ 23 }$$

$$\begin{eqnarray}&&{t}_{0}^{2}\ddot{E}(t)+a(t)E(t)=-a(t)D(t{)}^{2},\end{eqnarray} \tag{ 24 }$$

where $0\le t\lt \infty$ is dimensional physical time. The dimensionless function (scale factor) $a(t)$ is determined from the first Friedmann equation. In our investigations we consider the ΛCDM model, however we compare its predictions with those obtained in the classical CDM model.

To obtain the dependence ${L}_{0}(t)$, we use substitution $\lambda \to \lambda /{f}_{i}(\tau )$, ($\tau =t/{t}_{0}$) which yields (11)

$$\begin{eqnarray}&&H(\lambda,\tau )=\displaystyle \frac{2I(\lambda /{f}_{i}(\tau ))}{\pi \lambda },\quad i=1,2.\end{eqnarray} \tag{ 25 }$$

To derive ${f}_{\mathrm{1,2}}(\tau )$ in a particular model (the ΛCDM model in our case), it is necessary to calculate $a(t)$ from the first Friedmann equation, see Stephanovich & Godłowski (2015) for details

$$\begin{eqnarray}&&\displaystyle \frac{{da}}{{dt}}={H}_{0}\sqrt{{{\rm{\Omega }}}_{{\rm{\Lambda }}}{a}^{2}+\displaystyle \frac{1-{{\rm{\Omega }}}_{{\rm{\Lambda }}}}{a}}.\end{eqnarray} \tag{ 26 }$$

The solution of Equation (26) has the form

$$\begin{eqnarray}a(t) & = & \alpha {\sinh }^{2/3}(t/{t}_{0}),\\ \alpha & = & {\left(\displaystyle \frac{1-{{\rm{\Omega }}}_{{\rm{\Lambda }}}}{{{\rm{\Omega }}}_{{\rm{\Lambda }}}}\right)}^{1/3},\\ {t}_{0} & = & \displaystyle \frac{2}{3{H}_{0}\sqrt{{{\rm{\Omega }}}_{{\rm{\Lambda }}}}},\end{eqnarray} \tag{ 27 }$$

where ${{\rm{\Omega }}}_{{\rm{\Lambda }}}={\rm{\Lambda }}/(3{H}_{0}^{2})$ is the cosmological constant or so-called vacuum density, Λ is the cosmological constant and H₀ is the Hubble constant.

Having the function $a(t)$, we can solve Equation (23) numerically for $D(\tau )$ and then determine the function ${f}_{1}(\tau )={a}^{2}(\tau ){D}^{\prime}(\tau )\,({D}^{\prime}={dD}/d\tau )$. Accordingly, in the nonlinear regime, the function ${f}_{2}(\tau )={E}^{\prime}(\tau )$ could be calculated numerically from Equation (24).

One should note that functions ${f}_{2}(\tau )$, which are related to the second perturbative corrections, are negative. For instance, in the Einstein-de Sitter model ${f}_{1}(\tau )=(2/3)\tau$ and ${f}_{2}(\tau )=(-4/7){\tau }^{1/3}\lt 0$ (Doroshkevich 1970; Catelan &g Theuns 1996b). The same result (${f}_{2}(\tau )\lt 0$) can be obtained numerically for the ΛCDM model.

The dependences $H(\lambda,\tau )$ (25) for CDM (with above analytical expressions for ${f}_{i}(\tau )$) and ΛCDM models are shown in Figure 3. It is easy to observe that as time increases, the distribution function decreases, while its peak grows to infinity at $t\to 0$. As time grows, the whole distribution function "blurs" as its maximum shifts towards large t. It is also easy to notice that "blurring" of the distribution function at large times is much faster for the ΛCDM model. Also, both in linear and nonlinear regimes, $H(\lambda,\tau )$ increases with time. We emphasize once more that in the ΛCDM model this growth is much faster than in the CDM model. This is a consequence of the fact that functions ${f}_{i}(\tau )$ enter the exponent in the integrand (25). The comparison of the right and left panels of Figure 3 shows that the behavior of $H(\lambda,\tau )$ is qualitatively similar in linear and nonlinear regimes of fluctuation growth. This leads to the conclusion that even the linear regime gives a qualitatively correct approximation to the function $H(\lambda,\tau )$.

Zoom In Zoom Out Reset image size

The above results lead to the conclusion that angular momentum of galaxy clusters should increase in time. This hypothesis could be tested theoretically. This is because the limited speed of light causes the age of astronomical objects with different redshifts z to be different. So, assuming that galaxy clusters form in the same time instant, we expect that clusters with higher redshift z are younger. This means that galaxy alignment should decrease with z. Our preliminary analysis of the sample of 247 Abell cluster shows that in the case of χ² and K-S tests (Godłowski 2012; Aryal et al. 2013), the analyzed statistics decrease with z. Unfortunately this effect is not significant since parameter $t=a/\sigma (a)$ is less than 1.

One should note, however, that the basic catalog of galaxies is complete up to m = 18.3 mag, which means that the redshift of the most distant cluster $z\lt 0.12$. As a result, it is very difficult to detect such a subtle effect for a small cluster sample. Moreover, although the vast majority of clusters do not rotate (Regos & Geller 1989; Hwang & Lee 2007), this is not completely true for all clusters (Hwang & Lee 2007). Hwang & Lee (2007) study the dispersions and velocity gradients in 899 Abell clusters. They have found possible evidence for rotation in only six of them, i.e. less than 1%. A later sample of rotating clusters was studied by Aryal et al. (2013). The random orientation of galaxy angular momentum vectors in the analyzed clusters was found. Similarly, Yadav et al. (2017) found no preferred alignments of angular momenta vectors for galaxies in a sample of six dynamically unstable clusters. The presence of such cluster types, even relatively small, could introduce additional difficulties in observational investigation of the time evolution of cluster angular momenta. So, a larger sample of a cluster that extends to z is necessary to make an unambiguous conclusion regarding the above effect.

Our calculations demonstrate that although the gravitational interaction between stellar components (including dark matter halos) has a long-range multipole character, the observations (see below) give some confirmations that there is an additional short-range (like $\sim \exp (-r/{r}_{{\rm{c}}})$ with range ${r}_{{\rm{c}}}$) interaction. As a result, if the distance r between two objects (say galaxies) is smaller than ${r}_{{\rm{c}}}$, they are correlated, which means that their orbital momenta are aligned. This assumption works for dense (rich) galaxy clusters, which, by this virtue, have a high degree of orbital momentum alignment. For sparse (poor) clusters the situation is opposite. For such type of clusters with intergalaxy distance $r\gt {r}_{{\rm{c}}}$, the long-range multipole interaction prevails so that there is no alignment of the orbital momenta. The above statistical method accounts for this situation if we add the (empirical) short-range interaction term to the initial potential, Equation (2). In the analyzed case we obtain that the distribution function of random fields would depend on the average angular momentum ${L}_{\max }\equiv {L}_{\mathrm{av}}$ (see Stephanovich 1997; Semenov & Stephanovich 2002) and as a result we obtain the self-consistent equation for ${L}_{\mathrm{av}}$

$$\begin{eqnarray}&&{L}_{\mathrm{av}}=\displaystyle \int L(E)f(E,{L}_{\mathrm{av}}){d}^{3}E,\end{eqnarray} \tag{ 28 }$$

where $f(E,{L}_{\mathrm{av}})$ is the distribution function of gravitational field E, depending on ${L}_{\mathrm{av}}$ as a parameter. This function substitutes the expression in Equation (6) when a possible short-range interaction term is included. One should note that in the case of finite ${r}_{{\rm{c}}}$, the distribution function decays at $E\to \infty$ faster than in Equation (6) so that the integral in Equation (28) converges. As total interaction potential contains both luminous and dark matter components, Equation (28) allows us to ask the question about alignment of sub-dominant galaxies, even though the majority of cluster angular momentum is related to the smooth dark matter halo component. For instance, in the halo model (Schneider & Bridle 2010), when the luminous matter of galaxies is embedded in a dark matter halo, this halo by virtue of its mass may mediate the intergalaxy interaction, adding possible short-range terms to it. The self-consistent Equation (28) also permits including temperature into consideration (Semenov & Stephanovich 2003, 2002) and studying the time evolution of galaxies and their clusters (with respect to dark matter halos) within the ΛCDM model. Also, the combination of stochastic models (Garbaczewski & Stephanovich 2009; Garbaczewski et al. 2011) of primordial dynamics along with those of ΛCDM, most probably, would permit answering (at least theoretically) the question of whether the galaxies are initially aligned at the time of their formation, or if such an alignment is generated in some merger events, and how dark matter halos influence (mediate) this alignment.

Here we also show that there are different possible relations between angular momentum and the mass (richness) of the cluster. Note that the ${M}^{5/3}$ - power law for such a dependence as well as reasonable values for parameter λ in Equation (11) had been obtained by Heavens & Peacock (1988) followed by Schäfer & Merkel (2012).

Figure 2 reports our preliminary results of the dependence between analyzed statistics obtained for the sample of 247 rich Abell clusters (Godłowski et al. 2010). We conclude here that our comparison of the cases when the statistics grow as N, ${N}^{4/3}$ and ${N}^{5/3}$ does not permit establishing unambiguous correspondence of different dependences between angular momentum and richness of the structure. However, such unambiguous discrimination would be possible if a larger statistical sample of galaxy clusters is available. Moreover, we show that angular momentum of galaxies should increase with time. The latter fact follows from Equations (12) – (14) for the CDM model and from Figure 3 for the ΛCDM model. The physical mechanism of that has been discussed in details in our previous paper (Stephanovich & Godłowski 2015). It is related to the growing time evolution of scale factor $a(t)$ both in CDM (Doroshkevich 1970) and ΛCDM models, see Equation (27) for details. This means that the above theoretically predicted effect could be tested by observations as galaxy angular momentum should decrease with redshift z. Once more, the enlarged sample containing clusters with much higher z is necessary for such studies.

To summarize, in the present paper we analyze theoretically the observational dependences of the galaxies and their cluster angular momenta on their mass (richness). To do so, we use the method, introduced in our previous paper (Stephanovich & Godłowski 2015). Observational data are in agreement with our theoretical results and mainly Equations (15) and (13) where we have shown that under reasonable assumptions about cluster morphology the angular momentum of galaxy structures increases with their richness. The solution of Equation (28) will permit establishing a relation between the characteristics of possible short-range intergalaxy interaction and characterize their spin alignment.

We emphasize, however, that the above observational results about lack of alignment of galaxies for poor clusters, as well as evidence for such an alignment in rich galaxy clusters (Godłowski et al. 2005; Aryal et al. 2007, see also Godłowski 2011 for later improved analysis) clearly show that angular momentum of galaxy groups and clusters increases with their richness. The problem of cluster angular momenta in the context of their mutual interactions as well as those with dark matter halos has been discussed by Hahn et al. (2007) based on the results of computer simulations. The presence of a threshold value for cluster mass (that is to say richness) has been noticed in these simulations. This threshold value is related to mutual alignment of clusters and dark matter halo axes. As we have shown above, this fact can be explained by our model.

We finally note that direct computer simulations of stellar ensembles are still quite computationally expensive to simulate realistic (i.e. sufficiently large) parts of the Universe. Hence it seems to be a good idea to put some effort into developing new theoretical models for galaxy alignment with respect to dark matter halos and (possible) merger into larger structures like superclusters. Since galaxy morphology plays an important role in this behavior, our approach, linking the galaxy shapes with their characteristic distribution (especially in view that it permits calculating non-Gaussian pdfs), will improve the overall understanding, which can additionally be tested against observed galaxy shape distributions and alignments.

Here we present some more details about our model, based on the Hamiltonian (1). In this Hamiltonian, the explicit expression for the i-th galaxy quadrupolar moment Q_i has the form (Poisson 1998)

$$\begin{eqnarray}&&{Q}_{i}=\displaystyle {\int }_{{V}_{i}}{\rho }_{i}({\boldsymbol{x}})|{\boldsymbol{x}}{|}^{2}{P}_{2}({\boldsymbol{s}}\cdot {\boldsymbol{x}}){d}^{3}x,\end{eqnarray} \tag{ A.1 }$$

where ${P}_{2}(x)=(3{x}^{2}-1)/2$ is the corresponding Legendre polynomial (Abramovits & Stigan 1979), V_i is a volume of the i-th galaxy and ${\rho }_{i}({\boldsymbol{x}})$ is its mass density.

The geometry of the problem under consideration is shown in Figure A.1. It is seen first that the origin is not related to any specific galaxy or other astronomical object. Rather, it is situated at an arbitrary point in the Universe. Although ${{\boldsymbol{r}}}_{{ij}}$ is directed from one galaxy (in our case j) towards another (in our case i), it is by no means bound to these galaxies. It simply means the difference in their radius vectors, which connect the coordinate origin and position of each galaxy.

Zoom In Zoom Out Reset image size

The Hamiltonian in Equation (1) can be identically rewritten through the interaction energy

$$\begin{eqnarray}{\mathcal{H}} & = & -{{GM}}^{2}\displaystyle \mathop{\sum }\limits_{i}{p}_{i}{m}_{i}{W}_{i},\\ {W}_{i} & = & W({{\boldsymbol{r}}}_{i})=\displaystyle \mathop{\sum }\limits_{j}{m}_{j}V({{\boldsymbol{r}}}_{{ij}})\\ & \equiv & \displaystyle \mathop{\sum }\limits_{j}{m}_{j}V({{\boldsymbol{r}}}_{j}-{{\boldsymbol{r}}}_{i}).\end{eqnarray} \tag{ A.2 }$$

The interaction energy W_i is the energy exerted by the rest of the galaxy ensemble (due to intergalaxy interaction) on the galaxy at point i. We can see that after summation (actually integration, see below) over ${{\boldsymbol{r}}}_{j}$ the relative intergalaxy distance ${{\boldsymbol{r}}}_{{ij}}$ has actually disappeared.

The gradient of energy Equation (A.2) is indeed the gravity field, which acts on the i-th galaxy (or other astronomical object) from the other j objects in the ensemble

$$\begin{eqnarray}{{\boldsymbol{E}}}_{\mathrm{quad}}({{\boldsymbol{r}}}_{i}) & \equiv & {{\boldsymbol{E}}}_{\mathrm{quad},i}=\displaystyle \mathop{\sum }\limits_{j}{m}_{j}{\rm{\nabla }}V({{\boldsymbol{r}}}_{j}-{{\boldsymbol{r}}}_{i})\\ & = & {{\boldsymbol{i}}}_{r}{E}_{0}\displaystyle \mathop{\sum }\limits_{j}{m}_{j}\displaystyle \frac{3\,{\cos }^{2}{\theta }_{{ij}}-1}{{r}_{{ij}}^{4}},\end{eqnarray} \tag{ A.3 }$$

which is the expression in Equation (2) from the text, rewritten explicitly in terms of vectors ${{\boldsymbol{r}}}_{i}$ and ${{\boldsymbol{r}}}_{j}$.

Having the expression in Equation (A.3), we can explicitly write the distribution function of random quadrupolar fields, Equation (3) from the text

$$\begin{eqnarray}f({\boldsymbol{E}}) & = & \overline{\delta ({\boldsymbol{E}}-{{\boldsymbol{E}}}_{i})}\equiv \overline{\delta ({\boldsymbol{E}}-{{\boldsymbol{E}}}_{\mathrm{quad}}({{\boldsymbol{r}}}_{i}))}\\ & = & \overline{\delta \left({\boldsymbol{E}}-{{\boldsymbol{i}}}_{r}{E}_{0}\displaystyle \mathop{\sum }\limits_{j}{m}_{j}\displaystyle \frac{3\,{\cos }^{2}{\theta }_{{ij}}-1}{{r}_{{ij}}^{4}}\right)},\end{eqnarray} \tag{ A.4 }$$

where bar means averaging over random spatial configurations of galaxies and other astronomical objects.

In performing the actual averagings in the expression Equation (A.4) (see Fig. A.2), with respect to the fact that the number of galaxies is infinite and their "elementary separations" ${{\boldsymbol{r}}}_{{ij}}$ become very small, we can change summations in Equation (A.3) and Equation (A.4) to integrations using the expression for gravity field ${{\boldsymbol{E}}}_{i}$ in the form given by Equation (2). Further averagings in Equation (A.4) are prescribed in the text, see also Stephanovich & Godłowski (2015).

Zoom In Zoom Out Reset image size