text:
# Anharmonic double-phonon excitations in the interacting boson model

## I Introduction

Vibrational degrees of freedom in atomic nuclei can be described in terms of phonon excitations that arise from nuclear shape oscillations. Vibrations of nuclei with ellipsoidal symmetry can be of two types : $`\beta `$ vibrations which preserve axial symmetry and give rise to a band with $`K=0`$, and $`\gamma `$ vibrations which break axial symmetry and yield a $`K=2`$ band, where $`K`$ is the projection of the angular momentum on the axis of symmetry. At the experimental level, $`\gamma `$ bands have been identified in many well-deformed nuclei; in contrast, the identification of $`\beta `$ bands is still full of questions and difficulties. This is mainly because, when the energy surface has a well-deformed minimum in $`\beta `$ but is rather flat in $`\gamma `$, the $`\beta `$-band increases in excitation energy and approaches the energy region where other degrees of freedom are important. In that case band mixing may occur and can give rise to non-pure structures with decay patterns difficult to identify as those of a $`\beta `$ band .

Since single-$`\gamma `$ excitations are very well established, it is natural to look for double-$`\gamma `$ vibrations and to develop models that can deal with such multi-phonon excitations. Double-$`\gamma `$ excitations correspond to $`K^\pi =0^+`$ and $`K^\pi =4^+`$ bands which are the anti-parallel and parallel combinations of single-$`\gamma `$ phonons, respectively. The experimental identification of two-$`\gamma `$ states in deformed nuclei is difficult because their expected excitation energy is around the pairing gap and hence they can mix strongly with two-quasiparticle excitations. However, recent experimental improvements in nuclear spectroscopy following Coulomb excitation , inelastic neutron scattering , and thermal-neutron capture  have made possible the study of highly excited low-spin states. Many states have been proposed as possible candidates of double-$`\gamma `$ vibrations. There is, however, some controversy about their interpretation. The author of reference  claims that some of the presumed double-$`\gamma `$ states can be interpreted as single hexadecapole-phonon excitations. In fact, to identify the band-head of a double-$`\gamma `$ band it is not sufficient to analyze just $`B(E2)`$ values; data from single-nucleon transfer reactions, $`\beta `$-decay studies, and inelastic scattering experiments must be considered as well. One of the key properties to disregard a band as a double-$`\gamma `$ band is the fact that its members, in first order, cannot be populated in single-nucleon transfer reactions. Many examples of $`K^\pi =4^+`$ states that are identified as double-$`\gamma `$ excitations but are strongly populated in single-nucleon transfer reactions, can be found in the literature: 158Gd, 162Dy, 172Yb, 176,178Hf, and 190,192Os. Since the double-phonon character of the states in question is in doubt, they are not considered here. However, some candidates seem to have a genuine double-phonon nature. Such is the case with 164Dy and 166-68Er. In particular, in reference  a $`K^\pi =4^+`$ state in 164Dy at $`2.173`$ MeV is found to exhibit all properties of a double-$`\gamma `$ band. In references  the observation is reported of $`K^\pi =0^+`$ and $`K^\pi =4^+`$ double-$`\gamma `$ states in 166Er, at energies of $`1.949`$ MeV and $`2.029`$ MeV, respectively. Finally, in 168Er a $`K^\pi =4^+`$ double-$`\gamma `$ excitation is identified at an energy of $`2.055`$ MeV .

One of the most striking features of the observed double-$`\gamma `$ bands is their high anharmonicity, i.e. the ratio of double-$`\gamma `$ over single-$`\gamma `$ energy is different from $`2`$ and ranges from $`2.5`$ to $`2.8`$. This information is very important since it provides a stringent test of nuclear models. The nuclei 164Dy and 166,168Er have been interpreted in the context of many different models such as the quasi-phonon model , the geometrical model , the multi-phonon model , the self-consistent collective-coordinate method , and the $`sdg`$-IBM , and it is now of interest to revisit these models in connection with anharmonic vibrational behavior.

In this paper the simplest version of the interacting boson model (IBM)  is extended by adding to the usual Hamiltonian higher-order interactions between the bosons with the purpose of creating a framework that accommodates the high anharmonicities observed in 164Dy and 166,168Er. The structure of the paper is as follows. First, the IBM is reviewed with special reference to its harmonic character. In section III the inclusion of three-body terms in the Hamiltonian is discussed. The introduction of four-body terms is presented in section IV and some analytic results are pointed out. In section V a detailed study of possible four-body terms is carried out and realistic calculations for 164Dy and 166,168Er are presented. Finally, in section VI the conclusions of this work are made.

## II The IBM-1 as a harmonic model

The IBM describes low-lying collective excitations in even–even nuclei in terms of monopole ($`s`$) and quadrupole ($`d`$) bosons . The boson number that corresponds to a given nucleus equals half the number of valence nucleons ($`N=n/2`$). The rotationally invariant and number-conserving boson Hamiltonian usually includes up to two-body interactions between the bosons although higher-order terms can be added in principle. The most general two-body IBM Hamiltonian can be written in a multipole expansion as

$$\widehat{H}=\epsilon _s\widehat{n}_s+\epsilon _d\widehat{n}_d+\kappa _0\widehat{P}^{}\widehat{P}+\kappa _1\widehat{L}\widehat{L}+\kappa _2\widehat{Q}\widehat{Q}+\kappa _3\widehat{T}_3\widehat{T}_3+\kappa _4\widehat{T}_4\widehat{T}_4,$$
(1)
where $`\widehat{n}_s`$ and $`\widehat{n}_d`$ are the $`s`$\- and $`d`$-boson number operators, respectively, and

$`\widehat{P}^{}`$ $`=`$ $`{\displaystyle \frac{1}{2}}d^{}d^{}{\displaystyle \frac{1}{2}}s^{}s^{},`$ (2)
$`\widehat{L}`$ $`=`$ $`\sqrt{10}(d^{}\times \stackrel{~}{d})^{(1)},`$ (3)
$`\widehat{Q}`$ $`=`$ $`s^{}\stackrel{~}{d}+d^{}\stackrel{~}{s}+\chi (d^{}\times \stackrel{~}{d})^{(2)},`$ (4)
$`\widehat{T}_3`$ $`=`$ $`(d^{}\times \stackrel{~}{d})^{(3)},`$ (5)
$`\widehat{T}_4`$ $`=`$ $`(d^{}\times \stackrel{~}{d})^{(3)}.`$ (6)
The symbol $``$ represents the scalar product; in this paper the scalar product of two operators with angular momentum $`L`$ is defined as $`\widehat{T}_L\widehat{T}_L=_M(1)^M\widehat{T}_{LM}\widehat{T}_{LM}`$ where $`\widehat{T}_{LM}`$ corresponds to the $`M`$ component of the operator $`\widehat{T}_L`$. In the previous equations the operator $`\stackrel{~}{\gamma }_\mathrm{}m=(1)^m\gamma _\mathrm{}m`$ (where $`\gamma `$ refers to $`s`$ or $`d`$) is introduced so that the annihilation operator verifies the appropriate properties under spatial rotations.

It is not a priori clear to what extend $`\beta `$ and $`\gamma `$ vibrations are anharmonic in the IBM even if one just considers the Hamiltonian (1) with up to two-body interactions. A partial analysis of this problem was given in reference . There, the authors find that the IBM in its simplest version is a harmonic model in the limit of infinite boson number $`N`$ and even for finite $`N`$ the model cannot accommodate large anharmonicity if one considers up to two-body interactions; only the interplay between one+two-body terms and higher-order interactions can induce, in principle, a sizeable anharmonicity in the double-phonon excitations. The reason why one-body terms and two-body interactions cannot create a large anharmonicity can be understood as follows. If one considers a Hamiltonian with one parameter that controls the ratio of the strength of the one-body energies and the two-body interactions,

$$\widehat{H}=(1\xi )\left(\epsilon _s\widehat{n}_s+\epsilon _d\widehat{n}_d\right)+\xi \left(\kappa _0\widehat{P}^{}\widehat{P}+\kappa _1\widehat{L}\widehat{L}+\kappa _2\widehat{Q}\widehat{Q}+\kappa _3\widehat{T}_3\widehat{T}_3+\kappa _4\widehat{T}_4\widehat{T}_4\right),$$
(7)
where $`\xi `$ ranges from $`0`$ to $`1`$, one finds two ‘phases’ separated by a critical value, $`\xi _\mathrm{c}`$: a first phase where the one-body term plays the main role ($`\xi <\xi _\mathrm{c}`$) and a second phase where the two-body interaction is the driving force ($`\xi >\xi _\mathrm{c}`$). The crucial point is that the separation between the two phases is very sharp  and essentially no interplay between one- and two-body terms can be found. Since, to a good approximation, the force is either one body or two body but not both, harmonic behavior cannot be avoided.

The inclusion of high-order interactions in a system with a high boson number also leads to a harmonic description. Only for finite boson number the interplay between one+two-body terms and higher-order interactions can induce an anharmonicity in the double-phonon excitations that is comparable to the observed one. These ideas will be used as a guideline in the following sections.

To carry out a quantitative study of anharmonicities, it is convenient to define a ratio between single- and double-phonon excitation energies. Because the experimental situation for $`\beta `$ excitations is not clear we concentrate on $`\gamma `$ vibrations and define energy ratios for $`\gamma `$ phonons only (although similar definitions can be given for $`\beta `$ phonons):

$$R_0^\gamma \frac{E_\mathrm{x}(0_{\gamma \gamma }^+)}{E_\mathrm{x}(2_\gamma ^+)E_\mathrm{x}(2_1^+)},R_4^\gamma \frac{E_\mathrm{x}(4_{\gamma \gamma }^+)E_\mathrm{x}(4_1^+)}{E_\mathrm{x}(2_\gamma ^+)E_\mathrm{x}(2_1^+)},$$
(8)
where $`0_{\gamma \gamma }^+`$ and $`4_{\gamma \gamma }^+`$ are the band heads of the $`K^\pi =0^+`$ and $`K^\pi =4^+`$ double-$`\gamma `$ bands, respectively, and $`E_\mathrm{x}`$ stands for excitation energy. This particular definition removes any rotational influence.

## III Three-body Hamiltonians

Let us consider in the following a Hamiltonian that includes a quadrupole term, a rotational $`\widehat{L}^2`$ term, and three-body interactions between the $`d`$ bosons,

$$\widehat{H}=\kappa \widehat{Q}\widehat{Q}+\kappa ^{}\widehat{L}\widehat{L}+\underset{kl}{}\theta _l\left((d^{}\times d^{})^{(k)}\times d^{}\right)^{(l)}\left((\stackrel{~}{d}\times \stackrel{~}{d})^{(k)}\times \stackrel{~}{d}\right)^{(l)},$$
(9)
where $`\kappa =\kappa _2`$, $`\kappa ^{}=\kappa _1`$. Five independent three-body $`d`$-boson interactions exist which have $`l=0`$, 2, 3, 4, and 6. Interactions with the same $`l`$ but different $`k`$ are not independent but differ by a normalization factor only . The combinations $`(k,l)=(2,0)`$, (0,2), (2,3), (2,4), and (4,6) are chosen here.

The Hamiltonian (9) is certainly not the most general one+two+three-body Hamiltonian that can be considered. Notably, a vibrational term $`ϵ_d\widehat{n}_d`$ which dominates in spherical nuclei is omitted since it is thought of lesser importance in the deformed nuclei considered here. And, of all possible three-body interactions (seventeen terms), only those between the $`d`$ bosons are retained here since they are the more efficient terms to produce anharmonicity .

For the discussion of the anharmonicities of $`\gamma `$ vibrations we study the behavior of the energy ratios (8) as a function of the ratio $`\theta _l/\kappa `$. The identification of the states $`0_{\gamma \gamma }^+`$ and $`4_{\gamma \gamma }^+`$ is based on the $`B(E2)`$ values for decay into the single-$`\gamma `$ states. In figure 1Note that this figure differs from figure 2 in reference  in some scale factors due to an error in the definition of $`\theta _l`$. the influence of the various three-body interactions is shown for a typical value of $`\chi `$ ($`\chi =0.5`$) and for $`N=15`$ bosons. It is seen that a $`\gamma `$-vibrational anharmonic behavior is obtained which can be different for the $`K^\pi =0^+`$ and $`K^\pi =4^+`$ bands (e.g. positive for the former while negative for the latter.) Care has been taken to plot results only up to values of $`\theta _l`$ that do not drastically alter the character of rotational spectrum; beyond these values, the three-body interaction, being of highest order in the Hamiltonian (9), becomes dominant. Also shown in figure 1 are the ratios $`R_K^\gamma `$ as observed in 166Er , $`R_0^\gamma =2.76`$ and $`R_4^\gamma =2.50`$. Figure 2 shows the experimental spectrum of 166Er  and compares it to the eigenspectrum of Hamiltonian (9) with an $`l=4`$ three-body interaction. The parameters are $`\kappa =23.8`$ keV, $`\chi =0.55`$, $`\kappa ^{}=1.9`$ keV, and $`\theta _4=31.3`$ keV Note that the value $`\theta _4`$ is $`1/3`$ of that given in reference  which has an error in its definition. The results shown in that paper are correct after a simple rescaling of the parameter., with boson number $`N=15`$. With these values the calculated excitation energies of the double-$`\gamma `$ band heads are 1926 keV and 1972 keV for the $`K^\pi =0^+`$ and $`K^\pi =4^+`$ levels, respectively, leading to the ratios $`R_0^\gamma =2.82`$ and $`R_4^\gamma =2.45`$, in excellent agreement with observation. Note, however, that although all $`\gamma `$-band heads are well reproduced by the calculation, problems arise for the moments of inertia, in particular of the $`\gamma `$ band. In next sections we will come back on the moments of inertia to see that three-body Hamiltonians provide a very poor description of them.

## IV $`SU(3)`$ Hamiltonians with up to four-body interactions

In the previous section and in reference  a very good description of double-phonon excitation energy has been obtained, but at expense of spoiling the moments of inertia of ground and $`\gamma `$ bands. These drawbacks seem to be a general feature of three-body Hamiltonians. In this and the following section it is shown that the drawbacks of three-body interactions can be overcome by going to the next order.

Since $`\gamma `$ anharmonicity has been observed exclusively in well-deformed nuclei, it is appropriate to consider the problem in the SU(3) limit of the IBM which is suited to deal with nuclei in this mass region . Therefore, in a first approach, the $`SU(3)`$ limit is used and later, in the next section, these results are used as a guidance in more realistic calculations.

Let us consider the following Hamiltonian:

$`\widehat{H}`$ $`=`$ $`a\widehat{C}_2[SU(3)]+b_1\widehat{C}_3[SU(3)]+b_2\widehat{N}\widehat{C}_2[SU(3)]`$ (10)
$`+`$ $`c_1\widehat{C}_2[SU(3)]^2+c_2\widehat{N}\widehat{C}_3[SU(3)]+c_3\widehat{N}^2\widehat{C}_2[SU(3)],`$ (11)
where $`\widehat{C}_n[SU(3)]`$ stands for the Casimir operator of order $`n`$ of $`SU(3)`$. Note the inclusion of cubic terms for completeness in the $`SU(3)`$ analysis. Since only the spectrum of a single nucleus is of interest here, the number of bosons can be fixed in every case and the Hamiltonian (10) can be simplified by combining terms into a single one, leaving a Hamiltonian with three coefficients $`a`$, $`b`$, and $`c`$,

$`\widehat{H}`$ $`=`$ $`a\widehat{C}_2[SU(3)]+b\widehat{C}_3[SU(3)]+c\widehat{C}_2[SU(3)]^2,`$ (12)
with eigenvalues

$`(\lambda ,\mu )|\widehat{H}|(\lambda ,\mu )`$ $`=`$ $`a(\lambda ^2+\mu ^2+\lambda \mu +3\lambda +3\mu )`$ (13)
$`+`$ $`b(\lambda \mu )(2\lambda +\mu +3)(\lambda +2\mu +3)`$ (14)
$`+`$ $`c(\lambda ^2+\mu ^2+\lambda \mu +3\lambda +3\mu )^2.`$ (15)
No quartic Casimir operator exists for $`SU(3)`$ because the number of independent Casimir operators equals the number of labels that characterize an irreducible representation. The Hamiltonian (12) has no rotational term $`\widehat{L}^2`$ since of primary interest, at this point, is the description of band-head energies of single- and double-$`\gamma `$ excitations.

The definition of the energy ratios (8) must now be adapted to incorporate the symmetry labeling of the states. In the $`SU(3)`$ limit the $`\gamma `$ band belongs to the $`SU(3)`$ representation $`(2N4,2)`$ where $`N`$ is the number of bosons. The double-$`\gamma `$ band with $`K=4`$ is contained in the $`(2N8,4)`$ representation; the double-$`\gamma `$ band with $`K=0`$ is predominantly contained in the $`(2N6,0)`$ representation, although an important component is in $`(2N8,4)`$ . In this section the energy ratios (8) are thus defined as follows:

$$R_0^{SU(3)\gamma }\frac{E_\mathrm{x}(2N6,0)}{E_\mathrm{x}(2N4,2)},R_4^{SU(3)\gamma }\frac{E_\mathrm{x}(2N8,4)}{E_\mathrm{x}(2N4,2)}.$$
(16)
Because no rotational term is included, the ratios (16) can be compared directly with equations (8). In the following the effect of the different terms in equation (12) on the degree of anharmonicity of the two-phonon excitation energy is analyzed.

* $`a<0`$, $`b=0`$, $`c=0`$.

This corresponds to the simplest version of IBM; the energy ratios become

$`R_0^{SU(3)\gamma }={\displaystyle \frac{24N18}{12N6}},R_4^{SU(3)\gamma }={\displaystyle \frac{24N36}{12N6}}.`$ (17)
An almost pure harmonic $`\gamma `$-vibrational spectrum is found since the energy ratios (17) are only slightly lower than $`2`$. This will be referred to as negative anharmonicity as opposed to the positive anharmonicity for energy ratios above $`2`$.
* $`a<0`$, $`b0`$, $`c=0`$.

In this case the Hamiltonian (12) is a combination of two- and three-body terms. For given values of $`a`$ and $`b`$ and for a high enough boson number, only the three-body part of the Hamiltonian is dominant and a harmonic spectrum is recovered. For obtaining an anharmonic spectrum the values of the Hamiltonian parameters are very constrained once the number of bosons has been fixed, as the following analysis shows.

Without lost of generality $`a`$ can be fixed to $`a=1`$. The energy ratios are then

$`R_0^{SU(3)\gamma }`$ $`=`$ $`{\displaystyle \frac{34N+b(24N^236N+27)}{(2N1)(6bN+9b1)}},`$ (18)
$`R_4^{SU(3)\gamma }`$ $`=`$ $`2{\displaystyle \frac{2N3}{2N1}}.`$ (19)
To keep the energy of the $`\gamma `$ excitation positive, the value of $`b`$ has an upper limit

$$b<\frac{1}{6N+9}.$$
(20)
The value $`b=1/(6N+9)`$ leads to a divergence in $`R_0^{SU(3)\gamma }`$ and around this value anharmonic behavior is found. From equation (19) one observes that the behavior of the $`(2N8,4)`$ representation is completely harmonic and does not depend on $`b`$. On the other hand, from equation (18) one sees that a wide range of anharmonic ratios is found for the $`(2N6,0)`$ representation. As an illustration, in figure 3 equation (18) is represented as a function of $`b`$, for three values of $`N`$ ($`5`$, $`10`$, and $`15`$). Only positive values of $`b`$ are plotted because for the negative ones $`R_0^{SU(3)\gamma }`$ decreases smoothly to the asymptotic values $`1.27`$, $`1.58`$, and $`1.70`$ for $`N=5`$, $`10`$, and $`15`$, respectively.

The conclusion is that a Hamiltonian with $`c=0`$ does not agree with the experimental situation observed in the mass region of well-deformed nuclei, where the $`\gamma _{K=4}^2`$ state is highly anharmonic.
* $`a<0`$, $`b=0`$, $`c0`$.

In this case the Hamiltonian (12) is a combination of two- and four-body terms. As in the previous case, anharmonicity requires very constrained Hamiltonian parameters once the number of bosons is fixed.

Again, without loss of generality, we fix $`a=1`$. The energy ratios then read as follows:

$`R_0^{SU(3)\gamma }`$ $`=`$ $`{\displaystyle \frac{(4N3)(2c(4N^26N+9)1)}{(2N1)(8cN^2+6c1)}},`$ (21)
$`R_4^{SU(3)\gamma }`$ $`=`$ $`{\displaystyle \frac{2(2N3)(4c(2N^23N+9)1)}{(2N1)(8cN^2+6c1)}}.`$ (22)
To keep the energy of the $`\gamma `$ excitation positive, the value of $`c`$ has the upper limit

$$c<\frac{1}{8N^2+6}.$$
(23)
The value $`c=1/(8N^2+6)`$ produces a divergence in the two energy ratios and in its neighborhood highly anharmonic behavior is found. In figures 4 and 5 are plotted the ratios $`R_0^{SU(3)\gamma }`$ and $`R_4^{SU(3)\gamma }`$, respectively. Again, for negatives values of $`c`$, $`R_0^{SU(3)\gamma }`$ goes asymptotically to the values $`1.45`$, $`1.69`$, and $`1.78`$ for $`N=5`$, $`10`$, and $`15`$, respectively, while $`R_4^{SU(3)\gamma }`$ goes to $`1.33`$, $`1.59`$, and $`1.71`$. These negative anharmonicities have no phenomenological interest.

In this case the experimental situation can be nicely described. Both $`\gamma _{K=0}^2`$ and $`\gamma _{K=4}^2`$ states can be accommodated in a anharmonic description. The conclusion is thus that a Hamiltonian with $`\widehat{C}[SU(3)]`$ and $`\widehat{C}[SU(3)]^2`$ terms seems to be a good starting point to treat the $`\gamma `$ anharmonicity in deformed nuclei.

The description presented here only provides band heads and, to recover a rotational structure, a $`\widehat{L}^2`$ term must be included. The rotational structure is the same in every band because in the $`SU(3)`$ limit no mixing exists between rotational and vibrational degrees of freedom. In next subsection these results are illustrated with a schematic calculation for energies and transition probabilities.

### A A schematic application

Let us consider the case of 166Er which is, as already mentioned, one of particular interest because both double-$`\gamma `$ excitations (with $`K=0`$ and $`K=4`$) have been identified . To carry out the schematic calculation, we use the Hamiltonian (12) with $`b=0`$. The experimental values for the single and double-$`\gamma `$ energy ratios are $`R_0^\gamma =2.76`$ and $`R_4^\gamma =2.50`$ and can be compared directly with the expressions (21-22) leading two values for $`c/a`$, namely $`5\times 10^4`$ and $`1\times 10^4`$. Both solutions are fairly close and any value in between them will correctly describe the anharmonicity of the $`K=0`$ and $`K=4`$ bands. The value of $`a`$ and the strength of the rotational term are fixed from the excitation energies of the $`2_1^+`$ and $`2_2^+`$ levels. After these simple considerations one arrives at the following Hamiltonian:

$$\widehat{H}=13.43\widehat{L}\widehat{L}20.84\widehat{C}_2[SU(3)]+9.29610^3\widehat{C}_2[SU(3)]^2,$$
(24)
where the coefficients are given in keV. The theoretical and experimental spectra are compared in figure 6 and a very good agreement is obtained. However, due to the simplicity of the calculation, $`\gamma `$ and $`\beta `$ bands are degenerate in energy which is not the case experimentally.

To complete the description, $`E2`$ transition probabilities must be computed also. The calculation of $`B(E2)`$ values in the $`SU(3)`$ limit, as in other situations where degeneracies occur, must be treated with care and an appropriate basis must be chosen for states with the same energy. A natural way to do this work is to slightly lift the degeneracy of the $`SU(3)`$ Hamiltonian. The levels can be split using in the Hamiltonian a value $`\chi =1.30`$ which is very close to its $`SU(3)`$ value. With this $`SU(3)`$ breaking the degeneracy is lifted in a natural way because the $`\beta `$ band is pushed up in energy, as is observed. One may expect that with this small change in $`\chi `$ the $`SU(3)`$ spectrum will keep its properties. The $`E2`$ transition operator is:

$$\widehat{T}(E2)=e_{\text{eff}}(s^{}\stackrel{~}{d}+d^{}\stackrel{~}{s}+\chi (d^{}\times \stackrel{~}{d})^{(2)}).$$
(25)
The value of $`\chi `$ that best reproduces the data is $`\chi =0.26`$. Note that $`\chi `$ in the $`\widehat{T}(E2)`$ operator and in the Hamiltonian is different. The effective charge is fixed to reproduce $`B(E2;2_1^+0_1^+)`$: $`e_{\text{eff}}^2=(1.97)^2`$ W.u.. In table I theoretical and experimental transition rates involving the ground and the $`\gamma `$ bands are compared.

This simple analysis suggests that a four-body operator of the type in (24) provides a good description of both single- and double-$`\gamma `$ bands. In the next section this schematic analysis is extended to non-$`SU(3)`$ situations.

## V General Hamiltonians with up to four-body interactions

A general Hamiltonian with all possible three- and four-body terms can, in principle, be constructed but the number of parameters is so high that a study, even a schematic one, of the effect of the different terms on energy spectra and electromagnetic transitions is impossible. Schematic IBM Hamiltonians have been used for many years and, in particular, the quadrupole-quadrupole interaction has been very successful in describing a wide variety of nuclear spectra . On the other hand, in the previous section it was shown that an expansion in terms of Casimir operators, which are mainly related to quadrupole operators, leads to a satisfactory description of ground, single- and double-$`\gamma `$ bands. It is worth noting that such an expansion in terms of Casimir operators has been successfully used in Molecular Physics where spectroscopic data provide many anharmonic states . This leads us to propose a Hamiltonian as a quadrupole expansion that includes up to four-body terms. An alternative Hamiltonian can be based on an expansion in terms of pseudo-Casimir operators, which we define here as operators that become a true Casimir operator only for a particular choice of one structure parameter. For example, the $`\widehat{C}_2[SU(3)]`$ operator is related to the quadrupole operator through $`2\widehat{Q}\widehat{Q}+\frac{3}{4}\widehat{L}^2`$ where $`\widehat{Q}=s^{}\stackrel{~}{d}+d^{}\stackrel{~}{s}\pm \sqrt{7}/2(d^{}\times \stackrel{~}{d})^{(2)}`$. If the value $`\pm \sqrt{7}/2`$ is changed to $`\chi `$, a new operator $`\widehat{C}_2[SU(3)]_\chi `$ is obtained which we refer to as a pseudo-Casimir operator. It is a Casimir operator only for $`\chi =\pm \sqrt{7}/2`$.

Guided by the results of the previous section, two possible Hamiltonians that include up to four-body interactions, can be proposed, one based on a quadrupole expansion and the other on a pseudo-Casimir expansion:

$`\widehat{H}_Q`$ $`=`$ $`\kappa ^{}\widehat{L}\widehat{L}+a\widehat{Q}\widehat{Q}+b(\widehat{Q}\times \widehat{Q}\times \widehat{Q})^{(0)}+c(\widehat{Q}\widehat{Q})(\widehat{Q}\widehat{Q}),`$ (26)
$`\widehat{H}_{pC}`$ $`=`$ $`\kappa ^{}\widehat{L}\widehat{L}+a\widehat{C}_2[SU(3)]_\chi +b\widehat{C}_3[SU(3)]_\chi +c\widehat{C}_2[SU(3)]_\chi ^2,`$ (27)
where $`\widehat{Q}=\widehat{Q}(\chi )`$ and

$`\widehat{C}_2[SU(3)]_\chi `$ $`=`$ $`2\widehat{Q}(\chi )\widehat{Q}(\chi )+{\displaystyle \frac{3}{4}}\widehat{L}^2,`$ (28)
$`\widehat{C}_3[SU(3)]_\chi `$ $`=`$ $`4\sqrt{35}(\widehat{Q}(\chi )\times \widehat{Q}(\chi )\times \widehat{Q}(\chi ))^{(0)}{\displaystyle \frac{9}{2}}\sqrt{15}(\widehat{L}\times \widehat{L}\times \widehat{Q}(\chi ))^{(0)}.`$ (29)
For simplicity and taking into account the analysis done in the previous section, in the following $`b=0`$ is taken in equations (26) and (27).

### A Double-$`\gamma `$ band heads

There are three nuclei that have double-$`\gamma `$ bands identified without ambiguity: 164Dy, 166Er, and 168Er . In this section the band heads of these nuclei are studied using the Hamiltonians (26-27) to get an improved description of the anharmonicity phenomenon.

The number of bosons for 164Dy and 168Er is $`N=16`$ while for 166Er it is $`N=15`$. In the different calculations shown in this section the parameters of the Hamiltonian have been chosen as to reproduce as well as possible not only the heads of single- and double-$`\gamma `$ bands but also the structure of the bands. Because these calculations are schematic and do not try to be the best answer, the parameters of the Hamiltonian, for simplicity, will not be given fully and only in the case of final spectra will be shown. For comparison, the calculation with three-body terms (see section III and reference ) is also included.

In figures 7, 8, and 9 the heads of single- and double-phonon bands are shown. In each figure six panels are included: a) experimental data, b) calculation with Hamiltonian (26) and $`\chi =\sqrt{7}/2`$, c) calculation with Hamiltonian (26) and $`\chi =0.55`$, d) calculation with Hamiltonian (27) and $`\chi =\sqrt{7}/2`$, e) calculation with Hamiltonian (27) and $`\chi =0.55`$, and f) calculation with a three-body term $`(d^{}\times d^{}\times d^{})^{(4)}(\stackrel{~}{d}\times \stackrel{~}{d}\times \stackrel{~}{d})^{(4)}`$. In each panel, from left to right are represented the ground state, the $`\gamma `$ band head, the $`\beta `$ band head, the double-$`\gamma `$ $`K=0`$ band head, and the double-$`\gamma `$ $`K=4`$ band head. Due to the controversy on the nature of the $`\beta `$ band, several candidates for the latter have been included. For the same reason in the figures will be used the label “$`\beta `$”. In 164Dy and 168Er no information on the double-$`\gamma `$ $`K=0`$ band exists. The value $`\chi =0.55`$ is chosen as an alternative to the $`SU(3)`$ value because it describes very well the $`E2`$ transition probabilities in this mass region. Also, the same value of $`\chi `$ is taken in the Hamiltonian and in the transition operator, in line with the Consistent-Q Formalism (CQF) . In section V C this ansatz is used in the complete analysis of spectra and $`B(E2)`$ transitions for 164Dy, 166Er, and 168Er.

The most striking feature of figures 79 is that in all calculations the position of the double-$`\gamma `$ band heads and the degree of anharmonicity is well reproduced. Thus, the energies of the different band heads only are not sufficient to completely determine the Hamiltonian. Nevertheless, not all possible terms are able to create sufficient anharmonicity in the double-$`\gamma `$ bands. For example, in the case of three-body terms, a phenomenological study of the most relevant type of term shows that only $`(d^{}\times d^{}\times d^{})^{(4)}(\stackrel{~}{d}\times \stackrel{~}{d}\times \stackrel{~}{d})^{(4)}`$ is able to produce the required anharmonicity (see section III and reference ). On the other hand, only a few four-body interactions have been explored here and it cannot be excluded that other four-body terms can produce the appropriate degree of anharmonicity.

In order to decide which Hamiltonian is more appropriate, a description should be attempted not only of band heads but also of the structure of the bands and of $`E2`$ transition probabilities.

### B Moments of inertia

The study of the moments of inertia of the lowest bands is a very sensitive way to test the different calculations shown in figures 7, 8, and 9. Particular attention will be paid to the dynamic moment of inertia which can be obtained from the relation between angular momentum and $`\gamma `$-ray energy  and can be approximated by,

$$2\mathrm{}^2\frac{dJ}{dE_\gamma },$$
(30)
where $`J`$ is given dimensionless. In a plot of $`\gamma `$-ray energy versus $`J`$, $``$ will be the slope. Equation (30) can be used to study the structure of the rotational bands in comparison with experimental results.

In figure 10 the moments of inertia of the ground and $`\gamma `$-bands in the nuclei 164Dy, 166Er, and 168Er are compared to those obtained with the Hamiltonian of section III. The parameters for 164Dy are $`\kappa =24.2`$ keV, $`\chi =0.55`$, $`\kappa ^{}=6.0`$ keV, and $`\theta _4=51.0`$ keV, for 166Er are indicated in section III and for 168Er are $`\kappa =24.1`$ keV, $`\chi =0.55`$, $`\kappa ^{}=1.9`$ keV, and $`\theta _4=31.6`$ keV. The predicted moments of inertia disagree completely with the almost pure rotational structure observed experimentally. So, although this simple description based on one single three-body term describes well the anharmonic position of the double-phonon band heads, it fails in the moments of inertia of ground and $`\gamma `$-bands. A more realistic description needs others three-body terms. A recent analysis  shows, however, that Hamiltonians with two different three-body terms cannot get the correct moment of inertia.

The results with the Hamiltonian (26) are shown in figure 11. Two different calculations are shown which correspond to panels $`b`$ and $`c`$, respectively, of figures 79. In this case a better agreement is obtained but still some discrepancies remain, especially when a realistic value for $`\chi `$ is taken ($`\chi =0.55`$). Finally, in figure 12 the results with the Hamiltonian (27) are compared with the data. Again, two different calculations are shown which correspond to panels $`d`$ and $`e`$, respectively, of figures 79. Here, good agreement is found both for $`\chi =\sqrt{7}/2`$ and $`\chi =0.55`$.

The different results with the Hamiltonians (26) and (27), can be understood qualitatively by analyzing the structure of $`(\widehat{Q}(\chi )\widehat{Q}(\chi ))(\widehat{Q}(\chi )\widehat{Q}(\chi ))`$. With just two-body terms the Hamiltonians (26) and (27) are equivalent. This is different when up to four-body terms are included because $`(\widehat{Q}(\chi )\widehat{Q}(\chi ))(\widehat{Q}(\chi )\widehat{Q}(\chi ))`$ does not only contribute to $`\widehat{C}_2[SU(3)]_\chi ^2`$ and $`\widehat{L}^2`$ but also to $`(\widehat{L}^2)^2`$. This can be clarified with the equation

$`(\widehat{Q}(\chi )\widehat{Q}(\chi ))(\widehat{Q}(\chi )\widehat{Q}(\chi ))={\displaystyle \frac{1}{4}}\widehat{C}_2[SU(3)]_\chi ^2{\displaystyle \frac{3}{8}}\widehat{C}_2[SU(3)]_\chi ^2\widehat{L}^2+{\displaystyle \frac{9}{64}}(\widehat{L}^2)^2.`$ (31)
As a consequence, $`(\widehat{Q}(\chi )\widehat{Q}(\chi ))(\widehat{Q}(\chi )\widehat{Q}(\chi ))`$ substantially modifies the rotational structure of a band even in the case of pure $`SU(3)`$. This is how can be qualitatively understood that a description in terms of pseudo-Casimir operators is the most appropriate for dealing with anharmonic vibrations. In the next section a complete analysis is given of the nuclei under study in the framework of CQF using the Hamiltonian (27).

### C Realistic calculations

The complete calculated spectra of 164Dy, 166Er, and 168Er and the most relevant $`E2`$ transition probabilities are presented in this section. They are compared with existing data. In each calculation the same value of $`\chi `$ has been used both in the Hamiltonian (27) and in the electromagnetic operator (25). Finally, the effective charge in the transition operator (25), $`e_{\text{eff}}`$, was fixed for each nucleus to reproduce the $`B(E2;2_1^+0_1^+)`$ value.

The experimental and calculated spectra for 164Dy, 166Er, and 168Er are shown in figures 13, 14, and 15, respectively. The parameters used in the calculations are listed in table II. The parameters that yield the best fit to the energy spectra are very similar in the three cases which is consistent with the analysis carried out in the preceding sections. The overall description of the energies is satisfactory. The calculated low-lying bands are in good agreement with their experimental counterparts (see also figure 12) while the double-$`\gamma `$ band-head energies are close to the experimental values. The calculated ratios (8) are: $`R_0^\gamma =3.31`$ and $`R_4^\gamma =2.81`$ for 164Dy, $`R_0^\gamma =3.08`$ and $`R_4^\gamma =2.43`$ for 166Er, and $`R_0^\gamma =2.93`$ and $`R_4^\gamma =2.46`$ for 168Er, to be compared with the experimental ones: $`R_4^\gamma =2.84`$ for 164Dy, $`R_0^\gamma =2.82`$ and $`R_4^\gamma =2.45`$ for 166Er, and $`R_4^\gamma =2.50`$ for 168Er. Only for the $`\gamma _{K=0}^2`$ vibration the experimental and theoretical results are slightly different in the sense that this framework overestimates the anharmonic behavior for the $`\gamma _{K=0}^2`$ band. This can be corrected by increasing the value of $`|\chi |`$ in the Hamiltonian which, however, will introduce one more parameter because in the electromagnetic operator a different value of $`\chi `$ must be used.

For the calculation of $`E2`$ transition probabilities the same $`\chi `$ values as in the Hamiltonian are adopted. The effective charges are $`e_{\text{eff}}^2=(1.66)^2`$ W.u., $`e_{\text{eff}}^2=(1.83)^2`$ W.u., and $`e_{\text{eff}}^2=(1.67)^2`$ W.u., for 164Dy, 166Er, and 168Er respectively. In tables III, IV, and V the observed $`B(E2)`$ values and ratios concerning $`\gamma `$-vibrational states are compared with the theoretical results. In general, a good overall agreement is obtained in the three cases under study.

## VI Conclusions

In this paper the problem of anharmonicity in the $`\beta `$ and $`\gamma `$ vibrations of deformed nuclei was addressed in the context of the interacting boson model. The occurrence or not of anharmonicity was shown to be related to the order of the interactions between the bosons and the conclusions of the analysis can be summarized as follows. If the Hamiltonian includes up to two-body interactions, no sizeable anharmonicity can be obtained and the observed behavior cannot be obtained. The origin of this behavior is related to the existence of a first-order phase transition between rotational and vibrational nuclei  which excludes any interplay between one- and two-body terms, necessary to obtain anharmonic spectra. For up to three-body interactions that preserve $`SU(3)`$ symmetry it can be shown that the observed anharmonicity cannot be fully reproduced. Furthermore, extensive numerical calculations indicate that even a general IBM Hamiltonian that includes up to three-body interactions has difficulty in reproducing all observed aspects of ground, single- and double-$`\gamma `$ bands. However, due to the large parameter space of three-body interactions which is difficult to search exhaustively, this cannot be considered as a firm conclusion and alternative approaches (e.g. mean-field) should be tried to tackle the same problem. Finally, it was shown that a simple parametrization of the IBM Hamiltonian that includes up to four-body interactions can account for all observed properties of the three deformed nuclei with firmly established double-$`\gamma `$ vibrations. In particular, the introduction of $`SU(3)`$ pseudo-Casimir operators allows to describe the double-phonon states while keeping correct the properties of low-lying bands.

It should be emphasized that, in spite of its fourth-order character, the Hamiltonian considered here is only slightly more complex than the usual IBM Hamiltonian (one more parameter) and is a straightforward extension of the Consistent-Q Formalism that was previously successfully applied to many nuclei. Also clear from our study is the need for more experimental information about the double-phonon vibrations in deformed nuclei beyond the three cases known at present: our analysis indeed has shown that this information represents a challenging test of any theoretical description of deformed nuclei.

## Acknowledgements

We are grateful to C.E. Alonso, K. Heyde and, F. Iachello for valuable comments. Two of the authors (PVI and JEGR) wish to thank the Institute for Nuclear Theory, University of Washington, where this work was initiated. One of the authors (JEGR) thanks the F.W.O. for financial support. This work was supported in part by the Spanish DGICYT under project number PB98-1111.


format: text

text:
# Coherent radio emission from the magnetic chemically peculiar star CU Virginis

## 1 Introduction

Magnetic chemically peculiar (MCP) stars can produce radio emission at centimetric wavelengths. This radio emission is consistent with gyrosynchrotron emission from continuous ejected, mildly relativistic nonthermal particles trapped in the magnetosphere (Drake et al. drake (1987)).

MCP stars are characterized by mainly dipolar magnetic fields, whose strength steeply decreases with the stellar distance, and thus each radio frequency is expected to be emitted in a well localized stellar shell of the circumstellar region. Since the dipole axis is tilted with respect to the rotation axis, Leone (leo91 (1991)) suggested that the observed radio emission from MCP stars should be periodically variable. Combining their own observations with Drake et al. (drake (1987)) and Phillips & Lestrade (phil (1988)) data, Leone & Umana (leo93 (1993)) have shown that the 6 cm density fluxes of HD 37017 and HD 37479 vary with the stellar rotation period. The coincidence of radio maxima with the extrema of the effective magnetic field lead Leone & Umana to suggest that radio emitting regions are located above the magnetic poles. The variability of the 6 cm emission with the rotational period and has been also observed for the MCP star HD 133880 by Lim et al. (lim (1996)).

We have monitored along the rotational period CU~Vir (= HD 124224) with the Very Large Array (VLA) at four frequencies. This star is particularly suitable for studying the radio emission from MCP stars because of its small distance (80 pc), very short rotational period (0.52 days) and reversing magnetic field. Here we present the behavior of the 1.4 GHz radiation whose properties are compatible with cyclotron maser emission.

## 2 Observations and data reduction

We have monitored at 1.4, 5, 8.4 and 15 GHz the B9p Si star CU~Vir (= HD 124224) over three separate days using the VLA 111The Very Large Array is a facility of the National Radio Astronomy Observatory which is operated by Associated Universities, Inc. under cooperative agreement with the National Science Foundation. The observations have been carried out on June 1, 6 and 11, 1998 from 23 to 07 UT with all the available telescopes at the beginning of the A to B reconfiguration. For each frequency, we used the standard observing mode with two independent 50 MHz bands in Right and Left Circular Polarizations (RCP and LCP), adopting a 10 sec integration time. At 1.4 GHz the two bands are separated by 80 MHz, being centered at 1385 and 1465 MHz.

A typical observing cycle consisted of 10-min on source preceded and followed by 2-min on the phase calibrator 1354$-$021, which is only $`6\mathrm{°}`$ far from CU~Vir. The four frequencies have been observed alternately with all the telescopes, so that at each frequency CU~Vir has been observed for about $`1^\mathrm{h}20^\mathrm{m}`$ over a total time of $`8^\mathrm{h}`$. The sequence of the frequencies during the three observing runs have been organized to get the best possible sampling of the rotational phases, avoiding redundancies. To get a reliable flux scale, the amplitude calibrator 3C286 was observed at the beginning and the end of each run ($`40\mathrm{°}`$ and $`35\mathrm{°}`$ of elevation respectively).

Data were calibrated and mapped by using the standard procedures of the Astronomical Image Processing System (A.I.P.S.). CU~Vir was found unresolved on VLA baselines at all the observed frequencies. Its position was determined by using the task JMFIT, and found to be coincident with the position given in the Hipparcos Catalogue (Perryman et al. hip (1997)): $`\alpha (2000)=14^\mathrm{h}12^\mathrm{m}15\stackrel{s}{.}8`$ and $`\delta (2000)=+2\mathrm{°}24\mathrm{}34\stackrel{}{.}0`$. The temporal variation of the Stokes I and V parameters was determined with the task DFTPL. This task performs the direct Fourier transform of the visibilities as a function of time for an arbitrary position in the map. The results of this task can be affected by the sidelobes of any other strong source in the field. However, since the only source close to our target is very weak, having a flux density of only 3 mJy at 1.4 GHz, the possible confusion is negligible. Heliocentric correction was then applied to the times of observation. Fig. 1 shows the radio light curve at 1.4 GHz during the three days of observation.

## 3 Stellar parameters

Spectral, photometric and magnetic variations of MCP stars are characterized by a single period. Stibbs (stibbs (1950)) proposed that these stars present a dipolar magnetic field, whose symmetry axis is tilted with respect to the rotational axis, and a non homogeneous distribution of chemical elements on the stellar surface. Thus the variability period is the rotational one.

The rotational period of CU~Vir was firstly determined by Deutsch (deutsch (1952)) equal to 0.52067 days. Several authors have later re-determined this period with values ranging from 0.52067 to 0.5207 days, see Catalano & Renson (cat (1998)) for a complete list of references. Recently, Pyper et al. (pyper (1998)) analyzed all the photometric data covering 40 years. They suggested that the rotational period of CU~Vir increased abruptly between 1983 and 1987 of a factor $`\mathrm{\Delta }P/P4.9\times 10^5`$. This unexpected jump of the rotational period was discussed by Stȩpień (stepien (1998)). He suggested that the spin down occurs in the outher stellar envelope alone, due to a change of the moment of inertia of the envelope itself because of the internal magnetic field strength. The new ephemeris given by Pyper et al. (pyper (1998)), referred to light minimum, are:

$`HJD=2435178.6417+\{\begin{array}{cc}0^d.52067780E\hfill & \hfill \text{JD }<\text{ 2\hspace{0.17em}446\hspace{0.17em}000}\\ 0^d.52070308E\hfill & \hfill \text{JD }>\text{ 2\hspace{0.17em}446\hspace{0.17em}000}\end{array}`$

Using the projected rotational velocity $`v_\mathrm{e}\mathrm{sin}i=146\pm 2\mathrm{km}\mathrm{s}^1`$ (Hatzes hatzes (1997)) and the stellar radius $`R_{}=2.2R_{\mathrm{}}`$ (North north (1998)), in the relation

$$\mathrm{sin}i=P_{\mathrm{days}}\frac{v_\mathrm{e}\mathrm{sin}i}{50.6}(R_{}/R_{\mathrm{}})^1,$$

valid for a rigid rotator, we measure an inclination of the rotational axis with respect to the line of sight $`i=43\mathrm{°}\pm 7\mathrm{°}`$.

Borra & Landstreet (borra (1980)) found that the effective magnetic field ($`B_{\mathrm{eff}}`$), measured from 1976 to 1978, is variable with the period established by Winzer (winz (1974)) from photometry. $`B_{\mathrm{eff}}`$ is the average of the component of the local magnetic field along the line of sight over the whole stellar disk. From the value of the inclination $`i`$ and Borra & Landstreet (borra (1980)) measurements of $`B_{\mathrm{eff}}`$, we obtain an obliquity of the dipole $`\beta =74\mathrm{°}\pm 3\mathrm{°}`$ and a polar magnetic field equal to $`B_\mathrm{p}=3000\pm 200`$ gauss. Pyper et al. (pyper (1998)) give new values of $`B_{\mathrm{eff}}`$; however, we prefer to infer the magnetic geometry from Borra & Landstreet measurements which are based on hydrogen lines. In fact, Pyper and co-workers measurements are based on the 634.7 and 637.1 nm lines of silicon and this element is not homogeneously distributed on the surface of CU~Vir (Kusching et al. kusching (1999)). Even if we cannot exclude a more complex topology of the magnetic field, like a dipole+quadrupole or a decentered dipole, as proposed by Hatzes (hatzes (1997)), the assumption of a simple dipole at an height above the photosphere where the radio emission is generated seems reasonable, as the strength of the quadrupole field decreases one order of magnitude faster than the dipole field.

A summary of the adopted and derived parameters for CU~Vir is reported in Table 1.

## 4 The 1.4 GHz emission

The observed 1.4 GHz emission is shown in Fig. 2 as a function of the phase during the three days of observation separately. In Fig. 3 a summary of the behavior of other three frequencies as a function of the phase for the whole observing run is reported. The smooth variations of the radio flux and of the circular polarization at the higher frequencies will be discussed elsewhere, together with a numerical model for the continuum radio emission from MCP stars.

### 4.1 Correlation with the magnetic field

Since Leone & Umana (leo93 (1993)) reported a correlation between radio emission and effective magnetic field in HD 37017 and HD 37479, we look first at a similar correlation for our data. Fig. 4, upper panel, shows the 1.4 GHz emission as a function of the rotational phase. This emission is characterized by two or three components: a basal flux of $`23`$ mJy and very large increments of the flux around phases 0.35-0.45 and 0.75-0.85, where it goes up to 15 mJy, and secondary peaks of up to 7 mJy, further discussion of which is deferred to Section 5.2.2. Only in two moments, at phase 0.85 and 0.97, the intensity drops below 1 mJy. Phasing the magnetic data from Borra & Landstreet (borra (1980)) and Pyper et al. (pyper (1998)) (Fig. 4, lower panel), we found that those peaks of the radio emission coincide approximately with the the null magnetic field, i.e. when the axis of the dipole is almost perpendicular to the line of sight. The peak visible at phase 0.35-0.45 is defined by the observations on June 2 (maximum) and 6 (rising phase), and the 0.75-0.85 one on June 2 (rising phase) and 11 (maximum) (see also Fig. 2). According to the accuracy of the period given by Pyper et al. (pyper (1998)), our data are phased with respect to the magnetic data better than $`3\times 10^3`$.

It is worthy to note that the flux density of the peaks is about five times larger than those previously reported in the literature at this frequency for MCP stars. Previous observations of CU~Vir reported by Leone et al. (leo96 (1996)) show that the 1.4 GHz flux is 2.6 mJy, in agreement with the ”out of peaks” emission here reported.

The high flux increment, that occurs at the particular orientation of the magnetosphere, indicates an emission mechanism that is not explainable with the emission models up to now proposed. In the following, we will analyze further characteristics of this emission.

### 4.2 Polarization

A further exceptional aspect of the 1.4 GHz is its high degree of circular polarization. Fig. 5, where the Stokes V parameter (V = 1/2(RCP$``$LCP)) is plotted versus the rotational phase, shows that the increment of the flux occurs only in the right-hand circular polarization. Six peaks, denoted with the letters from a to f in Fig. 5, are visible, the largest and broad emission being detected in peaks a and d, that occur in coincidence with the null magnetic field. Outside the peaks, the Stokes parameter V is statistically null. The percentage of polarization $`\pi _\mathrm{c}`$ goes up to 80 % during the main peaks a and d. If we subtract the continuum radio emission, that contributes for about 3 mJy, we get $`\pi _\mathrm{c}100`$ %.

In order to verify if the behavior of the right-hand circular polarization was not due to an instrumental problem, we looked for other sources in the field of our star to be monitored for the polarization at the same times of our observations. Leone et al. (leo96 (1996)) found in the 6 cm VLA frames centered at the CU~Vir position a radio source with coordinates: $`\alpha (2000)=14^\mathrm{h}12^\mathrm{m}24\stackrel{s}{.}96`$ and $`\delta (2000)=+2\mathrm{°}22\mathrm{}04\stackrel{}{.}7`$, that they used to check possible instrumental effects. This field source is visible also in the 1.4 GHz map showed in Fig. 6 and has a flux density (Stokes I) of about 3 mJy, comparable with CU~Vir outside the peaks. By comparing the polarizations (Stokes V) of the field source and of CU~Vir, we can rule out any instrumental effects (Fig. 7).

### 4.3 Directivity of the radio emission

Since the polarized peaks have been observed in the three different days of observation, we can conclude that this radiation is stable at least in a period of weeks. The observed variations can be due to the different inclination that the oblique dipole forms with the line of sight as the star rotates. In a particular geometric configuration, highly beamed radiation is emitted toward the Earth, producing the observed peaks. As already pointed out, the main peaks of the right-hand circular polarization corresponds almost to the zero of the average magnetic field over the surface of the star. This means that the maxima of this emission occur when the magnetic axis is almost perpendicular to the line of sight.

Borra & Landstreet (borra (1980)) found that the magnetic curve is delayed with respect to the light curve. Assuming that the two main peaks occur in two symmetric orientations of the magnetosphere, we can compute this shift with good accuracy: from our data, peaks $`a`$ and $`d`$ show their maxima at phases $`\varphi _\mathrm{a}=0.37`$ and $`\varphi _\mathrm{d}=0.80`$ respectively, and therefore

$$\varphi _0=(\varphi _\mathrm{a}+\varphi _\mathrm{d})/20.5=0.08.$$

Once the geometry of the magnetosphere of the star has been defined by the values of the inclination $`i`$ and the obliquity $`\beta `$ (see Table 1), we can find the angle $`\theta _\mathrm{M}`$ that the line of sight forms with respect to the axis of the dipole, defined by:

$$\mathrm{cos}\theta _\mathrm{M}=\mathrm{sin}\beta \mathrm{sin}i\mathrm{cos}(\varphi \varphi _0)+\mathrm{cos}\beta \mathrm{cos}i.$$

The error of $`\theta _\mathrm{M}`$ estimated from this relation and from the uncertainty in the angles $`i`$ and $`\beta `$ is about $`4\mathrm{°}`$ when $`\theta _\mathrm{M}90\mathrm{°}`$ and $`8\mathrm{°}`$ when $`\theta _\mathrm{M}30\mathrm{°}`$ and $`120\mathrm{°}`$. Fig. 8 shows the change of the orientation of the star, and so of the angle $`\theta _\mathrm{M}`$, with the rotation.

Fig. 9 shows the behavior of the polarized component of the emission as a function of the angle $`\theta _\mathrm{M}`$ in a polar view. The main peaks are emitted at an angle of about $`85\mathrm{°}`$ with respect to the axis of the dipole (indicated as the arrow B), and have a high degree of directivity, with an half power beam width of about $`5\mathrm{°}`$. For the less intense peaks $`b,c,e`$ and $`f`$, the beams have a narrower width of about $`1\mathrm{°}`$. Fig. 10 shows the details of the main peaks as a function of the angle $`\theta _\mathrm{M}`$. It is important to note that the magnetic longitude of the line of sight during the two main peaks are different, since in these configurations the star shows opposite hemispheres (see Fig. 8). This means that the observed phenomenon is not related to an active longitude, depending only on $`\theta _\mathrm{M}`$. It seems to be similar to the pulses observed in pulsars, the magnetic fields of which are, like MCP stars, typically characterized by oblique dipole geometries.

## 5 Discussion

To explain the radio continuum and the X-ray emissions from young magnetic B stars and Bp-Ap stars, André et al. (andre (1988)) proposed a model where the stellar wind plasma flows out near the magnetic poles along the field lines. Far from the star, such a wind draws the field lines near the equator into current sheets which should be location of particle acceleration. Mildly relativistic electrons return to regions near the star by traveling along the magnetic field lines, emitting gyrosynchrotron radiation. This model has been applied by Linsky et al. (linsky (1992)) to explain the flat spectra of the MCP stars.

The characteristics of the polarized component of the radio emission at 1.4 GHz from CU~Vir, that show high degrees of circular polarization and high directivity, can be explained in terms of coherent emission. The two major mechanisms that have been suggested to explain coherent radio emission are plasma radiation due to Langmuir waves and cyclotron maser. Plasma radiation has been invoked in several cases, such as solar microbursts at 1.4 GHz (e.g. Bastian bastian (1991)), solar millisecond spikes (e.g. Wentzel wentzel (1993)) and radio bursts from the flare star AD~Leo (Abada-Simon et al. abada (1997)). The theory of Electron Cyclotron Maser Emission (ECME) (e.g. Wu & Lee wu (1979) or Melrose & Dulk meldulk (1982)) seems to be the favourable emission mechanism to explain coherent radiation from the magnetosphere of Jupiter, solar spike bursts and flare stars, like M dwarfs and RS CVn binary systems.

### 5.1 Plasma Radiation

Plasma radiation is a two-stage process where longitudinal waves in the plasma (Langmuir waves) are first generated and later their energy is converted into radiation (e.g. Dulk dulk (1985)). The frequency of the radiation is the plasma frequency ($`\nu _\mathrm{P}\mathrm{9\hspace{0.17em}000}\sqrt{N_\mathrm{e}}`$ Hz, with $`N_\mathrm{e}`$ the plasma density number in $`\mathrm{cm}^3`$) or its second harmonic. To be observed at 1.4 GHz, $`N_\mathrm{e}2.4\times 10^{10}`$ or $`6\times 10^9\mathrm{cm}^3`$ respectively for $`s=1,2`$. Plasma radiation can occur when the magnetic field is relatively weak ($`\nu _\mathrm{P}\nu _\mathrm{B}`$, where $`\nu _\mathrm{B}=2.8\times 10^6B`$, with $`B`$ in gauss). This gives an upper limit for the magnetic field $`B`$ in the region where the radiation is generated: $`B<500`$ ($`s=1`$) or $`B<250`$ gauss ($`s=2`$). Assuming a dipolar topology of the magnetic field, $`B=B_{}(R_{}/R)^3`$, with $`B_{}=\mathrm{3\hspace{0.17em}000}`$ gauss, we get that the region where plasma emission occur must be located at $`R>1.8R_{}`$ ($`s=1`$) or $`R>2.3R_{}`$ ($`s=2`$).

However, the theory of plasma radiation does not foresee any high directivity. For this reason this mechanism is not suitable to explain our observations.

### 5.2 Electron Cyclotron Maser Emission

Following the theory of the ECME, electrons reflected by the magnetic mirrors can develop a pitch angle anisotropy (or loss cone anisotropy), becoming candidates for cyclotron maser emission if the local plasma frequency is relatively small ($`\nu _\mathrm{P}\nu _\mathrm{B}`$). The frequency of the maser emission is given by $`\nu s\nu _\mathrm{B}`$, being $`s`$ the harmonic number and $`\nu _\mathrm{B}`$ the gyrofrequency. The faster growth rate is for the first few harmonic number ($`s4`$) in the extraordinary mode (x-mode). The radiation generated where the magnetic field intensity is $`B_0`$ has a frequency $`\nu _0s\nu _B`$; when crossing a more external layer with $`B_1=\frac{s}{s+1}B_0`$, it can be suppressed by the gyromagnetic absorption of the thermal plasma. Melrose & Dulk (meldulk (1982)) show that the $`1^{\mathrm{th}}`$ harmonic is generally suppressed, while the $`2^{\mathrm{th}}`$ can escape from layers absorbing at higher harmonic numbers. The ECME is confined in a hollow cone of half-angle $`\theta `$ with respect to the line of the magnetic field, with $`\mathrm{cos}\theta =v/c`$ and $`v`$ the speed of the emitting electrons. The thickness of the hollow cone is $`\mathrm{\Delta }\theta v/c`$. If the maser is emitted in a region of constant magnetic field, the relative bandwidth is $`\mathrm{\Delta }\nu /\nu =\frac{\nu s\nu _\mathrm{B}}{s\nu _\mathrm{B}}\mathrm{cos}^2\theta `$. However, if the emission comes from a layer where the magnetic field ranges from $`B_\mathrm{a}`$ to $`B_\mathrm{b}`$, then the observed bandwidth will be much larger, ranging the radiation from $`s\nu _\mathrm{a}`$ to $`s\nu _\mathrm{b}`$. In this case, the angle $`\theta `$ and the thickness of the hollow cone of radiation $`\mathrm{\Delta }\theta `$ do not change, depending only on $`v/c`$.

#### 5.2.1 The main peaks

We observe that the two main peaks $`a`$ and $`d`$ are beamed at an angle $`\theta _\mathrm{M}85\mathrm{°}`$ with respect to the axis of the dipole, and have a full width half maximum (FWHM) of about $`5\mathrm{°}`$ (Fig. 10). This strongly suggests that we are in presence of ECME.

We propose the following scenario: the electrons accelerated in the current sheets out of the Alfvén radius flow toward the photosphere close to the magnetic pole; they are eventually mirrored back by the increasing magnetic field, traveling along field lines almost parallel to the axis of the dipole. After the reflection, they develop a loss cone anisotropy because of the interaction with the thermal plasma, leading to electron cyclotron maser emission. In this hypothesis the angle $`\theta _\mathrm{M}`$ is just the angle $`\theta `$ of the ECME theory. Since $`\mathrm{\Delta }\theta v/c`$ and $`\mathrm{cos}\theta v/c`$, we get $`v/c0.09`$.

If the magnetic field of CU~Vir is a dipole, we should expect a symmetry in the ECME with the stellar phase. The absence of any beamed emission at $`\theta _\mathrm{M}95\mathrm{°}`$ (i.e. at $`85\mathrm{°}`$ from the direction of the south pole) suggests an asymmetry of the magnetosphere. The presence of a quadrupole component, that has been observed in other MCP stars like HD~32633 and HD~175362 (Mathys mathys (1991)), or of a decentered dipole, as recently proposed for CU~Vir (Hatzes hatzes (1997)), could explain the observed asymmetry of the ECME. In fact, the presence of a quadrupole, as for example shown by Michaud et al. (michaud (1981)) in their figure 1d, can inhibit the wind, and so the radio emission, from the magnetic south pole. The interpretation of an asymmetry in the wind is not new. In fact Brown et al. (brown (1985)), from the behavior of the UV lines of the MCP star HD~21699, inferred that the wind flows from ”only one of the magnetic poles”.

#### 5.2.2 The secondary peaks

The secondary peaks $`b`$, $`c`$, $`e`$ and $`f`$ are, like the main peaks $`a`$ and $`d`$, also circularly polarized. The maximum flux density that they show is about 3 mJy in the Stokes V, as it is possible to see inspecting Fig. 12. Peak $`b`$ lasts for 2 minutes (half power to half power) and a further rise is possible. It has been detect on June 6 at phase 0.5; at the same phase, on June 11, no flux enhancement has been observed, meaning that probably peak $`b`$ is a transient phenomenon. Peak $`c`$ lasts about 3–4 minutes and is followed by the rise of the main peak $`d`$. Peak $`e`$ is very short in duration, about one minute, and peak $`f`$ is about 4 minutes. While the main peaks $`a`$ and $`d`$ are detected in all the three days of observation, we cannot say if $`b`$, $`c`$, $`e`$ and $`f`$ are sporadic impulsive emissions or if they are stable as the main ones. Further observations are needed to clarify this point.

#### 5.2.3 The bandwidth

Our observations have been performed at two bands of 50 MHz separated by 80 MHz. No difference of flux between the two bands has been found. So, the bandwidth of the masing radiation is $`\mathrm{\Delta }\nu 80\mathrm{MHz}`$. Since in the ECME theory $`\mathrm{\Delta }\nu /\nu \mathrm{cos}^2\theta `$, with $`\theta =85\mathrm{°}`$ we expect a bandwidth of about 10 MHz. This apparent incongruence can be explained if the region where the maser emission is generated covers a wide range of magnetic field strength. The observed radiation is the envelope of a continuous series of maser spots along the field lines, the higher frequency being emitted in regions closer to the star, according to $`\nu s\nu _\mathrm{B}`$. In a dipole $`B=B_{}(R_{}/R)^3`$; with $`B_{}=\mathrm{3\hspace{0.17em}000}`$ gauss, $`\nu 8.4s(R_{}/R)^3`$ GHz. To be observed at 1.4 GHz, the maser spots are located at $`R1.8s^{1/3}`$, with $`s4`$. On the contrary, we did not observed any coherent emission at 5 GHz, as we will discuss in a following paper. So, no condition for the maser mechanism is expected close to the star. The maser mechanism is efficient at a distance $`R>1.5s^{1/3}`$. Probably the electrons are thermalized close to the star, due to a higher density.

#### 5.2.4 Why only right hand circular polarization?

The Stokes parameter V is always positive, that means the radiation is right hand polarized. The theory of the ECME foresees that the radiation is almost entirely polarized in X-mode, as observed for the auroral kilometric radiation (AKR) and for the Jupiter’s decametric emission (DAM) (Wu & Lee wu (1979), Melrose melrose (1976)). In fact, in the x-mode the sense of rotation of the electric vector is the same as the helicity of the emitting particles. Electrons moving in a magnetic field directed toward us are seen to rotate in counter clockwise, as the right hand circular polarization. In a perfectly symmetric configuration, the electrons mirrored outward move in the same direction as the magnetic field lines in the north hemisphere, in opposite direction in the south hemisphere, emitting respectively in right and left hand polarization.

Our data show no LCP enhancement at any rotational phase. This suggests that there is no condition for cyclotron maser emission at 20 cm in the magnetic south hemisphere. Again, this can be imputed to an asymmetry in the magnetosphere of CU~Vir.

## 6 Conclusions

The discovery of the coherent radio emission from CU~Vir, that can be explained in terms of Electron Cyclotron Maser Emission, opens new prospectives in the field of the magnetic chemically peculiar stars.

An important question that can be put is: is the coherent radio emission at low frequency a common characteristic of the magnetic chemically peculiar stars? Low frequency observations of a sample of MCP stars that present a reversal magnetic field have to be carried out in order to answer this question. Such a survey will certainly put more stringent parameters to draw the physical mechanisms that act in the magnetosphere of the MCP stars.

Further observations on CU~Vir are needed to confirm that the coherent radio emission here reported is stable in time. In particular, a continuous coverage with the rotational phase to better define the angular distribution of the coherent emission is crucial. Moreover, further observations covering several rotational periods are also needed to study the nature of the secondary peaks and to determine whether or not the phases at which they occur are truly random.

The bandwidth of the radiation is found to be larger than 80 MHz. Observations at frequency lower than 1400 MHz are needed to better define the spectrum.

###### Acknowledgements.

We thank the referee, Dr. S. Drake, for comments and suggestions which enabled us to improve this paper.


text:
# Rolling up of graphite sheet: Energetics of shell formation

## I  Energy of curved graphite surface

Three phenomenological parameters were used  to calculate the additional energy of CNC formation, as compared with the known specific energy of an infinite graphite sheet (graphene). The energy of any carbon cluster with a curved surface is decomposed into additive terms. Each term is given by some characteristic energy (the model parameter) multiplied by a number depending solely on the cluster geometry. This essential simplification of the calculation is based on natural reasoning. There are 3 sources bringing at least 3 parameters into the energy (excepting the ”zero energy” of an atom belonging to an infinite graphite sheet; this term is excluded from further expessions). The first additional term comes from dangling bonds. The corresponding energy parameter, $`E_b`$, is the dangling bond energy. An extra energy due to a change of the degree of a hybridisation of an electron orbit on the curved surface is described using the second parameter, $`E_c`$. The last parameter, $`E_5`$, is the total energy of twelve pentagonal rings. This parameter reflects that an electron belonging to the pentagon must have energy different from that in a graphite-like hexagon. Hence, it partly takes into account non-equivalent bonds. Note that this is a topologically determined donation to the total energy of any closed cluster .

In general, the curvature of a surface is a function of a point. The curvature is constant and equal to $`3/R^2`$ for a sphere of radius $`R`$, and $`9/8R^2`$ for a tube of the same radius. For the sake of clarity, we consider here only the surface with constant curvature. One easily obtains the formation energy of an infinite carbon tube. The specific energy of dangling bonds is negligible since it goes to zero with the tube size going to infinity. A tube has no pentagon, because a cylinder is topologically equivalent to a plane. Hence, the energy has only one term, the ”curvature energy” $`E=E_c\mathrm{\hspace{0.33em}9}N/8R^2`$, where $`N`$ is the number of atoms and $`R`$ is the dimensionless tube radius. Henceforth, all lengths will be measured in bond lengths ($`b1.4`$Å not supposed to vary according to the model). We compared the tube energy with the relevant results of other computations . This gave the first parameter of the model: $`E_c0.9`$ eV (see also ).

Most likely, the energy of a dangling bond in CNC does not differ too much from the value for graphite, and therefore, we took $`E_b=2.36`$ eV. The total energy associated with the dangling bonds is proportional to the cluster perimeter length (length of the open boundary of the CNC lattice) multiplied by the density of dangling bonds along the perimeter, which in turn depends on the local geometry. One can imagine, for example, a finite tube as a parallelogram, carved from hexagonal graphene, which has its opposite sides glued . This parallelogram has a variable number of sites on the open perimeter $`\zeta `$ (dangling bond density), depending on its orientation relative to the lattice vector. None the less the density varies in a narrow region. If we take the geometrical factor in $`\zeta `$ to be unity for the ”zigzag” tube, then it becomes $`2/\sqrt{3}`$ for the ”armchair” tube; one can easily calculate $`\zeta `$ for the bond density of a specific tube lying between these extreme examples.

It was shown in our previous papers , that there occurs competition  in the formation energy of a finite tube between two terms with $`E_c`$ and $`E_b`$ (the energy associated with dangling bonds can not be neglected for a tube of finite length). The tube energy reads as: $`E(N,R)=E_c\mathrm{\hspace{0.33em}9}N/8R^2+E_b\mathrm{\hspace{0.33em}4}\pi \zeta R/\sqrt{3}`$. It is advisable to make the total perimeter $`4\pi R`$ shorter in order to diminish the number of dangling bonds. The radius decrease costs an extra total energy owing to stronger curvature. We called the cluster having the minimal energy at a certain number of atoms ”an optimal cluster”. The optimal cluster is the configuration governed by energy consideration at any fixed $`N`$. (The nanotube energetics and the tube formation road from a flat graphene were discussed in details in Ref.). Thus, we constructed a phenomenological approach and predicted the optimal cluster shape, based on the energetics of any nanotube without performing a quantum-chemical calculation each time.

The third parameter of the model was determined by fitting the formation energy of a spherical cluster to the experimental value for $`C_{60}`$. The total energy of a sphere consists of two components . The first constant component is the energy of 12 pentagons. The sphere energy also includes the total curvature energy which is independent of the radius in this case (because the number of atoms being proportional to squared curvature, $`1/R^2`$). The second component is a correction excluding pentagonal bonds from the total count of curved bonds. Therefore, the energy is as follows:

$$E_{\mathrm{sph}}=\left(E_5+\frac{16\pi E_c}{\sqrt{3}}\right)\frac{N_s}{N}E_c$$
(1)
where we introduce a characteristic number $`N_s=2\times 60\times 16\pi /3\sqrt{3}1161`$. The number $`N_s`$ includes the number of bonds belonging to the pentagonal rings in the fullerene (cf. also ) and the (dimensionless) curvature of the surface times the (dimensionless) area of the sphere covered by hexagons. We took the experimental value of C60 formation energy from Ref. and fitted our last parameter $`E_517.7`$ eV, which takes into account all 12 pentagons of the spheroidal closed cluster. The sphere energy evidently remains positive (for the given model parameters $`E_c`$=0.9 eV and $`E_5`$) for the smallest spherical cluster considered here, $`C_{60}`$. Then, comparing the experimental formation energy of $`C_{70}`$ with our calculation, we find that the difference is about 1% and all three parameters are self–consistent.

## II How sphere is rolled up

A sphere has the minimal curvature energy between closed nanoclusters of any definite number of atoms, $`N=`$const. It seems interesting to compare the energies of closed clusters with the energy of an open fragment having dangling bonds. For an open CNC, the shorter the perimeter (number of dangling bonds), the stronger the curvature. The decreasing energy of dangling bonds makes favourable a small rolled–up cluster with a shorter perimeter. The curvature energy demands that the CNC be extended and flat. As a result of the energy competition, the system decreases its total energy via elimination of dangling bonds owing to the fact that the infinitely large sphere is a configuration resulting in the global minimum of the total CNC energy .

We consider below a sphere (of area $`S_o`$) with a round hole of angular size $`\mathrm{\Omega }`$ varying from 0 to $`4\pi `$, from the sphere to the round piece of the plane (see Fig.1). This means that we choose a flat CNC with the minimal perimeter (round) and a closed CNC with the minimal curvature (sphere) at fixed area or number of atoms. We name an intermediate cluster ”an open sphere”. The curvature energy decreases with the open sphere surface area as $`(1\mathrm{\Omega }/4\pi )S_o`$. This can be seen from the expression for the number of atoms in an ideal graphene sphere with a round hole: $`N=(1\mathrm{\Omega }/4\pi )16\pi R^2/3\sqrt{3}`$. It is also natural to feather out the dependence of the topological energy $`E_5(\mathrm{\Omega })`$ in the hole size $`\mathrm{\Omega }`$ (actually varying stepwise, each step correspondes to creating a pentagonal defect) and to substitute a linear dependence. Instead of choosing a specific way to create a defect and to place it in the cluster lattice, the simplest uniform distribution of pentagons over the open sphere surface is used here . This is nearly equivalent to placing pentagons as far apart as possible. The kinetics of the process  is not touched on in this paper, as well as the pentagon–pentagon distance optimised elsewhere .

In terms of the independent variables $`R`$ and $`\mathrm{\Omega }`$ (or $`N`$ and $`\mathrm{\Omega }`$ alternatively) the energy of the open sphere reads as:

$`\begin{array}{c}E(N,\mathrm{\Omega })=\left(1{\displaystyle \frac{\mathrm{\Omega }}{4\pi }}\right)\left(E_5+E_c\left[{\displaystyle \frac{16\pi }{\sqrt{3}}}{\displaystyle \frac{120}{R^2}}\right]\right)+{\displaystyle \frac{4\pi RE_b}{\sqrt{3}\zeta }}\sqrt{{\displaystyle \frac{\mathrm{\Omega }}{4\pi }}}\sqrt{1{\displaystyle \frac{\mathrm{\Omega }}{4\pi }}}=\\ =\left(1{\displaystyle \frac{\mathrm{\Omega }}{4\pi }}\right)\left[E_5+{\displaystyle \frac{16\pi E_c}{\sqrt{3}}}E_c{\displaystyle \frac{N_s}{N}}\left(1{\displaystyle \frac{\mathrm{\Omega }}{4\pi }}\right)\right]+2\left[E_5+{\displaystyle \frac{16\pi E_c}{\sqrt{3}}}\right]{\displaystyle \frac{1}{\sqrt{N_{\mathrm{th}}}}}\sqrt{{\displaystyle \frac{N\mathrm{\Omega }}{4\pi }}}\end{array}`$ (4)
where $`\zeta `$ is a geometrical multiplier of about unity , and we introduce $`N_{\mathrm{th}}`$, a ”threshold” size of planar cluster, given by: $`{\displaystyle \frac{\sqrt{\pi \sqrt{3}}}{\zeta }}{\displaystyle \frac{E_b}{2[E_5+\frac{16\pi }{\sqrt{3}}E_c]}}=\sqrt{N_{\mathrm{th}}^1}`$. Its meaning is discussed just below. The first term in Eq.(4) is the energy of an incomplete sphere. The second is the dangling bond energy. Of course, the accuracy of the model is poorer for small $`N`$ and strong curvature . As it is seen from Figure 2, the energy of the open sphere with hole size of neither 0 nor $`4\pi `$ has a maximum at fixed number of atoms when $`NN_{\mathrm{th}}254`$. Partial differentiation of Eq.(4) with respect to $`\mathrm{\Omega }`$ (at fixed $`N`$) gives an analytical expression for the size of the open sphere with the maximal energy. The function is given (with high accuracy for $`N>100`$) by:

$$\mathrm{\Omega }_B4\pi \left(1\left(1\frac{N}{N_{\mathrm{th}}}\right)\frac{N}{N+4N_s/\left[{\displaystyle \frac{E_5}{E_c}}+{\displaystyle \frac{16\pi }{\sqrt{3}}}\right]}\right).$$
(5)
The line given by Eq.(5) intersects the abscissa axis at $`\mathrm{\Omega }=4\pi `$ and $`N_{\mathrm{th}}`$, named ”the threshold flat cluster” in the following sense.

Let us consider a set of CNCs with $`N<N_{\mathrm{th}}`$ and $`4\pi <\mathrm{\Omega }<0`$: the energy of a cluster with $`N=\mathrm{const}`$ grows in moving away from both extremities of this interval and has a maximum at $`\mathrm{\Omega }=\mathrm{\Omega }_B`$. That is, going from $`\mathrm{\Omega }=4\pi `$ to $`\mathrm{\Omega }=0`$ (where the minimum of the energy occurs at fixed $`N`$), we pass through a barrier of height $`\delta E(N)=E(N,\mathrm{\Omega }_B)E(N,4\pi )`$ decreasing with $`N`$ and tending to zero at $`N_{\mathrm{th}}`$. Therefore, as follows from energy consideration solely, planar clusters with $`NN_{\mathrm{th}}`$ can be rolled up without any barrier. Below $`N_{\mathrm{th}}`$, the flat CNC is metastable to the rolling–up. Suppose we have initially some distribution of planar graphite fragments $`f(N)`$ (which is consistent with the two–step–model of carbon cluster formation proposed in  where nucleation and growth processes were separated). It is energetically favourable to increase the number of atoms, but this process depends on availability of the building material (carbon atom concentration). One can suppose that in the second step of synthesis the only conformation is possible between different clusters (even through some exchange of carbons).

We consider below that the change in number of atoms is suppressed. Suppose that the most energetically favourable process is rolling–up of a flat CNC into a sphere . Then clusters with $`N>N_{\mathrm{th}}`$ are rolled–up for energetical reason. At some finite temperature $`T`$, the probability to roll up a smaller cluster is exponentially small: $`W\mathrm{exp}(\delta E/T)`$, where $`\delta E`$ is the barrier height (shown in Figure 2) depending on $`N`$ and $`\mathrm{\Omega }`$. This is the difference between the energy of the given flat cluster and that of a ”barrier” cluster with the same number of atoms.

As a result of rolling–up, the initial size distribution of planar clusters is changed. We discussed above only a single path of conformation (rolling–up into the spherical-shell cluster). For this reason we can only speculate about the distribution of the rest of planar clusters. It has a maximum depending on temperature and characteristic number of atoms in the initial distribution $`f(N)`$, being determined by the conditions of two distinct synthesis regimes and, therefore, treated as independent .

In summary, we proposed a heuristic model for calculating the energy of a carbon nanocluster formation. Within the model, we use only three phenomenological parameters extracted from computer simulations and experimental data. The model allows evaluating the formation energies of various CNCs and finding the preferable cluster shape.

We compared the formation energies of the most preferable closed cluster (spherical cluster) and the most preferable flat fragment (round piece of graphene). We conclude that for energy reasons rolling–up of a flat fragment is always favourable, since the energy of a closed cluster is the lowest. We predict that a small enough flat graphite cluster can be metastable with respect to rolling into a shell. The energy barrier for rolling is calculated and shown to disappear for the clusters containing more than three hundred atoms.

Acknowledgments.  This work was partially supported by RFBR grants no. 96–15–96348 and 98–02–18117.

The supposed path of spherical-shell closed cluster formation: from a flat round fragment of graphene, $`\mathrm{\Omega }=4\pi `$ (right), via open spheres with $`0<\mathrm{\Omega }<4\pi `$ (center), to a shell, $`\mathrm{\Omega }=0`$ (left). The number of atoms is kept constant. Therefore, the linear size decreases.

The energy landscape of the process of rolling–up. The barrier height for rolling–up is shown. The energy of a flat round fragment of graphene is lower than that of an intermediate state, the open sphere, depending on hole size $`\mathrm{\Omega }`$ and number of atoms $`N`$. The barrier disappears at $`N>N_{\mathrm{th}}`$, as described in the text.


text:
# Death of Stellar Baryonic Dark Matter Candidates

## 1 Introduction

My basic conclusions of this talk are the following: I. It is looking very likely that 50-90% of our Galaxy is made of nonbaryonic dark matter. II. Regarding the 13-17 microlensing events interpreted as being in our Halo: these are not yet understood.

The nature of the dark matter in the universe and in our Galaxy is one of the great unanswered questions in astrophysics. It is clear from rotation curves of galaxies including our own that most of the matter in galaxies is not in the form of bright stars and instead consists of an unknown component of dark matter. Ten years ago there were two different camps of people on this subject. In the first camp, there were those who believed that the simplest solution would be baryonic dark matter. In particular, the most likely solution appeared to be stellar or substellar objects including faint stars, brown dwarfs, white dwarfs, or neutron stars. In the second camp, there were physicists (particularly motivated by particle physics) who believed there must be a dominant nonbaryonic contribution due to particles such as massive neutrinos, axions, or supersymmetric particles. The main point of my talk is to show that the objects preferred by the first camp, namely stellar baryonic dark matter candidates, are ruled out (see also my conference proceedings in \[confproc1 , confproc2 , confproc3 \] for a longer discussion). Thus nonbaryonic dark matter seems to be favored as explaining the mass of our Galaxy.

In 1986 Hegyi and Olive hegyi  ruled out many obvious candidates for baryons in the Halo of our Galaxy. They ruled out diffuse hot gas, cool neutral hydrogen, small lumps or snowballs of hydrogen, and rocks or dust. In the past decade, microlensing experiments including MACHO, EROS, and OGLE were designed to look for MACHOs, or Massive Compact Halo Objects, which are objects (probably baryonic) in the $`(10^71)M_{}`$ mass range. Instead of resolving the dark matter puzzle, these experiments have raised new issues. The most recent results from the MACHO experiment are discussed by E. Aubourg in this volume. These experiments (macho:2yr , macho:6yr , ansari , laserre ) have the very strong result that they have ruled out a significant component in the Halo of objects in the mass range $`(10^710^2)M_{}`$. These experiments have also found tens of objects not yet understood, with a best fit mass $`0.5M_{}`$, near the mass of a white dwarf. As a consequence, there has been a great deal of recent focus on a possible white dwarf component in the Halo.

In this talk, I will discuss work showing that stellar baryonic candidates for the dark matter are ruled out. The stellar candidates are: 1. Faint Stars. These are objects heavier than about $`0.1M_{}`$ that shine due to hydrogen burning in their cores. 2. Brown Dwarfs. These are objects lighter than about $`0.1M_{}`$ that do not have hydrogen burning in their cores; hence the easiest way to find these objects is to look for them gravitationally (such as with microlensing experiments). 3. White Dwarfs. These are objects with mass $`0.6M_{}`$ and are the remnants of $`(18)M_{}`$ stars. As mentioned above, these are the best fit to the MACHOs found by the microlensing experiments. 4. Neutron Stars. These are the 1.4$`M_{}`$ remnants of stars heavier than $`8M_{}`$.

As I will show in this talk, recent work has shown that none of these four candidates are found in sufficient abundance to explain the mass of our Galactic Halo.

Five years ago, many astronomers believed that the numbers of faint stars and brown dwarfs in the Halo could be quite substantial; in fact these appeared to be the most plausible candidates for the Halo dark matter. It appeared that, as one looked at lower and lower masses, the number of stars seemed to increase such that there could be very many low mass stars and substellar objects. Instead, my work of the last few years with David Graff as well as the work of other authors has shown that the first two candidates, faint stars and brown dwarfs, add up to less than 3% of the mass density of the Galactic Halo. Hubble Space Telescope data found faint stars, and one (bfgk , gf96a , gf96b , mcs , fbg ) can use these data to constrain the mass density of faint stars in the Halo. We showed that faint stars are seen to comprise roughly 1% of the Halo. Brown dwarfs are constrained by both the microlensing experiments (as discussed above) and by our work using a combination of parallax data and theory. We (gf96b ) looked at the faint stars in the parallax data of dahn  and used theory to extrapolate down into the brown dwarf regime. We found that brown dwarfs account for at most a few percent of the Halo. I don’t have time for further discussion of the constraints on these two classes of candidates. I recommend the reader to my previous conference proceeding (confproc3 ) for a fuller discussion.

## 2 White Dwarfs

Next I will proceed to discussion of white dwarfs as dark matter candidates. These are the most interesting stellar candidates currently as they have the masses best fit to the MACHOs seen by microlensing data (if one assumes that these MACHOs are indeed in the Galactic Halo). Is the dark matter made of white dwarfs? There are four problems and issues that need to be addressed: 1. Infrared Radiation 2. Initial Mass Function 3. Baryonic Mass Budget 4. Element Abundances (C, N, He4)

We will see that, for each of these topics, none of the expected signatures of a significant Halo white dwarf population is found.

### 2.1 Constraints from multi-TeV $`\gamma `$-rays seen by HEGRA

The mere existence of multi-TeV $`\gamma `$-rays seen in the HEGRA experiment places a powerful constraint on the allowed abundance of white dwarfs. This arises because the progenitors of the white dwarfs would produce infrared radiation that would prevent the $`\gamma `$-rays from getting here. The $`\gamma `$-rays and infrared photons would interact via $`\gamma \gamma e^+e^{}`$.

Multi-TeV $`\gamma `$-rays from the blazar Mkn 501 at a redshift z=0.034 are seen in the HEGRA detector. The cross section for (1-10)TeV $`\gamma `$-rays peaks at infrared photon energies of (0.03-3)eV. Photons in this energy range would be produced in abundance by the progenitor stars to white dwarfs and neutron stars. By requiring that the optical depth due to $`\gamma \gamma e^+e^{}`$ be less than one for a source at $`z=0.034`$ we limit the cosmological density of stellar remnants (Graff, Freese, Walker, and Pinsonneault gfwp ),

$$\mathrm{\Omega }_{\mathrm{WD}}(13)\times 10^3h^1.$$
(1)
This constraint is quite robust and model independent, as it applies to a variety of models for stellar physics, star formation rate and redshift, mass function, and clustering. In addition, we can be absolutely certain the main sequence progenitors of the white dwarfs produced light!

Note that recent direct observations of infrared light (gorjian , wright ) give comparable constraints on the white dwarf abundance.

### 2.2 Mass Budget Issues

First, I discuss the mass budget issues (based on work by Fields, Freese, and Graff ffg ) general to all Halo Machos, regardless of the type of object.

#### Contribution of Machos to the Mass Density of the Universe:

There is a potential problem in that too many baryons are tied up in Machos and their progenitors. We begin by estimating the contribution of Machos to the mass density of the universe. Microlensing results macho:6yr  predict that the total mass of Machos in the Galactic Halo out to 50 kpc is $`M_{\mathrm{Macho}}=(616)\times 10^{10}M_{}`$ (note that the new numbers are a factor of two lower than the previous estimates of macho:2yr  and in agreement with ansari ). Now one can obtain a “Macho-to-light” ratio for the Halo by dividing by the luminosity of the Milky Way (in the B-band), $`L_{MW}(1.32.5)\times 10^{10}L_{},`$ to obtain $`(M/L)_{\mathrm{Macho}}=(2.613)M_{}/L_{}.`$ ¿From the ESO Slice Project Redshift survey zuc , the luminosity density of the Universe in the $`B`$ band is $`_B=1.9\times 10^8hL_{}\mathrm{Mpc}^3.`$ If we assume that the $`M/L`$ which we defined for the Milky Way is typical of the Universe as a whole, then the universal mass density of Machos is

$$\mathrm{\Omega }_{\mathrm{Macho}}\rho _{\mathrm{Macho}}/\rho _c=(0.0020.01)h^1$$
(2)
where the critical density $`\rho _c3H_0^2/8\pi G=2.71\times 10^{11}h^2M_{}\mathrm{Mpc}^3`$.

We will now proceed to compare our $`\mathrm{\Omega }_{\mathrm{Macho}}`$ derived in Eq. (2) with the baryonic density in the universe, $`\mathrm{\Omega }_\mathrm{B}`$, as determined by primordial nucleosynthesis. To conservatively allow for the full range of possibilities, we will adopt $`\mathrm{\Omega }_\mathrm{B}=(0.0050.026)h^2.`$ Thus, if the Galactic halo Macho interpretation of the microlensing results is correct, Machos make up an important fraction of the baryonic matter of the Universe. Specifically, the central values give

$$\mathrm{\Omega }_{\mathrm{Macho}}/\mathrm{\Omega }_\mathrm{B}0.4h.$$
(3)
However, the lower limit on this fraction is considerably less restrictive,

$$\frac{\mathrm{\Omega }_{\mathrm{Macho}}}{\mathrm{\Omega }_\mathrm{B}}0.1h,$$
(4)
where we have used the lowest $`\mathrm{\Omega }_{\mathrm{Macho}}`$ and the highest $`\mathrm{\Omega }_\mathrm{B}`$.

#### Mass Budget constraints from Machos as Stellar Remnants: White Dwarfs or Neutron Stars

In general, white dwarfs, neutron stars, or black holes all came from significantly heavier progenitors. Hence, the excess mass left over from the progenitors must be added to the calculation of $`\mathrm{\Omega }_{\mathrm{Macho}}`$; the excess mass then leads to stronger constraints. Typically we find the contribution of Macho progenitors to the mass density of the universe to be

$$\mathrm{\Omega }_{\mathrm{prog}}=4\mathrm{\Omega }_{\mathrm{Macho}}=(0.0080.04)h^1.$$
(5)
The central values of all the numbers now imply

$$\mathrm{\Omega }_{\mathrm{prog}}>\mathrm{\Omega }_B,$$
(6)
which is obviously unacceptable. One is driven to the lowest values of $`\mathrm{\Omega }_{\mathrm{Macho}}`$ and highest value of $`\mathrm{\Omega }_B`$ to avoid this problem.

### 2.3 On Carbon and Nitrogen

The overproduction of carbon and/or nitrogen produced by white dwarf progenitors is one of the greatest difficulties faced by a white dwarf dark matter scenario, as first noted by Smecker and Wyse (wyse ) and Gibson and Mould (gm ). Stellar carbon yields for zero metallicity stars are quite uncertain. Still, according to the Van den Hoek & Groenewegen (1997) yields, a star of mass 2.5$`M_{}`$ will produce about twice the solar enrichment of carbon. However, stars in our galactic halo have carbon abundance in the range $`10^410^2`$ solar. Hence the ejecta of a large population of white dwarfs would have to be removed from the galaxy via a galactic wind.

However, carbon abundances in intermediate redshift Ly$`\alpha `$ forest lines have recently been measured to be quite low, at the $`10^2`$ solar level sc , for Ly$`\alpha `$ systems at $`z3`$ with column densities $`N3\times 10^{15}\mathrm{cm}^2`$ (for lower column densities, the mean C/H drops to $`10^{3.5}`$ solar lsbr .

In order to maintain carbon abundances as low as $`10^2`$ solar, only about $`10^2`$ of all baryons can have passed through the intermediate mass stars that were the predecessors of Machos (Fields, Freese, and Graff ffg ). Such a fraction can barely be accommodated by our results in section 4.1 for the remnant density predicted from our extrapolation of the Macho group results, and would be in conflict with $`\mathrm{\Omega }_{\mathrm{prog}}`$ in the case of a single burst of star formation. Note that stars heavier than 4$`M_{}`$ may replace the carbon overproduction problem with nitrogen overproduction vdhg  lattbooth ).

Using the yields described above, we calculated the C and N that would result from the stellar processing for a variety of initial mass functions for the white dwarf progenitors. We used a chemical evolution model based on a code described in Fields & Olive fo98  to obtain our numerical results. Our results are presented in the figure.

In the figure, we make the parameter choices that are in agreement with D and He4 measurements (see the discussion below) and are the least restrictive when comparing with the Ly$`\alpha `$ measurements. We take an initial mass function (IMF) sharply peaked at 2$`M_{}`$, so that there are very few progenitor stars heavier than 3$`M_{}`$ (this IMF is required by D and He4 measurements). In addition (see the figures in Fields, Freese, and Graff ffg2 ) we have considered a variety of other parameter choices. By comparing with the observations, we obtain the limit,

$$\mathrm{\Omega }_{\mathrm{WD}}h2\times 10^4.$$
(7)
As a caveat, note that it is possible that carbon never leaves the (zero metallicity) white dwarf progenitors, so that carbon overproduction is not a problem chabriernew .

### 2.4 Deuterium and Helium

Because of the uncertainty in the C and N yields from low-metallicity stars, we have also calculated the D and He4 abundances that would be produced by white dwarf progenitors. These are far less uncertain as they are produced farther out from the center of the star and do not have to be dredged up from the core. Panel a) in the figure displays our results. Also shown are the initial values from big bang nucleosynthesis and the (very generous) range of primordial values of D and He4 from observations. ¿From D and He alone, we can see that the white dwarf progenitor IMF must be peaked at low masses, $`2M_{}`$. We obtain

$$\mathrm{\Omega }_{\mathrm{WD}}0.003.$$
(8)
### 2.5 Is Dark Matter Made of White Dwarfs?

To reiterate, there are four major problems with a white dwarf Halo: 1) infrared radiation, 2) initial mass function, 3) baryonic mass budget, and 4) element abundances (C, N, He4). We have found that of the expected signature of a white dwarf population is found. Hence white dwarfs cannot explain the full dark matter of the Halo.

A second question remains: can the Macho data be explained by white dwarfs? If one compares $`\mathrm{\Omega }_{\mathrm{Macho}}`$ with the HEGRA limit presented above in eq.(1), all MACHOs can still be white dwarfs if one takes the extreme numbers. However, according to the limit in eq.(7), NOT all Machos can be made of white dwarfs, even with the most extreme numbers. Gates and Gyuk (gates ) have proposed the following explanation of the Macho data: white dwarfs in an enlarged protodisk can explain the Macho events with $`M=7\times 10^{10}M_{}`$ and $`\mathrm{\Omega }_{\mathrm{Macho}}3\times 103`$. These white dwarfs would comprise (3-4)% of the Halo density. Even such a small Halo fraction is hard to reconcile with eq.(7) above.

Recent observations (ibata1 ) (ibata2 ) have found evidence of direct optical detections of objects that may be Halo white dwarfs. The situation regarding these objects is unclear. The objects found in the Hubble Deep Field (ibata1 ) are in conflict with what was found earlier in the Luyten survey; if the new observations are correct, large numbers of these objects should have been found in the Luyten survey (flynn ). Regarding the two white dwarfs found in (ibata2 ), if one takes into account the Poisson statistics, these two white dwarf are consistent with white dwarfs from known stellar populations (graff3 ).

In any case the bulk of the dark matter has not yet been found.

## 3 Zero Macho Halo?

The possibility exists that the microlensing events that have been interpreted as being in the Halo of the Galaxy are in fact due to some other lensing population. One of the most difficult aspects of microlensing is the degeneracy of the interpretation of the data, so that it is currently impossible to determine whether the lenses lie in the Galactic Halo, or in the Disk of the Milky Way, or in the LMC. In particular, it is possible that the LMC is thicker than previously thought so that the observed events are due to self-lensing of the LMC. All these possibilities are being investigated. More data are required in order to identify where the lenses are.

## 4 Conclusions

Microlensing experiments have ruled out a large class of possible baryonic dark matter components. Microlensing experiments have ruled out objects in the mass range $`10^7M_{}`$ all the way up to $`10^2M_{}`$. In this talk I discussed the heavier possibilities in the range $`10^2M_{}`$ to a few $`M_{}`$. Brown dwarfs and faint stars are ruled out as significant dark matter components; they contribute no more than 1% of the Halo mass density. Stellar remnants are not able to explain the dark matter of the Galaxy either; none of the expected signatures of stellar remnants, i.e., infrared radiation, large baryonic mass budget, and C,N, and He4 abundances, are found observationally.

Hence, in conclusion, 1) Nonbaryonic dark matter in our Galaxy seems to be required, and 2) The nature of the Machos seen in microlensing experiments and interpreted as the dark matter in the Halo of our Galaxy remains a mystery. Are we driven to primordial black holes carr  jedam , nonbaryonic Machos (Machismos?), mirror matter Machos (mohap ) or perhaps a no-Macho Halo?

## 5 References