text:
# 1 Introduction - the Rotating Mass Matrix

## 1 Introduction - the Rotating Mass Matrix

By a rotating mass matrix, we mean one which undergoes unitary transformations through scale changes.

As most quantities in quantum field theory vary under changing scales, in particular the masses of particles or in other words the eigenvalues of mass matrices, it is no surprise that the orientations of mass matrices may also change.

Indeed, even in the Standard Model as usually conceived, nontrivial mixing between up and down fermion states will generally induce rotating mass matrices. That this is so can be seen from the following renormalization group equations:

$$16\pi ^2\frac{dU}{dt}=\frac{3}{2}(UU^{}DD^{})U+(\mathrm{\Sigma }_uA_u)U,$$
(1.1)
$$16\pi ^2\frac{dD}{dt}=\frac{3}{2}(DD^{}UU^{})D+(\mathrm{\Sigma }_dA_d)D.$$
(1.2)
satisfied by the mass matrices $`U`$ and $`D`$ of respectively the up and down fermions . Nontrivial mixing between up- and down-states means that the matrices $`U`$ and $`D`$ are related by a nondiagonal matrix $`D=VUV^{}`$ or $`U=V^{}DV`$, so that $`U`$ and $`D`$ cannot be simultaneously digonalized. Suppose now at some $`t`$ we diagonalize $`U`$ in (1.1) by an appropriate unitary transformation. All other terms in (1.1) are then diagonal but not the term $`DD^{}U`$. Hence, this term will necessarily de-diagonalize the matrix $`U`$ on running by (1.1) to a different $`t`$, or in other words the matrix $`U`$ will rotate with changing $`t`$ . Similar arguments hold also for the rotation of the matrix $`D`$. Given that quarks have long been known to have a nontrivial (CKM) mixing matrix , and that recent experiments on neutrino oscillation  strongly suggest a nontrivial (MNS) mixing  also for leptons, we have to conclude that the mass matrices for both quarks and leptons will rotate with changing scales even if there are no forces and interactions in nature other than those currently studied in the Standard Model.

However, there is clearly a possibility—perhaps one might even say a strong theoretical reason to suppose—that forces may exist in nature other than those now conventionally studied in the Standard Model, which can give further rotations to the fermion mass matrices. Indeed, in the usual formulation of the Standard Model, the fact that there are 3 generations of fermions and that they mix is taken as an input from experiment, with the result that the framework depends on a large number of empirical parameters, about three-quarters of which are traceable to the mystery of generations. One could thus hope that in future when the generation puzzle is solved, some or even most of the empirical parameters appearing in the present Standard Model will be predictable. One favourite and perhaps most natural assumption for theoretical attempts in this direction is that generations originate as a broken ‘horizontal’ gauge symmetry . If that is the case, then the interactions associated with this symmetry which mix the generation index are likely also to rotate the fermion mass matrices in generation space. This rotation would be over and above that driven by the nondiagonal CKM or MNS matrices via the mechanism described in the preceding paragraph.

Whatever its origin, however, the rotation of fermion mass matrices is theoretically of the same standing as the running of masses and coupling constants, and if testable by experiment would lend equal support to the quantum field theory framework. Moreover, if experiment shows that the rotation is different from simply that given by renormalization group equations (1.1) and (1.2) via the empirical mixing matrix, then it suggests that there are forces at work other than those currently studied in the Standard Model which are of the type considered in the preceding paragraph. This means that we would have opened up a new window for investigations furthering our long-standing quest for the origin of fermion generations. And since, as we shall see, mass matrix rotation leads to an entirely new and very distinctive category of physical phenomena, there is a fair chance of its effects being observed by experiment in the not too distant future. The purpose of the present paper is to explore such possibilities.

## 2 Specification of Fermion States

The fact that the mass matrix rotates poses immediately a question of physical interpetation at a rather basic level, namely the question of how to define the state vectors of the various fermion states, and hence the mixing matrices between them. For the familiar case of a non-rotating mass matrix, one defines the state vectors of the 3 generations as its eigenvectors, which are by assumption scale-independent and, since the matrix is hermitian, are also mutually orthogonal as they should be if they are to represent independent physical entities. The mixing matrix between the up and down fermion types, being the overlap matrix between 2 scale-independent orthonormal triads of state vectors, one for each type, is then also scale-independent and automatically unitary. Only the eigenvalues of a mass matrix then run, and the masses of the 3 states (generations) can be defined as the running values each evaluated at the scale equal to its value. However, when the mass matrix rotates, although one can still diagonalize the matrix at any scale by a triad of mutually orthogonal eigenvectors, this triad will be scale-dependent, and it is not immediately clear which vectors at which scale(s) are to be identified as the state vectors of the 3 generations. Although the ambiguity may not in all cases be numerically significant in view of present limitations in experimental accuracy, it still has to be resolved as a matter of principle against the day when better accuracy is achieved.

Now, the state vector of a physical state like the muon is normally, we think, taken as a scale-independent concept, or otherwise we are likely to meet with some awkwardness in its physical interpretation. For example, suppose we were to define the muon state as the eigenvector with the second highest eigenvalue of the charged lepton mass matrix at any scale. Then since the matrix rotates, the $`\mu `$ state vector will point at different scales in different directions in generation space. In that case, the muons obtained say from an energetic beam of charged $`\pi `$’s, which we know decay almost entirely into $`\mu `$’s at rest, will no longer appear as purely muons when it hits the target, but as a linear combination of $`e,\mu `$ and $`\tau `$, which is at variance with what is usually understood. Although one can in principle insist on defining fermion states as the eigenstates of the rotating mass matrix at every scale, and so long as one is consistent in their interpretation one would arrive in the end at the same physical results, such a procedure would seem to be rather inconvenient. In this paper therefore, we opt for a scale-independent definition of fermion state vectors, in which case it would be incumbent upon us to specify exactly the scale at which each of these vectors are to be defined.

If the rotation of the mass matrix is negligible to a certain approximation below some scale, then one may define to this approximation the state vectors as the orthonormal triad of eigenvetors of the mass matrix taken at that scale. This, though seldom stated explicitly is, we think, the tacit criterion usually adopted in the literature which we shall refer to in this paper as Fixed Scale Diagonalization FSD. It is, of course, at best an approximation for the mass matrix cannot stop rotating abruptly, and that approximation, as we shall see later, may not be all that good in certain circumstances. Besides, this criterion involves a degree of arbitrariness in choosing a certain fixed scale to effect the diagonalisation of the mass matrix, and begs the question of principle why it should be that particular scale and not some other.

In view of this, one may be tempted instead to do the following. As one defines the mass of a state as the eigenvalue evaluated at the scale equal to its value, one may try to define as its state vector the corresponding eigenvector also at the same scale. This would be more democratic and less arbitray than the FSD criterion but it will not work. Consider for example the $`U`$-type quarks. Following the above proposal, one would then define the $`t`$ state vector as the eigenvector with the largest eigenvalue $`m_3`$ of the matrix evaluated at the scale $`\mu =m_3(\mu )`$, the $`c`$ state vector as the eigenvector with the second largest eigenvalue $`m_2`$ at the scale $`\mu =m_2(\mu )`$, and the $`u`$ state vector as the eigenvector with the smallest eigenvalue $`m_1`$ at the scale $`\mu =m_1(\mu )`$. These 3 state vectors, however, will not be mutually orthogonal, for although the 3 eigenvectors of the hermitian mass matrix are mutually orthogonal when all evaluated at the same scale, they have no reason to be so when evaluated each at a vastly different scale since the matrix by assumption rotates with changing scales. It thus contradicts the assertion that they represent 3 independent physical entities; in particular, it would imply that the mixing (CKM or MNS) matrix, being a transformation matrix between the 2 triads of up and down state vectors, is not unitary as it ought to be.

There is, however, a working criterion for defining state vectors each at its own mass scale along the lines suggested in the above paragraph which takes account of the rotation between the different mass scales yet still gives mutually orthogonal state vectors and hence a unitary mixing matrix. This was first proposed by us  in connection with a scheme for quark and lepton mixing which we call the Dualized Standard Model (DSM), but the criterion can in fact be applied to define state vectors in any scheme with a rotating mass matrix . For 3 generations, the criterion goes as follow. We run the mass matrix $`m`$ down in scale until we have for its highest eigenvalue $`m_3`$ a solution to the equation $`\mu =m_3(\mu )`$. This value at this scale we define as the mass $`m_3`$, and the corresponding eigenvector the state vector $`𝐯_\mathrm{𝟑}`$ of the heaviest generation. Below that energy, the state 3 no longer exists as a physical state, and only the two lighter generations survive, whose state vectors have to be orthogonal to $`𝐯_\mathrm{𝟑}`$. We define then the mass matrix at energies below $`m_3`$ as the $`2\times 2`$ submatrix $`\widehat{m}`$ of $`m`$ in the subspace orthogonal to $`𝐯_\mathrm{𝟑}`$. To find now the mass and state vector for generation 2, we follow with $`\widehat{m}`$ the same procedure used above with $`m`$ for generation 3 and run $`\widehat{m}`$ down in scale until we find a solution to the equation $`\mu =\widehat{m}_2(\mu )`$, which value we call the mass $`m_2`$ and the corresponding eigenvector at that scale the state vector $`𝐯_\mathrm{𝟐}`$ of the generation 2. The state vector of the lightest generation 1 is now also defined as the vector orthogonal to both $`𝐯_\mathrm{𝟑}`$ and $`𝐯_\mathrm{𝟐}`$, while the mass of 1 is obtained by repeating the same procedure, namely by running down in scale the remaining $`1\times 1`$ submatrix, namely the expectation value $`𝐯_\mathrm{𝟏}|m|𝐯_\mathrm{𝟏}`$, until its value equals the scale. In this way, each mass is evaluated at its own appropriate scale while the physical state vectors of the 3 generations are mutually orthogonal as they should be. Furthermore, the mixing matrix, taken as the overlap matrix between the triads of state vectors so defined for the up and down states, will be both unitary and scale-independent. This criterion for defining state vectors, which is the only one we are aware of with the required properties, we shall refer to in future as Step-by-Step Diagonalization SSD.

## 3 Fermion Transmutation as Consequence

Whether one accepts the SSD proposal for defining physical fermion state vectors or chooses to ignore these subtleties and adopts the approximate criterion of FSD, the fact remains that the mass matrix will continue to rotate at high energy, and so long as the state vectors are defined as the eigenvectors diagonalizing the mass matrix at some prescribed scale(s), they will generally no longer diagonalize the mass matrix at some higher energy scales. Now, the fermion mass matrix appears for example in the fermion propagator of Feynman diagrams, and if it is nondiagonal, so also will be the transition amplitudes represented by these diagrams. Physically, this means that transitions will be induced in which a particle of one generation converts into a particle of another generation. For example,

the diagram of Figure 1, where a tilde denotes the occurrence of a nondiagonal fermion mass matrix, will have an off-diagonal element linking the $`e`$ state vector to the $`\mu `$ state vector, leading to a $`e^+e^{}\mu ^+e^{}`$ transition.

Such a conversion of $`e`$ to $`\mu `$ is quite distinctive and singular, differing from that induced, for example, by a flavour-changing neutral current (FCNC) exchanged between the initial $`e^+`$ and $`e^{}`$. The latter will normally produce in addition to the $`\mu ^+`$ at the upper vertex also a $`\mu ^{}`$ at the lower vertex, giving rise to a “double conversion”, although under certain circumstances where the exchanged FCNC bosons themselves mix, “single conversion” can also occur, as for example in the situation considered in . However, even in that last case, the conversion rate depends on the masses and couplings of the exchanged FCNC bosons which are still rather elusive quantitites, whereas the conversion effect considered here is due only to the rotation of the mass matrix, and may therefore be quite restricted. We propose therefore to call the latter “transmutation” to distinguish it from other conversion phenomena.

The effect of transmutation is so unique as a consequence of the rotating mass matrix that if it occurs with appreciable cross sections, it should be identifiable in experiment with little difficulty. At first sight, therefore, the prediction looks alarming, for the conservation of such quantities as muon number have been checked already to high accuracy, and their wholesale violation as suggested above seems hardly likely to survive the existing experimental bounds. On closer examination, however, the conclusion becomes less obvious, mainly for the following reason. At low energy, near the scale(s) at which the mass matrix is diagonalized to define the physical fermion states, as explained in the preceding section, the deviations of the mass matrix from diagonality are small, and so also will be the cross sections for transmutation. At high energy, on the other hand, where the off-diagonal mass matrix elements can become sizeable, it turns out that transmutation cross sections will be suppressed, as can be seen in the following example.

One of the simplest transmutation processes that one can think of is what one can call the photo-transmutation of leptons, which to leading order is given just by the Feynman diagrams in Figure 2

which are the same as for ordinary Compton scattering except that the fermion line here carries a rotating and generally non-diagonal mass matrix. Consider first the $`s`$-channel diagram (a) in which the lepton propagator $`1/(p\text{/}+k\text{/}m)`$ appears, carrying the non-diagonal mass matrix $`m`$ where transmutation originates. For $`s=(p+k)^2`$ large, we can expand the propagator in powers of $`m/\sqrt{s}`$, obtaining symbolically:

$$\frac{1}{(p\text{/}+k\text{/}m)}\frac{1}{(p\text{/}+k\text{/})}(1+\frac{m}{\sqrt{s}}),$$
(3.1)
where one sees that the off-diagonal terms leading to transmutation is suppressed by a factor of $`1/\sqrt{s}`$. Similarly, the lepton propagator in the diagram (b) for $`u`$ large can be approximated as:

$$\frac{1}{(p\text{/}k\text{/}^{}m)}\frac{1}{(p\text{/}k\text{/}^{})}(1+\frac{m}{\sqrt{|u|}}),$$
(3.2)
leading to transmutation effects suppressed by a factor $`1/\sqrt{|u|}`$. In the region where $`|u|`$ is small, off-diagonal contributions can remain sizeable, but the size of this region itself decreases as $`1/s`$ so that the integrated effect remains small.

For these reasons, transmutation effects arising from the rotating mass matrix are seldom very large and the question whether they can already be ruled out by existing data, and if not, how soon, and under what conditions, they will be verifiable, is not immediately obvious and can only be answered by detailed analysis. The answers will depend on the process investigated and on the speed with which the mass matrices rotate, and these in turn depend on the theory or model one is considering. We examine below and in a companion paper  some specific examples.

## 4 Two Rotation Mechanisms

In this section, we investigate the rotating fermion mass matrix in 2 specific schemes at, as it were, opposite ends of the spectrum so as hopefully to span most reasonable possibilities. At one end of the spectrum, we assume that there are no forces in nature other than those currently studied in the Standard Model so that all the rotation in the fermion mass matrix is driven by the electroweak Higgs terms via the empirical nondiagonal mixing (CKM or MNS) matrix. The generation phenomenon is thereby left unexplained. Fermion state vectors are to be specified by the FSD prescription (to contrast with the scheme below which only works with the SSD prescription) in which the mass matrix is diagonalized at some chosen (low) scale. In the case of leptons, since some mass and mixing parameters are but poorly known, these specifications have yet to be supplemented by further assumptions based on current popular views in a manner to be detailed later. We shall refer to this henceforth as the NSM scheme (for naively implemented Standard Model). At the other end of the spectrum, we consider a scheme we ourselves suggested , which in turn was abstracted from our so-called Dualized Standard Model (DSM)  although it need not strictly adhere to the duality concepts of that earlier work. In this scheme, generation itself originates from a broken $`SU(3)`$ gauge symmetry, the rotation is driven by Higgs bosons associated with the breaking of this symmetry, and even the nondiagonal mixing matrices as experimentally observed are themselves consequences of the mass matrix rotation. Here the SSD prescription for defining state vectors is essential. Given that in the DSM it is the rotation which gives rise to fermion mixing, and not the other way round as in the NSM case with the mixing driving the rotation, transmutation effects are usually, though not in all cases, larger in the DSM than in the NSM, and hence would be easier either to disprove or to confirm.

Our first task is to evaluate the rotating mass matrix over a range of scales, say, up to the energy achievable in the foreseeable future of around 100 TeV. For the NSM case, this in principle requires the solution of the equations (1.1) and (1.2) together with any other renormalization group equations they are coupled to, which though standard is a little involved. For our present purpose of a first exploration, however, we may approximate by linearizing the equations and keeping only those terms actually driving the rotation, namely the term $`DD^{}U`$ in (1.1) and the term $`UU^{}D`$ in (1.2). We need as input from experiment the masses of the fermion states and the mixing matrices measured at low energies. Assuming that the scale in the FSD prescription chosen for defining the fermion state vectors coincides with the scale at which the empirical mixing matrices are measured, the rotating mass matrices can then be evaluated by iterating the two equations. For quarks, both the masses and the CKM matrix have now been determined to sufficient accuracy for our purpose and pose no practical difficulty. This calculation has been performed but since it will not be of use in what follows in this paper, the result will not be presented. For leptons, however, the MNS matrix is still poorly known, while the Dirac masses of neutrinos on which their Higgs couplings depend are still almost entirely unconstrained. For these, therefore, we shall just insert for our exploration the most popular theoretical biases at present, namely for the MNS matrix the so-called bi-maximal mixing  version:

$$U_{MNS}=\left(\begin{array}{ccc}1/\sqrt{2}& 1/\sqrt{2}& 0\\ 1/2& 1/2& 1/\sqrt{2}\\ 1/2& 1/2& 1/\sqrt{2}\end{array}\right),$$
(4.1)
and for the Dirac mass $`m_3`$ of the heaviest neutrino a value around the mass of the $`t`$ quark. One shall need the Dirac mass $`m_2`$ of the second heaviest neutrino also, which is taken, for lack of any better choice, to be of the order of the charm mass. Since the equations have been linearized in approximation, the result of the rotation for charged leptons can be summarized simply by giving the rates of change of the off-diagonal elements of the (‘hermitized’) mass matrix with the logarithm of the energy, which are explicitly given as:

$$\frac{dL}{dt}=\frac{3}{128\pi ^2}\frac{1}{V^2}\left(\begin{array}{ccc}2m_2^2m_e& \frac{1}{\sqrt{2}}m_2^2(m_e+m_\mu )& \frac{1}{\sqrt{2}}m_2^2(m_e+m_\tau )\\ \frac{1}{\sqrt{2}}m_2^2(m_e+m_\mu )& 2m_3^2m_\mu & m_3^2(m_\mu +m_\tau )\\ \frac{1}{\sqrt{2}}m_2^2(m_e+m_\tau )& m_3^2(m_\mu +m_\tau )& 2m_3^2m_\tau \end{array}\right),$$
(4.2)
with $`V`$ being the vev of the electroweak Higgs, i.e. 246 GeV. Putting in then the suggested values of the masses, one obtains that $`\mu |m|\tau `$ changes by about $`5.5\times 10^3`$, $`e|m|\tau `$ by about $`1.77\times 10^7`$ and $`e|m|\mu `$ by about $`1.06\times 10^8`$ GeV per decade change in energy. These results will be useful later for estimating the rates of various transmutation effects implied by the NSM scheme. One notices that in (4.2), $`e|m|\tau `$ is proportional to $`m_2^2`$ only by virtue of the assumed zero in the bi-maximal mixing matrix. If that element $`U_{e3}`$ of the MNS matrix is nonzero, then $`e|m|\tau `$ would be proportional to $`m_3^2`$, which in turn can lead on iteration to terms of order $`m_3^4`$ in the element $`e|m|\mu `$, thus making both much larger than the above estimates, which represent therefore a sort of lower limit.

For the DSM case, the calculation, though more intricate, has been more explicitly specified. As the scheme depends on only 3 parameters  which have already been determined by fitting to some and then used to give very sensible predictions to others of the mass and mixing parameters of both quarks and leptons, no further information needs to be supplied. Indeed, in making these fits and predictions, the rotating mass matrices for quarks and leptons had perforce already been evaluated, though not previously presented. The result for the charged leptons which is of concern to us in this paper is now given in Figures 3.

To check our result, without going through the details of the fit in , one needs only to iterate the following equation for a vector $`(x^{},y^{},z^{})`$:

$$\frac{d}{dt}\left(\begin{array}{c}x^{}\\ y^{}\\ z^{}\end{array}\right)=\frac{5}{32\pi ^2}\rho ^2\left(\begin{array}{c}x_1^{}\\ y_1^{}\\ z_1^{}\end{array}\right),$$
(4.3)
with

$$x_1^{}=\frac{x^{}(x^2y^2)}{x^2+y^2}+\frac{x^{}(x^2z^2)}{x^2+z^2},\mathrm{cyclic},$$
(4.4)
starting from some initial value $`(x_I,y_I,z_I)`$ with a coupling $`\rho `$, all given in eq. 16 of . For an accuracy of about a percent, some 500 iterations are needed per decade change in energy, where the vector $`(x^{},y^{},z^{})`$ is to be normalized after every iteration. The resulting mass matrix is of the factorized form:

$$m^{}=m_T\left(\begin{array}{c}x^{}\\ y^{}\\ z^{}\end{array}\right)(x^{},y^{},z^{}),$$
(4.5)
where $`m_T`$ is a normalization depending on the fermion-type $`T`$, which may be taken as the mass of the heaviest generation state of that type. In (4.5), the mass matrix is given with respect to a fixed (gauge) basis in generation space, but the matrix with respect to the physical basis as presented in, for example, Figure 3 is also easily evaluated with the (physical) state vectors determined in accordance with the SSD prescription outlined in section 2. Thus, the physical state vectors of the charged leptons in the gauge basis are, explicitly:

$`𝐯_\tau `$ $`=`$ $`(0.9967,0.0760,0.0268)`$
$`𝐯_\mu `$ $`=`$ $`(0.0759,0.7741,0.6285)`$
$`𝐯_e`$ $`=`$ $`(0.0271,0.6285,0.7774).`$ (4.6)
For more details, the reader is referred to  and earlier references therein.

One particular feature of the DSM scheme is that all fermions lie on the same trajectory of the rotating vector $`(x^{},y^{},z^{})`$ so that at any scale the mass matrix for up and down fermion types have the same eigenvectors. This means that the $`U`$ and $`D`$ matrices in (1.1) and (1.2) can be simultaneously diagonalized, causing thus no de-diagonalization to each other via the electroweak Higgs loop effect. We recall that it was the latter which was what drove the rotation in the NSM case. Here instead all the rotation comes from the Higgs bosons associated with the broken ‘horizontal’ symmetry.

Comparing the actual values for the rotating mass matrices of the two schemes NSM and DSM, one finds the off-diagonal elements to be generally larger for the latter than the former, as we expected. Further, as can be seen in Figure 3, the mass matrix in the DSM scheme has quite an intricate structure with, for example, the off-diagonal elements passing through zero whenever the energy scale equals the mass of a physical fermion state, meaning that the mass matrix is there diagonal. This is a consequence of the SSD prescription by construction, and would apply also to the NSM had one chosen there also the SSD instead of the FSD prescription.

## 5 Examples of Transmutational Decays

Once given the rotating mass matrix as a function of the energy scale, it would appear that one has in principle already enough information to evaluate any implied transmutational effect. However, because of the mass matrix rotation leading to, among other things, some rather unfamiliar kinematics, adaptations to usual calculational procedures are required, which are not, at first encounter, entirely straight-forward. We find that this is true even in evaluating, for example, the simple Feynman diagrams of Figure 2 for the photo-transmutation of leptons that we have calculated. Therefore, in order to avoid clouding immediately the essentially simple general picture of transmutation with technical details, we relegate detailed calculations to other papers, e.g. on photo-transmutation to a companion paper , and limit ourselves here to qualitative estimates of widths for some illustrative examples in transmutational decays, which being single particle effects are simpler to analyse theoretically. Also, to avoid the complications of quark confinement, we shall deal here with only lepton transmutations, which are already of considerable physical interest, being constrained by some very stringent bounds in experiment.

(I) As a first example, consider the decays:

$$Z^0l_\alpha \overline{l}_\beta ,$$
(5.1)
with $`l_\alpha l_\beta `$ being different lepton states. As indicated already in Section 3, the amplitudes of such off-diagonal transmutational processes are expected to be suppressed with respect to those of their diagonal counterparts $`Z^0l_\alpha \overline{l}_\alpha `$ by a factor of the form $`\alpha |m|\beta /E`$, with $`E`$ being a measure of the typical energy carried by the fermion line, which we can take here to be of the order of (half) the $`Z^0`$ mass $`M_Z`$. The numerator $`\alpha |m|\beta `$ represents the off-diagonal mass matrix element linking the lepton states $`l_\alpha `$ and $`l_\beta `$ at the scale of the reaction, which is again $`M_Z`$ in the present case. Given any choice for $`l_\alpha `$ and $`l_\beta `$, $`\alpha |m|\beta `$ can be read off for the DSM in Figure 3. Thus, for example, for the decay $`Z^0\tau ^{}\mu ^+`$ the relevant matrix element $`\mu |m|\tau `$ at the scale $`M_Z=91`$GeV is seen to have value about $`0.11`$ GeV. For the NSM, the element $`\mu |m|\tau `$ was said in Section 3 to vary by about 0.5 percent per decade change in energy. According to the FSD prescription adopted there, we have to first fix a scale at which the mass matrix is diagonal. If this is taken, for lack of any obvious better choice, to be at the $`\tau `$ mass, we obtain for $`\mu |m|\tau `$ at $`M_Z`$ a value of about 0.01 GeV. To estimate the width for the decay $`Z^0\tau ^{}\mu ^+`$, the easiest way would be to compare it with either of the diagonal decays $`Z^0\tau ^{}\tau ^+`$ or $`Z^0\mu ^{}\mu ^+`$. Given that the masses of the final state leptons are so much smaller than the decaying $`Z^0`$ mass, one can neglect the small differences in kinematics between the different decay modes and write simply:

$$\frac{\mathrm{\Gamma }(Z^0\tau ^{}\mu ^+)}{\mathrm{\Gamma }(Z^0\tau ^{}\tau ^+)}\frac{|\mu |m|\tau |^2}{M_Z^2}.$$
(5.2)
Putting in the above estimates of the matrix elements and the empirical value 3.36 percent  for the branching ratio of the diagonal $`\tau ^{}\tau ^+`$ mode, one obtains then for the branching ratio of the transmutational mode $`\tau ^{}\mu ^+`$ an estimate of about $`4\times 10^8`$ for DSM and $`4\times 10^{10}`$ for NSM. Both of these estimates are much below the present empirical bound of $`1.2\times 10^5`$ listed in . Similarly, estimates for the branching ratios of the other modes are obtained, giving for the $`\tau ^{}e^+`$ mode a branching ratio of $`4\times 19^9`$ and for the $`\mu ^{}e^+`$ mode $`1.6\times 10^{11}`$ for DSM. Again these are way below the present empirical limits listed in  of respectively $`9.8\times 10^6`$ and $`1.7\times 10^6`$. The estimates in the NSM for these last two modes are minuscule due to the very small values for the mass matrix elements linking $`e`$ to $`\mu `$ and $`e`$ to $`\tau `$, as will also be true for the following examples. They shall therefore henceforth be ignored except, for a particular reason to be explained, for (IV) below. As one can easily see, the reason that the transmutation rates are all so small in $`Z^0`$ decay is that the off-diagonal effects are much suppressed by the high energies involved which is here of the order of the $`Z^0`$ mass.

(II) As a second example then, let us go to the other extreme in energy and consider the decay:

$$\pi ^0\mu ^{}e^+$$
(5.3)
which in $`\pi ^0`$-decay is the only kinematically accessible transmutational leptonic mode. The same arguments as before suggests that we compare this with the diagonal mode $`\pi ^0e^{}e^+`$ giving a suppression in the amplitude of the former with respect to the latter by a factor of the form $`e|m|\mu /m_\pi `$ where the off-diagonal mass matrix element is to be evaulated at the scale $`m_\pi `$. Here, the spins, parities and masses of the particles involved being such, a correction for the difference in kinematics between the two decays is warranted, which is made most easily by comparison with the analogous charged $`\pi `$ decays, giving thus:

$$\frac{\mathrm{\Gamma }(\pi ^0\mu ^{}e^+)}{\mathrm{\Gamma }(\pi ^0e^{}e^+)}\frac{|e|m|\mu |^2}{m_\pi ^2}\frac{\mathrm{\Gamma }(\pi ^\pm \mu ^\pm \nu _\mu )}{\mathrm{\Gamma }(\pi ^\pm e^\pm \nu _e)}.$$
(5.4)
The matrix element $`e|m|\mu `$ at scale $`m_{\pi ^0}=135`$ MeV obtained from the calculation giving Figure 3 is about $`3\times 10^4`$ GeV. Hence, putting in the empirical value of $`7.5\times 10^8`$ for the branching ratio of the diagonal mode $`\pi ^0e^{}e^+`$ and that of $`1.23\times 10^4`$ for the branching ratio of $`\mathrm{\Gamma }(\pi ^\pm e^\pm \nu _e)`$ , we have an estimate for the branching ratio of the transmuational mode $`\pi ^0\mu ^{}e^+`$ of about $`2.9\times 10^9`$. This is barely an order of magnitude below the experimental bound of $`1.72\times 10^8`$ . Indeed, one is saved here from a violation of the experimental limit only by the exceptionally small value of the off-diagonal mass matrix element due to the proximity of the reaction scale $`m_{\pi ^0}`$ to the muon mass where the mass matrix is constrained to be diagonal by the SSD prescription.

(III) As a third example, let us consider at an intermediate scale:

$$\psi l_\alpha \overline{l}_\beta .$$
(5.5)
Following the same procedure as above, one obtains:

$$\frac{\mathrm{\Gamma }(\psi \tau ^{}\mu ^+)}{\mathrm{\Gamma }(\psi \mu ^{}\mu ^+)}\frac{|\mu |m|\tau |^2}{m_\psi ^2}.$$
(5.6)
The matrix element $`\mu |m|\tau `$ for DSM as read from Figure 3 at the scale $`m_\psi `$ is about $`0.03`$ GeV. Combined with the empirical value of about 6 percent  for the branching ratio of $`\psi \mu ^{}\mu ^+`$, one obtains an estimate of about $`6\times 10^6`$ for the transmutational $`\tau ^{}\mu ^+`$ mode, which is larger than in either of the two earlier examples. The corresponding estimate for the NSM is about $`1\times 10^8`$. Although no empirical limit is given in  for this particular mode, the estimates seem close enough to the sensitivity of present experiment to make a search for it worthwhile. Similar arguments applied to the modes $`\psi \tau ^{}e^+`$ and $`\psi \mu ^{}e^+`$ yield the estimated branching ratios $`1.6\times 10^7`$ and $`6\times 10^{11}`$. For these modes again, no empirical limits are given in  but the estimates seem in any case beyond present experimental sensitivity. The same anlysis can also be applied to $`\mathrm{{\rm Y}}`$ decay giving for the branching ratios of the modes $`\tau ^{}\mu ^+,\tau ^{}e^+`$, and $`\mu ^{}e^+`$ respectively the following estimates: $`1.6\times 10^6,6\times 10^8`$, and $`7\times 10^{11}`$. Again, no empirical bounds for these modes are given in , but some, especially the $`\tau \mu `$ mode, will be worth searching for in future in e.g. $`B`$-factories.

(IV) Lastly, consider transmutational fermion decay, of which the prime examples would be:

$$\mu ^{}e^{}\gamma $$
(5.7)
and

$$\mu ^{}e^{}e^+e^{}$$
(5.8)
proceeding via the Feynman diagrams in Figure 4.

The bound as given in  for the branching ratio of the decay (5.7) is an impressive $`4.9\times 10^{11}`$, and of the decay (5.8) is $`1\times 10^{12}`$. And these coming on top of a main decay mode of $`\mu `$, namely $`\mu e\nu \overline{\nu }`$, proceeding by weak interactions with a width of only around $`3\times 10^{10}`$ eV, are a very stringent constraint indeed on any scheme which violates muon number. It is thus gratifying for us to note that transmutation in the DSM scheme survives this very dangerous-looking hurdle automatically. This comes about because the Feynman diagrams 4 for the transmutational decays (5.7) and (5.8) are to be evaluated at the scale of the $`\mu `$ mass, and at that scale the mass matrix in the DSM scheme is diagonal. That it is so by virtue of the SSD prescription for defining lepton states was noted already in preceding sections and can be seen explicitly in Figure 3 by the vanishing of the off-diagonal mass matrix element $`e|m|\mu `$ at the $`\mu `$ mass. Hence the DSM scheme actually predicts zero probability for this transmutational decay.

The same observation applies also to transmutational decays of $`\tau `$, for example:

$$\tau ^{}e^{}\gamma ,\mu ^{}\gamma ,$$
(5.9)
$$\tau ^{}e^{}e^+e^{},e^{}\mu ^+\mu ^{},\mu ^{}e^+e^{},\mu ^{}\mu ^+\mu ^{},$$
(5.10)
since for these decays the analogous Feynman diagrams are to be evaluated at the scale of the $`\tau `$ mass, and at that scale the mass matrix in the DSM scheme is again diagonal giving thus zero transitions. This is also fortunate, for although not quite as stringent as those on (5.7) and (5.8), the bounds on the branching ratios of (5.9) and (5.10) as given in  are still only of order $`10^6`$ for a total $`\tau `$ width of around $`3\times 10^3`$eV, which would be far from trivial to satisfy otherwise.

To illustrate how sensitive a test these lepton decays are to transmutation models, let us examine them in the NSM scenario. Here, as noted before, the rotation of the mass matrix is slower than in the DSM case giving thus usually smaller transmutation rates, especially for transmutations between $`\mu `$ and $`e`$. However, the mass matrix is required by the FSD prescription adopted to be diagonal only at some fixed mass scale. Suppose we choose this scale as the $`\tau `$ mass so as to make the transmutation rates vanish for both the decays (5.9) and (5.10) and so satisfy the experimental bounds on them automatically. As the energy lowers to the scale of the $`\mu `$ mass, however, the mass matrix will have rotated and acquired nonzero off-diagonal elements. Indeed, according to (4.2) the element $`e|m|\mu `$ would have acquired on running from the $`\tau `$ mass a value of around $`1.3\times 10^8`$ GeV, giving thus possibly nonzero rates to the decays (5.7) and (5.8). To estimates these rates following the previous procedure, let us compare for example (5.8) to the known main decay mode of $`\mu `$, namely:

$$\mu ^{}e^{}\overline{\nu }_e\nu _\mu ,$$
(5.11)
as given to first order by the Feynman diagram of Figure 5.

The mass of the electron being so much smaller than that of the muon, we can ignore the small difference in phase space between the two decays and write:

$$\frac{\mathrm{\Gamma }(\mu ^{}e^{}e^+e^{})}{\mathrm{\Gamma }(\mu ^{}e^{}\overline{\nu }_e\nu _\mu )}=\frac{|e|m|\mu |^2}{m_\mu ^2}\frac{e^2/m_\mu ^2}{g_2^2/M_W^2}$$
(5.12)
where we have taken for a typical momentum for the decay to be of order of the muon mass. Inserting the above NSM value for $`e|m|\mu `$ and the empirical values for the other quantities appearing in (5.12), one obtains a value of about $`1\times 10^3`$ for the ratio. This would be disastrous, being some 9 orders too large compared with the experimental bound given in  of order $`1\times 10^{12}`$. That this ratio turns out to be so large despite the very small NSM value for $`e|m|\mu `$ is because of the large $`W`$ mass occuring in the diagram for the normal decay mode of $`\mu `$, or in other words, because we are in (5.12) comparing what is basically, in spite of being off-diagonal, still a photon-mediated process to a weak decay. The above blatant violation of the experimental bound is of course not the fault of the Standard Model but only of the naive fashion that it has been implemented in the NSM scheme since it could be easily avoided by, among other means, choosing a smaller value for the Dirac mass $`m_2`$ of the second heaviest neutrino which is still at present entirely unconstrained by experiment. Nevertheless, this example shows how sensitive these lepton decays are as a test of transmutation models, and how delicately the DSM scheme has managed to survive it.

Obviously, what really mattered for passing that test was what we called the SSD prescription for defining the lepton state vectors which stipulates that the mass matrix is to be diagonal at the mass of every lepton state. Had one chosen in the NSM scheme to define fermion states in this way instead of the FSD prescription adopted, it would have passed the test also. With the FSD prescription, however, it would not have mattered if one had chosen the scale for diagonalizing the mass matrix to be $`m_\mu `$ instead of $`m_\tau `$, for although one would then pass the test for the decay $`\mu ^{}e^{}e^+e^{}`$, one would fail it for such decays as $`\tau ^{}\mu ^{}\mu ^+\mu ^{}`$. Indeed, running now from the $`\mu `$ mass, one would obtain a non-diagonal mass matrix at the $`\tau `$ mass with off-diagonal elements of the same order. The empirical bounds on transmutational $`\tau `$ decays are less stringent, e.g. the present bound on the branching ratio of $`\tau ^{}\mu ^{}\mu ^+\mu ^{}`$ is only $`1.9\times 10^6`$ , but that would still mean a violation of the bound by some 6 orders in magnitude by the above estimate. One could of course impose in the NSM a further condition on the FSD prescription stating that the mass matrix should stop rotating altogether once below the chosen diagonalization scale, but one does not know of a good theoretical reason for doing so.

Firm conclusions cannot be drawn until detailed calculations for the various decay widths have been performed, which to us seem quite feasible in the near future. Tentatively, however, one can claim that the DSM predictions for transmutation survive present experimental bounds in all cases so far studied, while the NSM (with FSD) survives except perhaps for the noted $`\mu `$ or $`\tau `$ transmutational decays. In most cases, the predicted rates are some way below present experimental limits, but for the decays $`\pi ^0e^\pm \mu ^{}`$, $`\psi \mu ^\pm \tau ^{}`$ and $`\mathrm{{\rm Y}}\mu ^\pm \tau ^{}`$ the estimated rates may soon be, if they are not already, within range of experimental sensitivity.

## 6 Remarks

The importance of possible fermion flavour violation has long been recognized and subjected to rigorous test by such experiments with ultra-high sensitivity as $`\mu e\gamma `$ decay. The negative outcome to-date of all these tests conveys the impression that such violations can exist, if at all, only to a very small extent. However, that conclusion was reached in the days before it was realized that the fermion mass matrix actually rotates with the energy scale so that the fact that flavour is conserved at one scale, or that the mass matrix is there diagonal in the flavour states, does not necessarily mean that it will remain so at some other scales. That being the case, the accuracy, though impressive, of such experiments as $`\mu e\gamma `$ searches, is not enough to rule out flavour violation in general but has to be supplemented by experiments at other energies. As our analysis above has shown, under some circumstances the mass matrix can acquire quite appreciable off-diagonal elements without being noticed in existing experiments. Hence, regardless of which scheme for mass matrix rotation one favours, it would be worthwhile for experimenters to perform routine checks under varying circumstances for the diagonality of fermion mass matrices by looking for effects of the type labelled in this paper as transmutations. With a rotating mass matrix, diagonality is no longer as sacrosanct as once believed.

As argued in the Introduction, that the mass matrix rotates forces on us a refinement of our definition of fermion flavour states. Given the sensitivity of such experiments as $`\mu ^{}e^{}e^+e^{}`$ searches, even the small rotations which may occur between the mass scales of two successive generations is enough to give huge violations of the experimental limits, as we have shown in the preceding section. It seems thus necessary either to define the fermion states in such a way as to ensure that the mass matrix is diagonal at the mass scales of all the fermion states, as was done in the SSD prescription, or else in a prescription like FSD to suppress altogether the rotations between fermion states, if not by some theoretical mechanism then by decree. For us, this conclusion has come as an agreeable surprise in that the SSD prescription which was originally invented just for internal consistency of our DSM scheme should turn out now to have the unexpected virtue of guaranteeing consistency with experiment on $`\mu e\gamma `$ and $`\mu ^{}e^{}e^+e^{}`$ decays, which would otherwise be a difficult hurdle to survive. As we shall see in , the same property of the SSD prescription will save one also from certain unwanted pole structures in the amplitudes for the photo-transmutation of leptons.

The rotation of its mass matrix being weaker, the NSM scheme will give generally smaller transmutation effects than the DSM, i.e. provided that it has been amended with something like the SSD prescription for defining fermion states to avoid the noted probable inconsistencies. Indeed, in all the examples so far studied NSM transmutation effects seem to be beyond present experimental sensitivity. That being the case, if in the routine check on transmutation bounds an experiment comes up with an effect, then it is probably due to a rotation of the mass matrix driven by some forces of nature other than those currently considered in the Standard Model. This would be exciting in any case, besides possibly opening up a new window for probing into the generation mystery, as we have already said in the Introduction.

For us advocates of the Dualized Standard Model (DSM), the results of this paper and its companion  are something of a relief, at least temporarily. Since the whole DSM scheme relies heavily on the concept of a rotating mass matrix driven by dual colour forces, on which all its seemingly sensible predictions on the mass and mixing patterns of fermions also depend, it would be a little disastrous if the relatively fast mass matrix rotation required for the scheme’s success predicts at the same time much larger transmutation effects than can be accommodated by present experiment. And, in contrast to flavour-changing neutral current (FCNC) effects  which also pose some hurdles for the scheme to overcome but depend on an unknown mass parameter, there is in transmutations no free parameter which could be adjusted as a loophole for escape. It is thus fortunate that a disaster does not seem so far to have happened although in some cases it has come rather dangerously close. Instead, one obtains as a bonus a whole new class of phenonmena, with quantitative predictions which, while remaining below present bounds, are still not that far off as to be entirely inaccessible to experiment. However, much more work will still be needed first, to make sure that disasters will not happen in circumstances not yet examined, and secondly, to identify specific cases where the scheme’s predictions can be checked by experiment, if not immediately, at least in the not too distant future.

In general terms, what we have done in this paper is to raise some basic questions of concept forced on us by the rotation of the fermion mass matrix and attempt to answer some of them. Although in answering these questions we may have raised as many more, we hope at least to have clarified a little some of their implications.


format: text

text:
# Quasar environments at 0.5≤𝑧≤0.8

## 1. Introduction

The differences and similarities between radio-loud and radio-quiet quasars (hereafter RLQs and RQQs, respectively) have kept astronomers busy for a long time. At essentially all wavelength ranges, except at radio wavelengths, the appearance of radio-loud and radio-quiet quasars is similar. RQQs appear compact with a weak radio component coinciding with the optical quasar nucleus, whereas RLQs have extended lobes of radio emission with hotspots at the outer edges of the radio structure. The radio-emitting lobes are being fed by powerful jets emerging from a bright, central core. RQQs can also have jet-like structures (Blundell & Beasley 1998), although with bulk kinetic powers $`10^3`$ times lower than for RLQs (Miller, Rawlings & Saunders 1993). This suggests that both quasar types have jet-producing central engines, but that the efficiency of the jet production mechanism is very different in the two cases.

One way to learn more about how the two quasar types are related is to study their host galaxies and their galactic environment. A longstanding belief that RQQs are hosted by spiral galaxies and RLQs by ellipticals is now being questioned. Recent studies (Dunlop et al. 1993; Taylor et al. 1996; McLure et al. 1999; Hughes et al. 2000) have found that powerful quasars at $`z>0.5`$, both RLQs and RQQs, seem to exist in galaxies above the break in the luminosity function at $`L^{}`$, but a clear picture has still not emerged. Some studies claim a high fraction of disk morphologies amongst the radio-quiets (e.g. Percival et al. 2000), whilst others suggest that nearly all quasars are in giant ellipticals (e.g. McLure et al. 1999).

The galaxy environment on scales larger than the host galaxy is also interesting. It may provide clues about quasar formation and evolution, since a period of quasar activity may be triggered by interactions and mergers (e.g. Stockton & MacKenty 1983; Ellingson, Green & Yee 1991). Also, a comparison of the galaxy environments around different types of AGN may help constrain the so-called ‘Unified Models’, e.g. the orientation-dependent unified scheme where a RLQ and a radio galaxy are believed to be intrinsically the same type of object, but viewed at different orientations to the line of sight (Barthel 1989). In the orientation-dependent unified scheme for RLQs and radio galaxies, one expects to find that the galaxy environments around the two are the same.

Radio galaxies and RLQs at $`0.5z0.8`$ are often found in environments with richer than average galaxy density (Yee & Green 1987; Hill & Lilly 1991; Ellingson, Yee & Green 1991; Wold et al. 2000), from groups of galaxies and poorer clusters to Abell Class 1 clusters or richer. This is perhaps not surprising since giant elliptical galaxies frequently reside in the centres of galaxy clusters. But what is the galaxy environment like for RQQs in this redshift range, and how does it compare to the environment around RLQs with comparable AGN luminosities?

In order to make a meaningful comparison one needs a method to quantify the galaxy environment that takes into account the depth of the survey, the angular coverage and also corrects for foreground and background galaxies. One such parameter that has been much used when quantifying galaxy environments around AGN is the amplitude of the spatial cross-correlation function, the ‘clustering amplitude’. Yee & López-Cruz (1999) find that the clustering amplitude is a robust estimator of galaxy richness in clusters.

The first attempt at comparing the environments around radio-loud and radio-quiet quasars by using the clustering amplitude was made by Yee & Green (1984). They found only a marginal difference between RLQ and RQQ fields, but later Yee & Green (1987) did deeper imaging of the same RQQ fields, and added more RLQ fields to the study, this time finding an even smaller difference, but they were unable to draw any firm conclusions due to the small number (seven) of RQQ fields.

The problem was later addressed by Ellingson et al. (1991b), who also used the clustering amplitude to quantify the environments around a sample of 32 RLQs and 33 RQQs at $`0.3<z<0.6`$. They found the RLQs to exist in richer than average galaxy environment, and frequently also in clusters as rich as Abell class 1. The RQQs, on the other hand, were found to be much less frequently situated in rich galaxy environments, suggesting that the two quasar types may be physically different objects. Since this study, there has not been much work aimed at comparing the Mpc-scale environments of the two quasar populations at moderate redshifts using the clustering amplitude.

We have therefore undertaken a study using the NOT to collect data on the fields around two samples of $`0.5z0.8`$ RLQs and RQQs that are matched in redshift and AGN luminosity. The redshift range was chosen to extend previous work which went up to $`z0.6`$, to as high a redshift as possible consistent with keeping the redshifted 4000 Å break shortward of the $`I`$-band. We selected quasars randomly from complete flux-limited samples spanning a wide range in both optical and radio luminosity (for the RLQs). By extending the redshift range of previous studies, whilst maintaining the luminosity range, we are able to better disentangle trends in the environmental richness due to cosmic evolution from those due to radio and/or optical luminosity. Our assumed cosmology has $`H_0=50`$ km s-1 Mpc-1, $`\mathrm{\Omega }_0=1`$ and $`\mathrm{\Lambda }=0`$.

## 2. The quasar samples

The RLQ sample was selected from two different radio-optical flux limited surveys, the Molonglo/APM Quasar Survey (Serjeant 1996; Maddox et al., in prep.; Serjeant et al., in prep.) and the 7C quasar survey (Riley et al. 1999) and consists of 21 radio-loud steep-spectrum quasars with redshifts $`0.5z0.82`$ covering a wide radio luminosity range of $`23.8\mathrm{log}\left(L_{408\mathrm{M}\mathrm{H}\mathrm{z}}/\mathrm{WHz}^1\mathrm{sr}^1\right)26.7`$.

The 20 RQQs were selected from three optical surveys with different flux limits in order to cover a wide range in quasar $`B`$-band luminosity within the given redshift range. Eight of the quasars are from the faint Durham/AAT UVX survey of Boyle et al. (1990) and ten of the quasars are from the intermediate luminosity Large Bright Quasar Survey by Hewett, Foltz & Chaffee (1995). There are also two high-luminosity quasars in the sample, selected from the Bright Quasar Survey (BQS) by Schmidt & Green (1983).

In Fig. 1 we show the distribution of the RLQs and the RQQs in the redshift-luminosity plane. There is a slight tendency for the quasars at bright $`M_B`$’s to appear at the highest redshifts. The correlation between redshift and luminosity is a well-known feature in flux-limited samples, but here it may instead be an artifact of the narrow redshift range of the sample. In our analysis we investigate if there are any correlations between the environmental richness and the quasar redshift and luminosity. It is therefore important that there is no underlying correlation between redshift and luminosity in the samples. We use Spearman’s partial rank correlation coefficients for this analysis, allowing the correlation coefficient between two variables (e.g. environmental richness and radio luminosity) to be determined when holding the third variable constant (redshift). For more details about the RLQ sample, see Wold et al. (2000). The RQQ sample will be presented by Wold et al. (in prep.).

## 3. Control fields

Since we are concerned with investigating if there is an excess of galaxies in the quasar fields, we aim to have a good determination of the foreground and background counts. For this purpose, we obtained several images of random fields in the sky at approximately the same galactic latitudes as the quasar fields. There are twelve different control fields of which five were imaged in two filters. In total, they cover 18 arcmin2 in $`V`$, 58.7 arcmin2 in $`R`$ and 73.7 arcmin2 in $`I`$. They were imaged along with sources in the quasar sample, so they have the same depth and were obtained in exactly the same manner as the quasar fields. This is important in order to obtain a robust estimate of the galaxy clustering in the quasar fields (Yee & López-Cruz 1999).

## 4. Observations

Most of the data were obtained using the High Resolution Adaptive Camera (HiRAC) at the NOT, equipped with either the 1k SiTe or the 2k Loral CCD giving a field of view of 3$`\times `$3 and 3.7$`\times `$3.7 arcmin, respectively. Some images were also obtained using the ALFOSC in imaging mode, with the 2k Loral CCD with pixel scale 0.189 arcsec. Typically, the integrations were divided into four exposures of 600 s each. Two RLQ fields were imaged with the 107-in telescope at the McDonald Observatory, and five RLQ fields with the HST (Serjeant, Rawlings & Lacy 1997).

The bulk of the data were obtained under photometric conditions, and the seeing FWHM is less than one arcsec in 11 out of the 14 RLQ fields imaged with the NOT, and in 15 out of the 20 RQQ fields. In several images ‘fuzz’ from the host galaxy is clearly visible around the stellar image of the quasar.

The filters were chosen such as to give preference to early-type galaxies with strong 4000 Å breaks at the quasar redshifts. For an early-type galaxy at $`z0.67`$ the 4000 Å break moves from $`R`$-band into $`I`$-band, so for the $`z0.67`$ quasars we used $`I`$-band imaging and for the $`z<0.67`$ quasars we used $`R`$-band. We also imaged 20 of the quasar fields in two filters, either $`V`$ and $`R`$, or $`R`$ and $`I`$ depending on the redshift of the quasar so as to straddle the rest frame 4000 Å break.

For photometry and object detection we processed the images in focas (Faint Object Classification and Analysis System). By performing completeness simulations in the images we find that the data are complete down to 24.0 in $`V`$, 23.5 in $`R`$ and 23.0 in $`I`$, with errors of $`\pm `$0.3 mag at the limits. All quasar fields, except two, lie at galactic latitudes $`|b|>42^{}`$ and have galactic reddening $`E\left(BV\right)<0.063`$. We corrected for galactic extinction using an electronic version of the maps by Burstein & Heiles (1982) and the Galactic extinction law by Cardelli, Clayton & Mathis (1989).

## 5. Results

To quantify the amount of excess galaxies in the quasars fields we use the amplitude, $`B_{\mathrm{gq}}`$, of the spatial galaxy–quasar cross-correlation function, $`\xi \left(r\right)=B_{\mathrm{gq}}r^\gamma `$, where $`\gamma =1.77`$. The amplitude is evaluated at a fixed radius of 0.5 Mpc at the quasar redshift, corresponding to $`1`$ arcmin at $`z=0.7`$ and has units of Mpc1.77. Longair & Seldner (1979) showed that $`B_{\mathrm{gq}}`$ can be found by first obtaining the amplitude of the angular cross-correlation function, which is directly proportional to the relative excess of galaxies. See Wold et al. (2000) for details of the analysis.

We counted the number of galaxies within the 0.5 Mpc radius in the quasar fields and averaged the counts for the $`R`$ and $`I`$-band data, i.e. for the $`z<0.67`$ and the $`z0.67`$ quasar fields. The average counts are shown in Fig. 2 where we also have plotted the average galaxy counts from the control images for comparison. The two left plots in Fig. 2 show the average $`R`$-band counts for the $`z<0.67`$ RLQ (13) and RQQ (7) fields, and the two plots to the right show the average $`I`$-band counts for the $`z0.67`$ $`I`$-band RLQ (6) and RQQ (13) fields. In the RLQ fields there appears to be a small excess of faint galaxies at $`R`$ and $`I>21`$, whereas the $`R`$-band RQQ fields show no excess. However, there is a clear excess of galaxies in the $`z0.67`$ RQQ fields at $`I>20`$. The errors in the background galaxy counts were calculated as 1.3$`\sqrt{N}`$ in order to take into account the non-random fluctuations in the counts due to the clustered nature of field galaxies.

We calculated the net excess of galaxies in each quasar field by subtracting the average background counts, and thereafter computed the clustering amplitude, $`B_{\mathrm{gq}}`$. In Fig. 3 we show $`B_{\mathrm{gq}}`$ for the RLQ and the RQQ fields as a function of redshift and quasar $`B`$ absolute magnitude. The dotted line across the plots show the value obtained by Davis & Peebles (1983) for the amplitude of the local galaxy-galaxy auto-correlation function, $`B_{\mathrm{gg}}=60`$ Mpc1.77. The mean clustering amplitudes for the RLQ and the RQQ samples are 213$`\pm `$66 and 189$`\pm `$83 Mpc1.77, respectively. We thus make two observations, first that the mean clustering amplitudes for the quasar fields are significantly larger than that of local galaxies, implying that the quasars exist in fields with richer than average galaxy density. Second, we note that the mean clustering amplitudes for the two samples are practically indistinguishable, i.e. on average, there is no difference in the galaxy environments on 0.5 Mpc scales for the RLQs and the RQQs. On average, both the RLQs and the RQQs seem to prefer environments similar to galaxy groups and poorer clusters of galaxies.

The quasar environments span a wide range, however. Some individual fields show no significant excess of galaxies, and other fields appear to be very rich, e.g. the field around the RQQ BQS 1538$`+`$477 with $`B_{\mathrm{gq}}`$ in the range 1100–1200 Mpc1.77 and a galaxy excess of $`40`$–50 galaxies. Wold et al. (2000) argue that an amplitude of $`740`$ Mpc1.77 corresponds to Abell richness class $`>`$$``$ 1, so the cluster candidate around BQS 1538$`+`$447 must qualify as a richness class 2 cluster. Two other fields with clustering amplitudes of 785$`\pm `$255 and 703$`\pm `$250 Mpc1.77 are probable Abell class 1 clusters. Plotting the galaxies in the four richest RQQ fields in a colour-magnitude diagram, reveals a hint of a red sequence at $`RI1.5`$–1.6, i.e. tentative evidence that these fields likely contain galaxy clusters at $`z0.7`$–0.8, see Fig. 4.

As seen in Fig. 3 there are no obvious trends in $`B_{\mathrm{gq}}`$ with either redshift or quasar luminosity. There is a hint that the low-$`z`$ RQQ fields have lower clustering amplitude than the higher redshift RQQ fields, but this is most likely an artifact of the narrow redshift range. The Spearman partial rank correlation coefficient giving the correlation between $`B_{\mathrm{gq}}`$ and redshift is 0.4 with a 1.5$`\sigma `$ significance.

## 6. A link between radio luminosity and environmental richness?

Looking at $`B_{\mathrm{gq}}`$ for the RLQs as a function of radio luminosity, reveals a weak, but significant, correlation between $`B_{\mathrm{gq}}`$ and radio luminosity with much scatter. This is shown in Fig. 5 where we also have plotted $`B_{\mathrm{gq}}`$’s for RLQ fields as found by Yee & Green (1987), Ellingson et al. (1991b) and Yee & Ellingson (1993). Here the, correlation coefficient between $`B_{\mathrm{gq}}`$ and $`L_{408\mathrm{M}\mathrm{H}\mathrm{z}}`$, holding redshift constant, is 0.4 with a 3.4$`\sigma `$ significance.

Does the correlation between radio luminosity and environmental richness imply that environment is the primary factor in controlling the radio luminosity of a RLQ? There are at least three ways in which such a situation could come about. In the first, the environment determines the bulk kinetic power in the radio jets. According to the relation between a galaxy’s black hole mass and the mass of the spheroidal component (Kormendy & Richstone 1995; Magorrian et al. 1998), giant elliptical hosts of radio galaxies and RLQs should have high black hole masses, $`10^8`$$`10^9`$ M. Assuming that the radio jets are powered by accretion and that the accretion rate is proportional to the black hole mass, galaxies with massive black holes will power more luminous radio sources. These massive galaxies will prefer richer environments and thus the correlation between radio luminosity and environmental richness may just reflect an increasing mass of the host.

The second possibility is that more fuel for the quasar and its radio jets is available in a group or cluster environment. A group or a poor cluster environment may be ideal for the fuelling of a black hole as encounters will be more common than in the field and will be of low enough relative velocity to disrupt the interacting galaxies and cause gas to flow into the centre (Ellingson et al. 1991a).

A third possibility is that the radio luminosity is almost independent of the bulk kinetic power in the radio jets and is instead largely determined by the density of the environment into which the source expands. Wold et al. (2000) constructed a simple model with these assumptions and found that the predicted relation between $`B_{\mathrm{gq}}`$ and radio luminosity was much too steep to fit the data, thereby ruling out as strong a $`B_{\mathrm{gq}}`$$`L_{408\mathrm{M}\mathrm{H}\mathrm{z}}`$ dependence as we would see if all RLQs had the same jet power and environment was entirely responsible for determining the radio luminosity. Nevertheless, as the luminosity function for the radio jet power is likely to be steeply declining at high powers, it seems not unlikely that selection effects could operate to produce some correlation between $`L_{408\mathrm{M}\mathrm{H}\mathrm{z}}`$ and $`B_{\mathrm{gq}}`$ without it being as strong as it would be in this rather extreme model in which the jet power is the same for all sources.

Given the large scatter in the correlation it is however quite possible that both environment and radio jet power play important roles in determining the radio luminosity. The relationship between radio sources and their environments must be complex, and the vast majority of radio sources may lie in some sort of cluster-like environment, from groups of only a few galaxies to clusters as rich as Abell class 1 or more.

Given this correlation, one might expect the RQQs, with radio luminosities typically two-three orders of magnitude lower than the RLQs, to be sited in poorer environments. Instead we find that the environments around RLQs and RQQs are indistinguishable. This is however fully consistent with the $`B_{\mathrm{gq}}`$-radio luminosity correlation. That the large-scale environments around RLQs and RQQs are similar suggests that the process that decides on radio loudness in a quasar is not dependent on the environment on Mpc scales, but may instead be found in the central regions of the host galaxy.

When an epoch of quasar activity is triggered in a galaxy, e.g. as a result of galaxy interactions and mergers, then some central process decides whether it becomes a RLQ or a RQQ. If a RLQ is born, the environment into which the source expands will to some extent determine its radio luminosity, and we observe the $`B_{\mathrm{gq}}`$-$`L_{408\mathrm{M}\mathrm{H}\mathrm{z}}`$ correlation. If instead a RQQ is born, we do not observe any correlations between quasar luminosity (in this case optical luminosity) and environmental richness, because this quasar does not have extended radio lobes that can interact with the surrounding galaxy environment. The radio jets and the radio lobes extending beyond the host galaxy thus work like sensors from which we can read off the physical conditions in the intergalactic medium.

## 7. The evolutionary states of AGN-selected clusters

The fuelling mechanism of luminous AGN is still not understood, although it is thought that companion galaxies in groups and clusters may be able to supply fuel either through mergers or via a cluster cooling flow (e.g. Hall, Ellingson & Green 1997). The merger and cooling flow models, however, make very different predictions for the evolutionary state of the cluster surrounding the AGN. If the AGN is fuelled by mergers and interactions, we might expect that the cluster is still forming by merging of sub-clumps, but in the cooling flow scenario the cluster may be well-established and virialized. We have therefore started to investigate the evolutionary states of AGN-selected galaxy clusters by weak lensing techniques.

Using deep images obtained with the ALFOSC on the NOT in sub-arcsec seeing, we have mapped the projected mass distribution in rich, X-ray luminous AGN-selected clusters, and preliminary results show that we have comfortable detections of the weak lensing signal ($`>3\sigma `$). These observations will also allow us to estimate the mass-to-light ratio for the clusters, and to investigate whether the mass-to-light ratio is different than for clusters selected solely on the basis of bright optical or X-ray emission from their baryonic matter component.

The weak lensing technique is being used more and more frequently since it is a powerful method for investigating galaxy clusters (see e.g. Dahle et al., these proceedings). The combination of the good seeing conditions at the NOT and the wide field imagers with good resolution currently available and underway at the NOT (the focalreducer FRED) makes the NOT a powerful instrument for doing weak lensing observations.

### Acknowledgments.

We are grateful to the staff at the NOT and the McDonald Observatory for help with the observations. We also thank the British Research Council and the Research Council of Norway for support. MW also wishes to thank the organizers of the meeting for support.

## References

Barthel, P.D., 1989, ApJ, 336, 606

Blundell, K.M., & Beasley, A.J., 1998, MNRAS, 299, 165

Boyle, B.J., Fong, R., Shanks, T., & Peterson, B.A., 1990, MNRAS, 243, 1

Burstein, D., & Heiles, C., 1982, AJ, 87, 1165

Cardelli, J.A., Clayton, G.C., & Mathis, J.S., 1989, ApJ, 345, 245

Clowe, D., Luppino, G.A., Kaiser, N., Henry, J.P., & Gioia, I.M., 1998, ApJL, 497, 61

Davis, M., & Peebles, P.J.E., 1983, ApJ, 267, 465

Dunlop, J.S., Taylor, G.L., Hughes, D.H., & Robson, E.I., 1993, MNRAS, 264, 455

Ellingson, E., Green, R.F., & Yee, H.K.C., 1991a, ApJ, 378, 476

Ellingson, E., Yee, H.K.C., & Green, R.F., 1991b, ApJ, 371, 49

Hall, P., Ellingson, E., & Green, R.F., 1997, AJ, 113, 1179

Hewett, P.C., Foltz, C.B., & Chaffee, F.H., 1995, AJ, 109, 1498

Hill G.J., & Lilly, S.J., 1991, ApJ, 367, 1

Hughes, D.H., Kukula, M.J., Dunlop, J.S., & Boroson, T., 2000, preprint (astro-ph/0002021)

Kormendy, J., & Richstone, D., 1995, ARA&A, 33, 581

Longair, M.S., & Seldner, M., 1979, MNRAS, 189, 433

Luppino, G.A., & Kaiser, N., 1997, ApJ, 475, 20

Magorrian, J., et al., 1998, ApJ, 115, 2285

McLure, R.J., Kukula, M.J., Dunlop, J.S., Baum, S.A., O’Dea, C.P., & Hughes, D.H., 1999, MNRAS, 308, 377

Miller, P., Rawlings, S., & Saunders, R., 1993, MNRAS, 263, 425

Percival, W.J., Miller, L., McLure, R.J., & Dunlop, J.S., 2000, preprint (astro-ph/0002199)

Riley, J.M., Rawlings, S., McMahon, R.G., Blundell, K.M., Miller, P., Lacy M., & Waldram, E.M., 1999, MNRAS, 307, 293

Schmidt, M., & Green, R.F., 1983, ApJ, 269, 352

Serjeant, S.B.G., 1996, DPhil Thesis, Univ. Oxford

Serjeant, S., Rawlings, S., & Lacy, M., 1997, in Quasar Hosts Conference, D.L. Clements, & I. Perez-Fournon, Springer–Verlag, Berlin, 188

Stockton, A., & MacKenty, J.W., 1983, Nat, 305, 678

Taylor, G.L., Dunlop, J.S., Hughes, D.H., & Robson, E.I., 1996, MNRAS, 283, 930

Wold, M., Lacy, M., Lilje, P.B., & Serjeant, S., 2000, MNRAS, in press (astro-ph/9912070)

Yee, H.K.C., & Ellingson, E., 1993, ApJ, 411, 43

Yee, H.K.C., & López–Cruz, O., 1999, AJ, 117, 1985

Yee, H.K.C., & Green, R.F., 1984, ApJ, 280, 79

Yee, H.K.C., & Green, R.F., 1987, ApJ, 319, 28


text:
# Chaos in a modified Hénon-Heiles system describing geodesics in gravitational waves

## 1 Introduction

Various open Hamiltonian systems in which unbounded orbits exhibit a chaotic behaviour have been studied (in particular in the context of chaotic scattering), see e.g. - and elsewhere. A standard method of investigation and description of these systems is to examine a typical set of initial conditions and observe which of the outcome possibilities is assigned to each initial data. The presence of fractal basin boundaries separating the possible escapes in the space of initial conditions is a suitable and invariant characterization of chaos -.

Among the simplest Hamiltonians of this type are those with two degrees of freedom and polynomial potentials, in particular the famous Hénon-Heiles (HH) system  given by

$$H=\frac{1}{2}(p_x^2+p_y^2)+V(x,y),$$
(1)
where

$$V=V_{HH}(x,y)\frac{1}{2}(x^2+y^2)+\frac{1}{3}(x^33xy^2).$$
(2)
In the literature, however, this is usually studied at energies below the critical escape value $`E_c=\frac{1}{6}`$, for which motions are bounded. These are regular for lowest $`E`$ but become chaotic at higher energies still obeying $`E<E_c`$. Chaos is also present above $`E_c`$, where unbounded orbits choose between three possible distinct escapes to infinity (see, e.g. ). The corresponding fractal basin boundaries have recently been investigated in .

In this paper a system with a potential similar to the standard HH model (2) is presented and studied. This naturally arises in the context of test particle motion in plane-fronted gravitational waves with parallel rays (pp-waves) . It is well-known that in the simplest (so called homogeneous) pp-waves the geodesic motion is regular but, as shown recently -, it is chaotic in the non-homogeneous case. It is the purpose of this paper to discuss in detail a combined effect of both homogeneous and non-homogeneous gravitational-wave contributions. In section 2 the geodesic equations of motion in pp-wave spacetimes are presented in a suitable form in which the relation to Hamiltonian dynamics is established. Specific energy manifolds are then described and classified in section 3. The onset of chaos investigated by numerical simulation of orbits in an interesting medium energy region is described in section 4. Finally, in section 5 the box-counting dimension of the fractal basin boundaries is calculated.

## 2 Geodesics in gravitational pp-waves with both homogeneous and non-homogeneous components

The line element of the well-known class of gravitational pp-wave solutions  in vacuum can be written

$$\mathrm{d}s^2=2\mathrm{d}\zeta \mathrm{d}\overline{\zeta }2\mathrm{d}u\mathrm{d}v(f+\overline{f})\mathrm{d}u^2,$$
(3)
where $`f(u,\zeta )`$ is an arbitrary complex function of the retarded time $`u`$ and the complex spatial coordinate $`\zeta `$. As shown in , , to determine the geodesic motion in these space-times it is crucial to solve the complex equation

$$\ddot{\zeta }+(U^2/2)\frac{\overline{f}}{\overline{\zeta }}=0,$$
(4)
for $`\zeta (\tau )`$, where a dot denotes the derivative with respect to the normalized ‘time’ parameter $`\tau `$, and $`U\dot{u}`$ is a constant. (The remaining functions $`u(\tau )`$ and $`v(\tau )`$ can be found by subsequent integrations.) Introducing two real coordinates by $`xe\zeta `$, $`ym\zeta `$, the Eq. (4) can be considered as a consequence of the Hamiltonian (1) with the potential $`V(x,y,u)=(U^2/2)ef`$.

The space-time (3) with the function $`f`$ linear in $`\zeta `$ represents the trivial flat Minkowski universe. The quadratic function $`f\zeta ^2`$ describes homogeneous pp-waves (see  and references therein), in which the particle motions are integrable. Geodesics in non-homogeneous pp-waves given by $`f\zeta ^n,n=3,4,5,\mathrm{}`$, have only recently been investigated , . The corresponding polynomial potential $`V(x,y)e\zeta ^n`$ is called an ‘$`n`$-saddle’ (a ‘monkey saddle’ for $`n=3`$). It has been found that motion in these gravitational waves with constant ‘amplitude’ are unbounded with $`n`$ possible escapes to $`|\zeta |\mathrm{}`$, and that the basin boundaries separating these outcomes are fractal. For non-constant amplitudes describing sandwich and impulsive gravitational waves a smearing of chaos is observed .

Now, it is natural and also physically interesting to investigate the combined effect of quadratic and higher polynomial terms in the function $`f`$. This is the main purpose of our work here. The simplest situation of this type occurs when $`f`$ has both quadratic and cubic terms,

$$f(\zeta )=\alpha \zeta ^2+\frac{2}{3}\beta \zeta ^3,$$
(5)
where $`\alpha `$ and $`\beta `$ are real positive constants. With this, (4) reduces to $`\ddot{\zeta }+(\alpha \overline{\zeta }+\beta \overline{\zeta }^2)U^2=0`$, in which the (non-trivial) constants $`\alpha `$, $`\beta `$ and $`U`$ can be removed by a simple transformation

$$\zeta (\tau )\frac{\beta }{\alpha }\zeta \left(\frac{\tau }{U\sqrt{\alpha }}\right).$$
(6)
In the real coordinates $`x`$ and $`y`$, the equations of motion then take the form

$$\ddot{x}+x+x^2y^2=0,\ddot{y}y2xy=0.$$
(7)
These equations can be considered to follow from the Hamiltonian (1) with the potential

$$V_m(x,y)=\frac{1}{2}(x^2y^2)+\frac{1}{3}(x^33xy^2),$$
(8)
which may be called a ‘modified Hénon-Heiles potential’ since this has almost the form of the standard Hénon-Heiles potential (2) with the exception of the negative sign of the $`y^2`$ term. However, as will be shown below, this ‘small’ difference results in a behaviour which is qualitatively different from that observed in the well-known HH system. Note that, from a previous investigation of a generalized HH system $`V_g(x,y)=\frac{1}{2}(Ax^2+By^2)+\frac{1}{3}(Cx^33Dxy^2)`$ ‘with adjustable coefficients’, it is known that the system (8) (corresponding to the choice $`A=C=D=1`$, $`B=1`$ in $`V_g`$) is not integrable .

## 3 Classification of possible orbits according to a conserved energy

First, we shall investigate specific properties of the polynomial potential (8) and compare them with those of the standard HH potential (2) and the monkey saddle potential $`V_{ms}(x,y)=\frac{1}{3}(x^33xy^2)`$, in which the quadratic terms are missing. There are two equilibrium saddle points of the potential $`V_m`$ at $`(x,y)=(0,0)`$ and $`(1,0)`$. This is different from the HH potential (2) for which there is an elliptic point at the origin (in addition to three saddle points), and from the monkey saddle potential which has only one unstable equilibrium point (at the origin). The equipotentials $`V_m(x,y)=E=\text{const.}`$ of the potential (8) are drawn in Fig. 1. Let us now briefly comment on various possible situations shown in the figure:

* For $`E<0`$ the energy manifold $`H=E`$ consists of three disjoint components (separated by a ‘forbidden’ grey region, see Fig. 1), each of which contains a single exit to infinity. The similar behaviour is observed for negative energies also in the monkey saddle potential $`V_{ms}`$ and in the HH system.
* For $`E=0`$ the two components of the energy manifold join at $`(0,0)`$, the third one remains separated.
* For $`0<E<\frac{1}{6}`$ the energy manifold consists of two disjoint components, one containing a single exit to infinity ($`x\mathrm{}`$), the second containing two distinct exits ($`y+\mathrm{},y\mathrm{}`$). This situation does not occur in the monkey saddle potential, in which all energy manifolds with $`E>0`$ are connected and contain three exits to infinity. The energy manifolds $`0<E<\frac{1}{6}`$ for the HH system consist of four disjoint components. One corresponds to bounded orbits and the remaining three are unbounded, each containing a single exit to infinity.
* For $`E=\frac{1}{6}`$ the two components of the energy manifold join at the point $`(1,0)`$.
* For $`E>\frac{1}{6}`$ the energy manifold is connected and contains three exits to infinity in this case. The topology of energy manifolds is the same for both the HH and the monkey saddle systems.

The orbits in the negative energy region cannot choose between different exits, so that no basin boundaries are present. On the other hand, the high energy region ($`E>\frac{1}{6}`$) provides three possible exits to infinity and numerical simulations show that the basin boundaries are fractal indicating the presence of chaos. These simulations give results very similar to those for the monkey saddle presented in , or for unbounded motions in the HH model . This is not surprising as the equipotentials of both the potential (8) and the standard HH potential (2) are well approximated by the monkey saddle equipotentials $`V_{ms}=E`$ for sufficiently high energies.

Thus, the most interesting energy region to be considered here is given by the condition $`0<E<\frac{1}{6}`$. Since the motion for $`E<0`$ is not chaotic while it is chaotic in the high energy region, chaos has to appear somewhere between $`E=0`$ and $`E=\frac{1}{6}`$. This is similar to the standard HH model, where the low energy orbits above $`E=0`$ are regular while they become chaotic for higher energies still below $`E=\frac{1}{6}`$. Such behaviour is not present in the ‘pure cubic’ monkey saddle potential $`V_{ms}`$ in which one can clearly distinguish only two energy regions, negative and positive. The existence of the medium energy region is caused by the presence of the quadratic terms in the potentials (2) and (8). In the following section we shall concentrate on the behaviour of orbits in this interval.

## 4 The onset of chaos in the medium energy region

For an investigation of orbits in the medium energy region $`0<E<\frac{1}{6}`$ by a fractal method (which is useful for the characterization of a sensitive dependence on initial conditions in systems with unbounded orbits) it is necessary to choose a suitable typical set of initial data. Here we may concentrate on motion in the region with two possible exits to infinity, $`y\pm \mathrm{}`$ (the right white region in Fig. 1 d, e, f). Let us consider orbits starting at $`\tau =0`$ from rest from the ‘right branch’ of the level curve $`B_r`$ (the branch of $`V_m(x,y)=E`$ with $`x>0`$), see Fig. 2 for typical trajectories of this type. This choice of initial conditions is also useful for localizing periodic orbits, which play an essential role in establishing the chaotic behaviour in both the monkey-saddle system and the standard HH system above the critical energy , . We can conveniently parametrize the starting points on $`B_r`$ using the polar coordinate angle $`\varphi `$ with respect to the origin (such that $`|\varphi |<\frac{\pi }{6}`$). Each orbit is then uniquely characterized by its energy $`E`$ and a starting angle $`\varphi `$.

In fact, we can concentrate on a smaller range of $`E`$ in the medium energy region. It is easy to prove that for energies $`0<E<\frac{1}{12}`$ the orbits starting from $`B_r`$ cannot cross the $`x`$-axis. Indeed, the second equation of (7) can be rewritten as $`\ddot{y}=y(1+2x)`$. Since $`x>\frac{1}{2}`$ in the considered component of the energy manifold for $`E<\frac{1}{12}`$, the acceleration $`\ddot{y}`$ has the same sign as $`y`$. The basin boundary between the two possible exits (‘up’ and ‘down’) will therefore be given by $`\varphi =0`$, as the orbits starting from $`B_r`$ with $`\varphi >0`$ will always escape directly to $`y\mathrm{}`$, while the orbits with $`\varphi <0`$ will always escape to $`y\mathrm{}`$. They never cross the $`x`$-axis, so the basin boundary is not fractal. The unique orbit with $`\varphi =0`$ is an unstable periodic orbit $`\mathrm{\Pi }_0`$ with trajectory along the $`x`$-axis.

Consequently, chaos has to appear in the energy region $`\frac{1}{12}<E<\frac{1}{6}`$. Numerical simulations show that behaviour of orbits with small $`\varphi `$ becomes much more complicated for energies $`E`$ above a certain critical value $`E_0[\frac{1}{12},\frac{1}{6}]`$. The value of $`E_0`$ will be discussed and specified later. Typical orbits for $`E=0.16`$ are drawn in Fig. 2. In addition to the unstable periodic orbit $`\mathrm{\Pi }_0`$ along the $`x`$-axis, there is an unstable periodic orbit $`\mathrm{\Pi }_3`$ with $`\varphi >0`$ (see the bold curves in Fig. 2). Another one $`\mathrm{\Pi }_2`$ with $`\varphi <0`$ (which is not shown) can be obtained by symmetry.

At this point it is convenient to introduce a simple notation. As mentioned above, there are two possible escapes to infinity for orbits in the region being discussed. Let us assign them a symbol $`j`$, where $`j=2`$ for the exit with $`y\mathrm{}`$, and $`j=3`$ for the exit with $`y\mathrm{}`$ ($`j=1`$ is reserved for the exit $`x\mathrm{}`$, which is present for $`E>\frac{1}{6}`$, to keep the notation consistent with  and ). Also, let us denote by $`\tau _s`$ the time $`\tau >0`$ at which the orbit reaches infinity (we put $`\tau _s=\mathrm{}`$ for orbits that never escape to infinity).

We shall now investigate the detailed structure of the orbits for the typical energy $`E=0.16`$. In Fig. 3 the functions $`j(\varphi )`$ and $`\tau _s(\varphi )`$ are plotted. (These are analogues of the scattering function and the time delay, respectively.) The basin boundaries separating possible outcomes in the initial condition state space are clearly fractal, which is confirmed here by zooming in up to the fourth level. This is an invariant evidence that the behaviour of orbits for energy $`E=0.16`$ is chaotic. At each level of the fractal in Fig. 3, there are three peaks of $`\tau _s(\varphi )`$. The left (middle, right) peak corresponds to orbits asymptotic (up to the given level) to the periodic orbit $`\mathrm{\Pi }_2`$ ($`\mathrm{\Pi }_0`$, $`\mathrm{\Pi }_3`$). Observe that a fractal structure occurs in a region in which there are only two possible outcomes.

Further observations analogous to those described in ,  could be made. Instead, we shall concentrate here on the feature which is not present in the monkey saddle system, namely the dependence of $`j(\varphi )`$ on the energy $`E`$. In addition to a single-parameter $`\varphi `$-dependence in Fig. 3, it is necessary to investigate the basins for the ‘up’ ($`j=3`$) and ‘down’ ($`j=2`$) exits in the two-parameter region $`\{(E,\varphi )|\frac{1}{12}<E<\frac{1}{6},|\varphi |<\varphi _{max}=\frac{\pi }{6}\}`$. The result for $`|\varphi |<0.1`$ is presented in Fig. 4, in which black color denotes $`j=3`$ and white corresponds to $`j=2`$. As expected, the only basin boundary for $`E`$ below the critical energy $`E_00.13`$ is the $`x`$-axis ($`\varphi =0`$). Above the critical energy this basin boundary bifurcates into three ‘branches’, the upper (with $`\varphi >0`$) corresponding to the unstable periodic orbit $`\mathrm{\Pi }_3`$, the middle corresponding to $`\mathrm{\Pi }_0`$ and the lower corresponding to $`\mathrm{\Pi }_2`$. Moreover, the basin boundaries above the critical energy $`E_0`$ are fractal. For example, the function $`j(\varphi )`$ shown in the previous Fig. 3 (which is clearly fractal) can be obtained as a section of Fig. 4 for $`E=0.16`$. Fig. 4 is also suitable for a better localization of the critical energy $`E_0`$ (from Fig. 3 it was only clear that $`E_0<0.16`$). More detailed numerical calculations (not presented here) led us to the value $`E_0=0.1364\pm 0.0001`$.

The chaotic part of the medium energy region ($`E_0<E<\frac{1}{6}`$) in Fig. 4 contains two interesting subregions which deserve special attention: energies just below $`E=\frac{1}{6}`$ and energies just above $`E=E_0`$. These subregions will be discussed below.

Fig. 5 shows an enlargement of Fig. 4 for energies close to $`E=\frac{1}{6}`$. It is obvious that there is another bifurcation at $`E_10.165`$ which is somewhat similar to that which occured at $`E=E_0`$. The fractal structure above $`E_1`$ becomes even more complex, roughly speaking with five peaks in $`\tau _s(\varphi )`$ at each fractal level, instead of three peaks shown in Fig. 3.

However, Fig. 5 indicates more: as the energy $`E`$ approaches $`\frac{1}{6}`$, at least one other bifurcation occurs. This naturally leads us to the following conjecture: There is a infinite number of bifurcations in the fractal basin boundary below $`E=\frac{1}{6}`$ at each of which the basin boundary along $`\varphi =0`$ bifurcates into three branches.

According to the conjecture, the density of bifurcation points grows to infinity as $`E\frac{1}{6}`$, where the third exit (to $`x\mathrm{}`$) is opened. The periodic orbit $`\mathrm{\Pi }_0`$ vanishes while a periodic orbit $`\mathrm{\Pi }_1`$ in the new ‘neck’ of the potential is created. Above $`E=\frac{1}{6}`$ there are three exits to infinity and three corresponding unstable periodic orbits $`\mathrm{\Pi }_1`$, $`\mathrm{\Pi }_2`$, $`\mathrm{\Pi }_3`$, similar to the structure for the monkey-saddle system111There are two other periodic orbits $`\mathrm{\Pi }_4`$ and $`\mathrm{\Pi }_5`$ in the monkey-saddle system, which are more difficult to localize. We do not discuss their analogues in our system. (see , ).

It would be interesting to investigate the transition region around the first bifurcation point $`E=E_0`$ in detail. However, it proved to be very demanding on the accuracy of numerical methods. Fig. 6 indicates why: typical orbits in this energy region reach infinity at $`\tau _s1`$, after ‘bouncing’ many times between the potential walls ($`\tau _s400`$ in Fig. 6, which is still moderate: we were able to follow number of orbits with $`\tau _s>10^4`$). The orbit shown in Fig. 6 started at $`\varphi =0.002`$ with $`E=0.13669`$ and performed 110 bounces before escaping to $`j=2`$. The numerical computation is extremely sensitive to the accuracy of the method. Such sensitivity is present close to any fractal basin boundary shown in Fig. 4, but it is usually restricted to a tiny range of initial conditions. However, it is present in the whole subregion bordered roughly by the rectangle $`E_0<E<E_0+0.0005`$, $`|\varphi |<0.003`$. This region near the first bifurcation and similar (but smaller) regions near the higher bifurcations represent orbits with the maximum level of unpredictability in the system under consideration.

## 5 The box-counting dimension

Finally, we present a more complete characterization of the topology of chaotic orbits in the medium energy region. In discussion so far, we have restricted our analysis to specific (yet important) orbits starting from rest from $`B_r`$. A completely general orbit can start at an arbitrary (not forbidden) point with given energy $`E`$ and direction of the initial impulse. For our purposes it is crucial to investigate the character of the basin boundaries, i.e. their regular vs. fractal nature. All unbounded orbits starting at points $`(x,y)`$ with $`y>0`$ (or $`y<0`$) which never cross the manifold section $`y=0`$ has the same certain exit $`j=3`$ (or $`j=2`$, respectively). Therefore, the basin boundaries can appear only among those orbits which cross $`y=0`$. In order to obtain the basin boundary portrait it is thus sufficient to systematically study orbits which start from the (negative) $`x`$-axis, i.e. from points $`P=(x0,y=0)`$, where $`x[x_{min}(E),0]`$, with initial impulse $`p_x,p_y`$. (A complementary discussion of orbits starting at $`x0,y=0`$ would be analogous.) These orbits belong to the energy manifold $`H=\frac{1}{2}(p_x^2+p_y^2)+V_m(x,0)=E`$. The set of initial conditions considered can thus be represented by a disc in the $`(p_x,p_y)`$ plane in which each point characterize a specific initial velocity $`\dot{x}=p_x,\dot{y}=p_y`$ of an orbit starting at $`(x,0)`$, where $`x0`$ is the solution of the equation $`x^2+\frac{2}{3}x^3=2Ep_x^2p_y^2`$. The circular boundary of the disc corresponds to geodesics starting at $`x=0`$ with maximum permitted impulse $`p_x^2+p_y^2=2E`$, whereas the center of the disc denotes $`p_x=0=p_y`$ and $`x_{min}(E)(1,0)`$. In Fig. 7 we present typical discs for $`E=0.16`$ (left) and $`E=0.166`$ (right) in which each point representing specific initial data has again been colored in black or white according to the exit (‘up’ and ‘down’, respectively) reached by the corresponding evolved orbit. Analogous discs in the monkey saddle potential can be found in , .

It is obvious that the discs contain fractal basin boundaries, patterns and bifurcations which resemble those shown in Fig. 4 and Fig. 5. Interestingly, the geometrical pictures are very similar despite the fact that they describe different situations. In Figs. 4 and 5 the main independent parameter is $`E`$, whereas $`E`$ is fixed in Fig. 7 with the main parameter being $`x`$.

The fractal structure of the basin boundaries can rigorously be characterized by the box-counting dimension , . A large number of randomly chosen initial conditions represented by points in the above discs is evolved. For each point $`P`$ its neighbourhood of radius $`\epsilon `$ is also inspected. If not all the points in the neighbourhood belong to the same basin as $`P`$ then the point is denoted as ‘uncertain’. The fraction $`P(\epsilon )`$ of uncertain initial conditions among all chosen points (we simulated 10000 points for each $`\epsilon `$) scales as $`P(\epsilon )\epsilon ^{2d}`$, where the parameter $`d`$ gives the box-counting dimension of the fractal basin boundaries of the discs shown in Fig. 7. Numerical simulation indicates that $`d=1.34\pm 0.07`$ for $`E=0.16`$, and $`d=1.43\pm 0.06`$ for $`E=0.166`$. In fact, the dimension grows with $`E`$ from $`d=1`$ for $`E<E_00.1364`$ to $`d1.50\pm 0.05`$ as $`E\frac{1}{6}`$.

## 6 Concluding remarks

We have discussed in detail the Hamiltonian system with the potential (8) which is a simple but interesting modification of the famous Hénon-Heiles potential (2). Interestingly, orbits in the system studied correspond to the geodesic motion of free test particles in the space-times describing plane-fronted gravitational waves with both homogeneous and non-homogeneous components in the structural function (5). We have demonstrated that for energies $`E>E_c0.1364`$ the motion becomes chaotic since the basin boundaries separating two possible exits to infinity have a fractal structure. For $`E>\frac{1}{6}`$ the third exit opens and the resulting chaotic behaviour has a similar structure as in the monkey saddle potential $`V_{ms}`$ described in , . This ‘pure cubic’ potential for the monkey saddle system is obtained by setting $`\alpha =0`$ in (5). On the other hand, if the cubic term in the function $`f`$ vanishes ($`\beta =0`$), the gravitational wave is homogeneous and the equations of geodesic motion become linear, i.e. integrable.

## Acknowledgments

We acknowledge the support of grants GACR-202/99/0261 from the Czech Republic and GAUK 141/2000 of the Charles University. We thank the developers of the software system FAMULUS which we used for computation and drawing of all the pictures. Also, we thank Jerry Griffiths for his help with the manuscript.

## Figure Captions

Fig. 1 The equipotentials $`V_m(x,y)=\frac{1}{2}(x^2y^2)+\frac{1}{3}(x^33xy^2)=E`$ for nine values of energy $`E`$ representing all typical situations. The forbidden regions are grey.

Fig. 2 An example of orbits starting from rest from the ‘right branch’ of the level curve $`B_r`$ with energy $`E=0.16`$. The initial positions on $`B_r`$ are parametrized by an angle $`\varphi `$. The orbits $`\mathrm{\Pi }_0`$ and $`\mathrm{\Pi }_3`$ are unstable periodic orbits.

Fig. 3 Plots of the function $`j(\varphi )`$ labelling the two possible exits to infinity and the function $`\tau _s(\varphi )`$ denoting the time at which the orbit starting at given $`\varphi `$ reaches infinity ($`E=0.16`$ here). The enlargements up to the fourth level confirm the fractal structure of basin boundaries.

Fig. 4 The basins of the two exits, ‘up’ (black) and ‘down’ (white), in the subset of initial conditions corresponding to orbits starting from $`B_r`$. These are uniquely characterized by energy $`E`$ and starting angle $`\varphi `$.

Fig. 5 Enlargement of Fig. 4 for energies close to $`E=\frac{1}{6}`$. The basin boundaries are fractal.

Fig. 6 Typical geodesics ($`\varphi =0.002`$, $`E=0.13669`$) ‘bounce’ many times between the walls before escaping to infinity.

Fig. 7 The discs $`(p_x,p_y)`$ for $`E=0.16`$ and $`E=0.166`$ show the topology of orbits in the corresponding energy manifold. Each point in the disc represents an orbit starting at specific point $`x<0`$, $`y=0`$ with a given velocity. Its colour denotes one of the two possible exits. Again, the fractal basin boundaries indicate chaos.