text:
# A CLEAN-based Method for Deconvolving Interstellar Pulse Broadening from Radio Pulses

## 1 Introduction

Multipath propagation in the interstellar medium (ISM) at radio wavelengths is caused by diffraction from irregularities in electron density on scales smaller than the Fresnel scale, $`\sqrt{\lambda D}10^{11}`$ cm (where $`\lambda `$ is the wavelength of observation and $`D`$ is the extent of the scattering medium). One consequence is that pulses from radio pulsars are distorted by the differences in arrival time of different portions of the scattered wavefronts. The resultant pulse broadening (e.g. Williamson, 1972) is usually quantified in terms of the pulse-broadening time, $`\tau _d`$, which is related to the scattering measure SM, where SM is the integrated value of $`C_n^2`$, the coefficient of the wavenumber spectrum of electron density irregularities. Pulse broadening is important information for quantifying the Galactic structure of electron density fluctuations. But it is also a hindrance in studying the intrinsic properties of pulsars that are manifested in the undistorted pulse shape.

For a Kolmogorov wavenumber spectrum with small inner scale (e.g. Rickett, 1990; Cordes & Lazio, 1991; Armstrong, Rickett & Spangler, 1995), the pulse broadening time scales as

$`\tau _d`$ $``$ $`1.10\mathrm{ms}W_\tau \mathrm{SM}^{6/5}\nu ^{22/5}D,`$ (1)
where $`\nu `$ is the frequency of observation (GHz), $`D`$ is the pulsar distance (kpc), SM has units kpc m-20/3, and $`W_\tau `$ is a geometric factor that accounts for the distribution of $`C_n^2`$ along the line of sight (e.g. Cordes & Rickett, 1998). Pulsars with larger dispersion measures ($`\mathrm{DM}`$) tend to have larger scattering measures. The effect is thus predominant in observations of distant pulsars and at lower frequencies, and it limits the utility of such scattered pulse profiles for the study of pulsar magnetospheres.

The observed pulse $`y(t)`$ can be modeled as a convolution of the intrinsic pulse $`x(t)`$ with the pulse-broadening function (PBF) $`g(t)`$ and the net instrumental resolution function, $`r(t)`$:

$`y(t)=x(t)g(t)r(t);`$ (2)
the PBF $`g(t)`$ is the response of the interstellar medium (ISM) to a delta function. The exact form of the PBF depends on the spatial distribution of scattering material along the line of sight and on its wavenumber spectrum (e.g. Williamson, 1972, 1973; Cordes & Rickett, 1998; Lambert & Rickett, 1999; Cordes & Lazio, 2001; Boldyrev & Gwinn, 2002).

In general, reconstruction of the original pulse shape from the scattered data is complicated by the fact that the PBF $`g(t)`$ is unknown. However, classes of possible PBFs corresponding to different distributions of scattering media are known.

Independent determination of the PBF requires knowledge of the wavenumber spectrum and distribution of scattering material along the line of sight. Alternatively, a priori knowledge of the intrinsic pulse shape could allow determination of the PBF. In general, however, both $`x(t)`$ and $`g(t)`$ are unknown; hence the utility of scattered pulse profiles is limited to measuring $`\tau _d`$ through a template fitting method, which assumes a form for the intrinsic pulse shape and a model for the ISM, rather than determining these quantities from the data. The assumed shapes can be incorrect, leading to estimates for $`\tau _d`$ which are only approximate, and have poorly quantified uncertainties.

In this paper, we describe a method for reconstruction of the intrinsic pulse shape as well as determination of the PBF. The basic idea is derived from the CLEAN algorithm extensively used in image reconstruction of synthesis imaging data (Högbom, 1974; Schwarz, 1978). Our one-dimensional application is similar to the usage of CLEAN in spectral analysis (Roberts, Lehar & Dreher, 1987). In §2, we outline the method and its practical realization, and compare it to other methods that have been used for this application. This CLEAN-based method is demonstrated using simulated data as well as some recent data from Arecibo observations, and the results are discussed in §3. Finally, we present our conclusions in §4.

## 2 Method and Practical Realization

Before outlining the CLEAN-based approach, we briefly describe other methods that have been used to estimate the pulse-broadening time from scattered pulse profiles. In general, a model is fitted to the observed pulse profile, $`y_o\text{ }(t)`$, where the model, $`y_m\text{ }(t)`$, is the convolution of a template profile, representing an estimate for the intrinsic pulse shape, with an impulse response, $`g(t)`$, that describes pulse broadening from the ISM for the particular line of sight.

### 2.1 Conventional Frequency-Extrapolation Approach

Typically, the intrinsic pulse shape at the frequency of interest is extrapolated from high-frequency observations where the pulse broadening is assumed to be negligible. In addition, the PBF shape is also assumed known, often as a one-sided exponential function,

$`g(t)=\tau _{d}^{}{}_{}{}^{1}\mathrm{exp}(t/\tau _d)U(t),`$ (3)
where $`U(t)`$ is the unit step function. As discussed below, this form applies to a thin scattering screen for which the density irregularities follow a square-law structure function. The pulse broadening time $`\tau _d`$ is estimated by minimizing $`\chi ^2`$ (as standardly defined) as a function of $`\tau _d`$, the sole parameter in this problem. Results from the method thus rely on use of an empirical pulse shape derived from radio frequencies quite different from that at which the scattering is prominent. However, it is well known that profile evolution with frequency can be quite complex for some pulsars (e.g. Hankins & Rickett, 1986), and extrapolation to the frequency of interest will be in error. It is also well known that the PBF is a one-sided exponential only for a medium described as a thin screen with a square-law structure function (e.g. Cordes & Rickett, 1998). While thin screens may exist along some lines of sight, a more reasonable default model is one where the medium is Kolmogorov (leading to a non-square law structure function) and fills a significant fraction of the line of sight. Despite such shortcomings, the method has been used to obtain estimates of $`\tau _d`$ for a large number of pulsars (e.g. Löhmer et al., 2001).

Variations on this approach include use of non-exponential pulse broadening functions and detailed estimation of the frequency dependence of the intrinsic pulse shape, for example by modeling individual components of the pulse shape and scaling each with frequency. The essence of the method is the same, however, in that measurements at other frequencies are used to estimate the intrinsic pulse shape at the frequency of interest.

### 2.2 Fourier Inversion

From Eq. 2, the Fourier transform of $`y(t)`$ is related to that of $`x(t)`$ as $`Y(f)=X(f)G(f)R(f)`$, where $`X,G`$ and $`R`$ are the respective Fourier transforms of the factors on the right-hand side of Eq. 2. Adopting a particular form for $`g(t)`$ and knowing $`r(t)`$, one can, in principle, perform the deconvolution by calculating $`X(f)=Y(f)/G(f)R(f)`$ and inverse-transforming to obtain the intrinsic pulse shape. This method was used by Weisberg et al. (1990) and Kuz’min & Izvekova (1993) with a one-sided exponential PBF. Kuz’min & Izvekova pointed out that the results depend on the assumed form of the PBF but did not present any results using alternative forms.

### 2.3 A CLEAN-based approach

In this paper, we propose an alternative method to realize both objectives of recovering the intrinsic pulse shape and determining the shape and characteristic time scale of the broadening function. The method is based on the CLEAN algorithm, an iterative method for deconvolving the instrumental point source response function from an interferometric image (Högbom, 1974; Schwarz, 1978). In the present case, the problem is much simplified as the data are one-dimensional. The method involves inverting the convolution in Eq. 2. By analogy with imaging applications, the measured, scattered pulse is equivalent to the “dirty map”; and the combination $`g(t)r(t)`$ (the PBF convolved with the instrumental response) is equivalent to the “dirty beam”. In synthesis imaging applications, CLEANing involves the iterative subtraction of scaled copies of the dirty beam from the dirty map at locations corresponding to real features until the residual map is indistinguishable from noise, followed by restoration of the subtracted CLEAN components. Here, we adapt the method to (one-dimensional) scattered pulse profiles.

While this method is inspired by CLEAN as used in synthesis imaging, we emphasize that there are some important differences between the two algorithms. Specifically, in imaging applications, the dirty beam is assumed to be perfectly known, and sets an upper limit on the achievable resolution. In the algorithm outlined here, we are performing a true deconvolution, and the resolution of the deconvolved pulse is determined by the response function of the instrument, rather than the width of the PBF. Additionally, the PBF is not known a priori, and we perform a search in both parameter space and function space for the best pulse-broadening function, as outlined below.

We consider the measured pulse to consist of delta functions, each convolved with the PBF, $`g(t)`$, and the resolution function, $`r(t)`$. In reality, radiometer noise adds to the measured quantity and the intrinsic pulse shape is itself perturbed by pulse-to-pulse fluctuations intrinsic to the pulsar. The CLEAN process identifies the amplitude and location of a CLEAN component (CC) by finding a suitable maximum in the residual profile in each iteration. Initially, the residual profile is the measured profile, $`y(t)`$. The amplitude of the CC is the identified maximum multiplied by a loop gain, $`\gamma `$. The CC $`y_c`$(t) is convolved with $`g(t)`$ and $`r(t)`$ and subtracted from the measured pulse profile to yield the residual pulse profile, $`\mathrm{\Delta }y(t)`$:

$$\mathrm{\Delta }y(t_i\text{ })=y(t_i\text{ })\left[y_c\text{ }(t)\left(g(t)r(t)\right)\right]_{t=t_i},i=1,\mathrm{},N,$$
(4)
where $`y_c\text{ }(t_i\text{ })=\gamma \text{ }\{\mathrm{max}[y(t)]\}\delta (tt_o\text{ })`$, and $`N`$ is the number of pulse phase bins. $`\mathrm{\Delta }y(t)`$ for a given iteration becomes $`y(t)`$ for the next iteration, and the process is repeated until the residual subtracted pulse, $`\mathrm{\Delta }y(t)`$, is comparable to the off-pulse rms noise level111The termination criterion is not well-defined for the CLEAN technique.. Upon termination of the iteration, $`n_c`$ CLEAN components are identified: $`C_j,t_j`$, $`j=1,\mathrm{},n_c`$ denote the amplitudes and times of these CLEAN components (i.e., the collection of delta functions at the end of the CLEAN process).

The restored pulse shape is built from the ensemble of CLEAN components, $`C_j`$, convolved with a ‘restoring function’ $`\rho (t)`$, details of which we will discuss later in this section:

$$y_r(t)=\underset{j=1}{\overset{n_c\text{ }}{}}C_j\delta (tt_j)\rho (t).$$
(5)
At the end of the CLEAN process the residual noise, $`\mathrm{\Delta }y(t)`$, is added to the convolved CCs to obtain the final restored pulse profile.

Examples of the CLEAN process are shown in Figure 1 for simulated data (left-hand panels) and real data (right-hand panels). The figure shows the measured pulse shapes, the chosen PBF (a one-sided exponential as in Eq. 3), the CCs, and the deconvolved profile. Compared to the traditional least-squares or frequency-extrapolation methods (§2.1), the approach to the pulse-broadening problem outlined here is much easier to apply because it is non-parameteric. Further, it does not make any assumption about the intrinsic pulse shape or its evolution with frequency, and thus can potentially yield useful information about the intrinsic pulse shape.

Below we elaborate on some aspects concerning the practical realization of the method.

#### 2.3.1 Loop Gain

Small values of the loop gain ($`\gamma `$) are employed in imaging applications. The performance of the method is expected to change only marginally for sufficiently small values (say, $`\gamma \text{ }0.1`$). Larger values tend to cause over-subtraction and consequent lack of convergence. We adopted $`\gamma `$=0.05 as a reasonable choice for our application. Thus, in each iteration, a delta function with an amplitude of 5% of the profile peak is convolved with the chosen PBF and the instrumental response function $`r(t)`$, and then subtracted from the peak of the residual pulse profile.

#### 2.3.2  Pulse Broadening Functions

Unlike radio synthesis imaging applications, where the dirty beam is often known exactly as the response of the instrument to a point source, the exact functional form for the PBF of the ISM is not known. It depends on both the distribution and the wavenumber spectrum of scattering material along the line-of-sight to the pulsar, neither of which is known in detail. The PBF is particularly sensitive to the spatial distribution of $`C_n^2`$. For instance, in the simplest case of a thin slab of scattering material located in between the pulsar and the Earth, the PBF is characterized by a sharp rise and exponential fall-off (eq. 3), whereas the PBFs of more realistic models (such as a thick slab or a uniform medium) are characterized by rounded edges and slower fall-offs, with different rise and peak times. For media with a square-law structure function (Cordes & Rickett, 1998; Lambert & Rickett, 1999) the PBF for a thin screen is given by Eq. 3 and for thick and uniform media the PBFs are given by Williamson (1972, 1973):

$`g(t)`$ $`=`$ $`\left(\pi \tau _d/4t^3\right)^{1/2}\mathrm{exp}(\pi ^2\tau _d/16t)`$ (6)
$`g(t)`$ $`=`$ $`\left(\pi ^5\tau _{d}^{}{}_{}{}^{3}/8t^5\right)^{1/2}\mathrm{exp}(\pi ^2\tau _d/4t)`$ (7)
where the time constant $`\tau _d`$ is the pulse-broadening time.

For other media, such as those with a Kolmogorov wavenumber spectrum having different inner scales (Lambert & Rickett, 1999) or media with non-Gaussian density fluctuations (Boldyrev & Gwinn, 2002) (but still having uniform statistics transverse to the line of sight), the PBFs are qualitatively similar to those for media with square-law structure functions in that thin screens show sharp rise times while extended media have slower rise times.

Nonuniformities of the scattering transverse to the line of sight, particularly any truncation of the scattering screen, produce PBFs truncated for times $`t>t_{\mathrm{max}}\text{ }`$, where $`t_{\mathrm{max}}`$ is the arrival time corresponding to the edge of the screen (Cordes & Lazio, 2001). As an example, a circular disk centered on the line of sight to a source will produce a PBF of the form

$$g(t)=\frac{\mathrm{exp}(t/\tau _d)}{\tau _d(1\mathrm{exp}(\zeta ))}U(t_{\mathrm{max}}\text{ }t),$$
(8)
where $`t_{\mathrm{max}}\text{ }=\zeta \tau _d`$, $`\zeta =\theta _{\mathrm{max}}^2\text{ }/2\sigma _\theta ^2\text{ }`$, where $`\theta _{\mathrm{max}}`$ and $`\sigma _\theta `$ are the maximum and RMS scattering angles, respectively. Cordes & Lazio (2001) also discuss PBFs for several other configurations. For example, in the case of scattering from a long, narrow filament, $`\theta _{\mathrm{max}}`$ will be much larger than $`\sigma _\theta `$, and consequently the scattered image becomes elongated, and the PBF tends toward the form $`t^{1/2}\mathrm{exp}(t)`$. For a large screen with elongated irregularities, the PBFs will be strongly frequency dependent, and will be characterized by multiple peaks that align at different frequencies.

In applying our algorithm, one can perform a grid search over various forms of the PBF as well as the relevant parameters and thus determine (a) the appropriate model of the PBF that best describe observations, and (b) the best fit $`\tau _d`$ that yields meaningful intrinsic pulse shapes. The figures of merit for determining the best fit case is discussed later in this section.

#### 2.3.3  Instrumental Response

The response function due to instrument can be described as the convolution of the functions that characterize various instrumental effects relevant in determining the effective pulse phase resolution, including (i) dispersion smearing ($`r_{dm}`$), (ii) profile binning ($`r_{pb}`$), (iii) the backend time averaging ($`r_{av}`$), and (iv) additional post-detection time averaging ($`r_{pd}`$). The response function of the instrument, or the effective resolution function of the data, is therefore given by

$$r(t)=r_{dm}(t)r_{pb}(t)r_{av}(t)r_{pd}(t).$$
(9)
The contributing factors in the above equation, in particular those which describe effects (i) to (iii), can be treated as approximately rectangular. Hence the response function will be approximately “trapezoidal” for cases where one or two factor(s) dominate the effective time resolution of the data. It will be more “Gaussian-like” due to the convolution of similar-sized rectangle functions when there are multiple contributions that determine the time resolution. We form the resolution function according to eq. 9, where the convolution is performed with a time resolution that is much larger than the narrowest factor; the resultant function is then re-sampled to the resolution of the actual data.

#### 2.3.4 Restoring function

The CLEANed pulse shape is obtained by convolving CCs with a restoring function, $`\rho (t)`$, equivalent to the CLEAN beam in imaging applications. The choice of a restoring beam is driven by several ideas (see Cornwell, Braun & Briggs (1999) for a detailed discussion). Following the analogy of imaging applications, where the CLEAN beam is approximated as a Gaussian fitted to the central part of the dirty beam, we adopt an equivalent Gaussian, $`r_g(t)`$, whose width reflects the effective resolution of the actual data (given by the width of $`r(t)`$ as discussed in the previous section). The choice of restoring function, $`\rho (t)=r_g(t)`$, de-emphasizes the higher frequency components that are sometimes spuriously generated by the CLEAN algorithm when directly using the response function $`r(t)`$ as the restoring function. The restoring function is normalized to unit area in order to ensure conservation of flux in the inversion process.

#### 2.3.5 Figures of Merit

We now discuss the figures of merit for determining the best fit case, i.e., when the observed pulse profile is optimally compensated for interstellar broadening. Assume that we know the shape of the PBF but not its characteristic scale, $`\tau _d`$. If we choose $`\tau _d`$ too large, we will tend to oversubtract pulsed flux and thus introduce negative going features in the deconvolved profile. Conversely, if too-small a value for $`\tau _d`$ is used, the deconvolved profile will display residual asymmetry due to un-deconvolved scattering. Therefore two figures of merit that are appropriate are: (1) positivity and (2) minimum asymmetry of the resultant pulse. In practice we numerically require the residual noise (i.e., measured$``$modeled pulse) to show zero mean and the same rms everywhere. We use all these criteria to choose the correct PBF and characteristic time scale, $`\tau _d`$.

We define a parameter ($`f_r`$) that measures positivity as

$$f_r=\frac{m}{N\sigma _{\mathrm{off}}^2\text{ }}\underset{i=1}{\overset{N}{}}[\mathrm{\Delta }y(t_i\text{ })]^2U_{\mathrm{\Delta }y}$$
(10)
where

$`U_{\mathrm{\Delta }y}=U\left(\mathrm{\Delta }y(t_i\text{ })x\sigma _{\mathrm{off}}\text{ }\right),`$ (11)
$`\mathrm{\Delta }y(t_i\text{ })`$ is the residual noise at the end of the CLEAN process, and $`m`$ denotes a ‘weight’ of order unity. The unit step function $`U`$ turns on when the residual $`\mathrm{\Delta }y`$ is significantly below the off-pulse noise level (i.e. when oversubraction is caused by using a PBF that is too broad). We adopted $`x=3/2`$ for our analysis so that distortions of the pulse more negative than $`3\sigma /2`$ yield a penalty manifested as a larger value of $`f_r`$. To minimize profile asymmetry, we minimize the skewness parameter of the CLEANed profile, $`y_r(t)`$, by explicitly minimizing the skewness ($`\mathrm{\Gamma }`$) of the CCs:

$`\mathrm{\Gamma }={\displaystyle \frac{t^3}{t^2^{3/2}}},`$ (12)
where

$`t^n`$ $`=`$ $`{\displaystyle \frac{{\displaystyle \underset{i=1}{\overset{n_c}{}}}(t_i\overline{t})^nC_i}{{\displaystyle \underset{i=1}{\overset{n_c}{}}}C_i}},`$ (13)
$`\overline{t}`$ $`=`$ $`{\displaystyle \frac{{\displaystyle \underset{i=1}{\overset{n_c}{}}}t_iC_i}{{\displaystyle \underset{i=1}{\overset{n_c}{}}}C_i}},`$ (14)
and $`C_i,t_i`$, $`i=1,\mathrm{},n_c`$ denote the amplitudes and times of the CLEAN components, as defined before.

In addition to the symmetry and positivity constraints, we also maximize the number of residual points in the on-pulse window that were consistent with the noise in an off-pulse window. Specifically we count the number of points $`N_f`$ that satisfy $`|y_i\text{ }y_{\mathrm{off}}\text{ }|3\sigma _{\mathrm{off}}\text{ }`$, where $`y_{\mathrm{off}}`$and $`\sigma _{\mathrm{off}}`$represent the mean and rms of the off-pulse region. Thus the figures of merit we employ include (a) maximum $`N_f`$, (b) minimum residual noise, (c) minimal skewness, and (d) minimum positivity parameter $`f_r`$.

Sample plots illustrating the figures of merit are shown in Figure 2 for the simulated and actual data of Figure 1. The best fit model will correspond to minimal values of both $`f_r`$ and $`\mathrm{\Gamma }`$, or equivalently to a minimum of the combined parameter, $`f_c=(\mathrm{\Gamma }+f_r)/2`$, as well as conforming to a maximum of $`N_f`$ and a minimum of the residual noise rms ($`\sigma _{\mathrm{offc}}`$) of the CLEANed pulse. Economy in the number of CCs required to fit the observed pulse profile may also provide a means to discriminate between different forms of the PBF (and hence, different distributions of scattering material).

As the deconvolution process removes smearing due to interstellar scattering, the resultant CLEANed pulse is narrower and the peak signal to noise (S/N) is higher, the increase being determined largely by the ratio of $`\tau _d`$ to pulse width. This is evident from examples of simulated and real data shown in Figs. 1 and 36; examination of Figure 1 suggests that the sum of CCs matches closely the area of the scattered pulse.

## 3 Further Examples

We illustrate the application of the above-described method using additional simulated and real data. In doing so, we show the sensitivity of the method for finding both the shape of the pulse broadening function and the appropriate $`\tau _d`$. We also demonstrate how the derived intrinsic pulse shapes vary with different PBFs.

### 3.1 Simulated Data

Pulse profiles are generated as single or multi-component Gaussians. Scattered profiles are then constructed for a given choice of PBF, and noise is added to the scattered profile to account for finite signal-to-noise ratios. The simulated scattered data are then input to our CLEAN algorithm to recover the intrinsic pulse shape as well as the best fit model PBF, and the results are compared with the original data and PBF used to generate the data. Figure 3 illustrates the case where the intrinsic pulse is a simple Gaussian and the ISM is described by a simple thin slab model. As evident from the figure, the intrinsic pulse obtained by application of the CLEAN algorithm is in very good agreement with the original data. A significant deviation from the best fit model PBF (a one-sided exponential with a characteristic timescale of 40 ms) leads to distortion of the resultant deconvolved pulse profile, leading to residual scattering (for PBFs narrower than the optimal case) or negative going features at trailing edge of the pulse (for PBFs wider than optimal case). This example thus demonstrates the ability of the technique to well characterize the PBF and recover the intrinsic pulse.

Fig 4 shows another example where the intrinsic pulse is modeled as a sum of Gaussians, to mimic a pulse profile that is composed of multiple components. As in the previous example, the ISM model assumed here is a thin slab of scattering material, which can be described by an exponential PBF with a characteristic time scale $`\tau _d=60`$ ms. The recovered pulse is in good agreement with the original data.

### 3.2 Data from Observations

Having demonstrated the feasibility of the technique using simulated data, we now discuss its application to some real data. The data used in this paper are from recent observations at Arecibo, where a large sample of high-DM pulsars has been studied at multiple frequency bands (Bhat et al., 2002). The sample includes several new discoveries from the Parkes multibeam survey (Manchester et al., 2001) that are visible at Arecibo. Details of the observations and results will be presented in a forthcoming paper. Figure 5 shows the 1475 MHz profile of PSR J1852+0031 (B1849+00), whose DM$`=680`$ $`\mathrm{pc}\mathrm{cm}^3`$ implies a distance D$`=8.4`$ kpc using the model of Taylor & Cordes (1993). We attempted deconvolving these data with two widely different types of PBFs: (a) a sharp-edged PBF (eq. 3) that characterizes a thin slab model, (b) a rounded PBF (eq. 7) for a uniform distribution of scattering material. As the figure shows, the results are strikingly different; the final CLEAN pulse appears approximately triple for an exponential PBF (with a characteristic time $`\tau _d`$=225$`\pm `$14 ms), whereas it is a distinct double for the case of a uniform medium ($`\tau _d`$=121$`\pm `$6 ms). At first glance, either of the models seem equally good. However, given that two-component profiles are more commonly observed at and near this radio frequency (which implies observations are more sensitive to a conal geometry than a core and cone case), we believe that a double profile is more likely. The double profile obtained with a rounded PBF is constituted from only 7 CCs, compared to 21 required to fit the data with an exponential PBF. As alluded to before, this economy in the number of CCs lends weight to the idea that the intrinsic profile is double. Very similar results are obtained for the same pulsar at 1175 MHz (Figure 6a), which implies that although the folded profile at 1475 MHz is an average from only $``$120 pulse periods, the result is not likely to be affected by inadequate profile stability.

More examples are shown in Figure 6, where the data span a wide range of signal to noise ratio and degree of scattering. The CLEAN method yields satisfactory and robust results even for cases where the scattering tail is not very prominent, as illustrated by the scattering tail of PSR J1902+0556, which extends only to a few percent of the pulse period. The 1175 MHz profile of PSR J1852+0031 illustrates the other extreme, where the scattering tail extends out to nearly half the pulse period. For the 430 MHz data on PSR J1901+0331, three different types of PBFs (eqs. 3, 7 and 7) were used for deconvolution. As evident from Figure 6c, the observations are better accounted for by a thin slab case with exponential PBF of $`\tau _d`$=60$`\pm `$3 ms or by a uniform medium with a PBF of characteristic time $`\tau _d`$=31$`\pm `$2 ms, rather than an intermediate case such as a thick slab (eq. 7). The residual asymmetry in the deconvolved pulse is most likely intrinsic to the pulsar beam, though it is possible that the true PBF for this line of sight differs from what we have considered. Detailed interpretations are deferred to a forthcoming paper where we will consider multifrequency data and additional forms for the PBF.

### 3.3 Estimates of Pulse-broadening Time

The CLEAN-based approach we outlined in this paper provides a general tool that can be applied to a wide variety of measured pulse shapes and degrees of scattering. For pulse shapes that show significant scattering tails, application of the frequency extrapolation approach as outlined in §2.1 can be expected to yield satisfactory results. We have applied the CLEAN and frequency-extrapolation methods to strongly scattered pulses, such as those in Figure 6, and the results are tabulated in Table 1. The data span a wide range of DM, signal to noise ratio and degree of scattering. In the frequency extrapolation method, the uncertainty in $`\tau _d`$ is estimated as the model parameter range where $`\chi ^2`$is unity above the minimum. For the CLEAN method, a similar estimate can be defined based on the figure of merit parameter $`f_r`$, which is essentially an indicator of the average deviation of the residuals from the noise level (normalized to the noise variance). For instance, if a model results in a value for $`f_r`$ that is unity above that for the best fit PBF, then the pulse is over-subtracted by $`\sigma _{\mathrm{off}}`$(on average). Note that the best fit model is characterized by a minimum in the combined figure of merit parameter $`f_c`$.

The CLEAN-based pulse-broadening measurements in Table 1 for PSRs J1852+0031 and J1853+0546 indicate scaling laws $`\tau _d\nu ^x`$ with $`x=3.45\pm 0.7`$ and $`x=3.4\pm 1.0`$, respectively. These scalings may be compared with those expected from a Kolmogorov wavenumber spectrum with a negligibly small inner scale ($`x=4.4`$) or from a medium with a square-law structure function ($`x=4`$). The empirical intervals for $`x`$ include both cases, though the centroids of the intervals are certainly biased below the expected values. We defer a detailed analysis to a later paper, but we point out here that the inferred scaling index for a given pulsar will depend on the PBF used and that the PBF shape may also be a function of frequency. Ideally, the data and the CLEAN procedure will indicate the best PBF through investigation of the figures of merit discussed above. Both the PBF and the scaling index of $`\tau _d`$ are useful information for understanding the nature and distribution of scattering along the line of sight (Cordes & Lazio, 2001; Löhmer et al., 2001; Boldyrev & Gwinn, 2002). For multifrequency data like these, the CLEAN-based approach we have outlined here can be extended to do a global fit (by assuming the same form for the PBF, and $`\tau _d`$ scaling with the frequency as a power law but with an unknown index) to yield the pulse-broadening times and the index of the power law, along with the intrinisic pulse shapes. A detailed treatment of this approach is deferred to a forthcoming paper where we report on multifrequency pulse-broadening measurements for a larger sample of distant pulsars.

## 4 Conclusions

We have developed and applied an algorithm for deconvolving interstellar pulse broadening from radio pulses. The technique realizes the two-fold objective of recovering the intrinsic pulse shape and determining the pulse-broadening function that describes the scattering process. The method parallels the CLEAN algorithm used in synthesis imaging applications, but differs from CLEAN in that the best pulse broadening function is not known in advance. Additionally, we perform a true deconvolution, where the resolution of the deconvolved pulse is limited by the instrumental resolution rather than the pulse broadening function.

Application of the technique is demonstrated through use of simulated data as well as data from Arecibo observations. Unlike the conventional frequency extrapolation approach, the method makes no prior assumptions about the intrinsic pulse shape, and can therefore be applied to a wider variety of pulse shapes and degrees of scattering.

The derived intrinsic pulse shape depends on the form and parameters of the PBF used for deconvolving, and further considerations have to be employed to select the best fit model. We employ figures of merit that quantify positivity and asymmetry of the CLEANed pulse, and the best fit PBF is taken to be the one that yields minimal residual scattering and skewness of the intrinsic pulse. As in any other nonlinear method, the uniqueness of the deconvolved result is not always assured: however, further considerations (such as data at different observing frequencies) can help resolve any such degeneracy. Intrinsic pulse shapes obtained by application of this method are likely to be closer to reality, and hence should prove to be useful for the study of pulsar emission properties. A paper is in preparation which applies our method to a large body of multi-frequency pulsar observations.

We thank the referee, S. R. Spangler, for several thoughtful comments that helped us to improve the paper. This work was supported by NSF Grant AST 9819931 and by the National Astronomy and Ionosphere Center, which operates the Arecibo Observatory under cooperative agreement with the NSF.


format: text

text:
# Formation of atomic tritium clusters and condensates

## Abstract

We present an extensive study of the static and dynamic properties of systems of spin-polarized tritium atoms. In particular, we calculate the two-body $`|F,m_F=|0,0`$ $`s`$-wave scattering length and show that it can be manipulated via a Feshbach resonance at a field strength of about 870 G. Such a resonance might be exploited to make and control a Bose-Einstein condensate of tritium in the $`|0,0`$ state. It is further shown that the quartet tritium trimer is the only bound hydrogen isotope and that its single vibrational bound state is a Borromean state. The ground state properties of larger spin-polarized tritium clusters are also presented and compared with those of helium clusters.

05.30.Jp,34.50,36.40.-c

In 1976 , Stwalley and Nosanow suggested, based on statistical arguments, that the system of spin-polarized bosonic tritium atoms behaves “very much like 4He”. To the best of our knowledge, their arguments have not yet been tested by a microscopic quantum mechanical treatment. A detailed theoretical study of spin-polarized tritium systems, namely spin-polarized atomic tritium clusters and optically-pumped tritium condensates, is the objective of the present work. In this Letter, we present results for both structural and scattering properties of tritium dimers, trimers, and clusters.

Pioneering experimental studies of the lowest quartet state of spin-polarized atomic trimers have been pursued recently for sodium and potassium . Study of these trimers, which were prepared on the surface of large 4He clusters, revealed that three-body effects are surprisingly important . We are not aware, though, of any experimental or theoretical studies of larger spin-polarized atomic cluster systems. Bosonic helium systems — i.e. liquid bulk 4He, two-dimensional 4He films, and finite size 4HeN clusters  — have, of course, been studied extensively. This Letter thus presents the first predictions for spin-polarized atomic clusters. In particular, we characterize spin-polarized tritium clusters \[in the following denoted by (T$``$)N\] with up to $`N=40`$ tritium atoms, and compare their energetic and structural properties with those of bosonic 4HeN clusters. We hope that this study will stimulate further experimental work. Of particular interest is the lowest quartet state of the tritium trimer, which we predict to be a Borromean or halo state, and tritium cluster formation in the presence of an external magnetic field.

We also point out the possibility for creating an optically-pumped gaseous tritium condensate. Bose-Einstein condensates (BECs) are, to first order, well characterized by the two-body $`s`$-wave scattering length between two atoms. It will be shown below that the triplet two-body $`s`$-wave scattering length $`a_t`$ of two tritium atoms is large and negative, implying an unstable condensate of spin-polarized tritium atoms. We find, though, that there is an unusually broad Feshbach resonance  for two high-field-seeking, $`|F,m_F=|0,0`$ tritium atoms ($`F`$ denotes the total angular momentum, and $`m_F`$ the magnetic quantum number of this state).

While condensation of atomic hydrogen was realized experimentally in 1998 , it was a difficult experiment — at least in part because of hydrogen’s small triplet scattering length that limits the utility of evaporative cooling. Nevertheless, owing to hydrogen’s simplicity, it remains an important species to study. For instance, properties such as the interatomic potential and spin relaxation rates can be obtained theoretically from first principles. Unfortunately, we find no Feshbach resonance at reasonable field strengths for hydrogen. Thus, the resonance for tritium may permit faster condensation of a hydrogen-like atom, and allow for the formation of a stable BEC of tritium atoms with controllable properties in an optical dipole trap . Formation of such a tritium condensate should enhance the lively interplay between theory and experiment.

The behaviour of atomic tritium clusters and condensates is primarily determined by the two-body interaction potential for two tritium atoms, which is identical to that for two H or D atoms, except for the isotope-dependent adiabatic correction. Since there are only two electrons, these dimers are among the few for which highly accurate ab initio potentials are available. In the following, we concentrate on the singlet ground state ($`S=0`$, $`X^1\mathrm{\Sigma }_g^+`$), and on the triplet ground state ($`S=1`$, $`b^3\mathrm{\Sigma }_u^+`$) of the tritium dimer.

To construct the two-body $`S=0`$ and $`S=1`$ Born-Oppenheimer interaction potentials for hydrogen and its heavier isotopes, highly accurate ab initio data for the short range part  that incorporate the mass-dependent adiabatic correction are connected smoothly with an analytical expression describing the long-range behavior . This procedure results in six potential curves describing H2, D2, and T2 in their $`S=0`$ and $`S=1`$ states, respectively, which are then used in the radial Schrödinger equation describing the relative motion of a particle with reduced mass $`m/2`$. The mass dependence of the singlet and triplet scattering lengths has been discussed extensively in the literature. Here, we estimate the uncertainty of our two-body scattering observables by solving the radial Schrödinger equation using both the reduced atomic mass and the reduced nuclear mass for each two-body potential described above.

Figure 1 compares the (T$``$)2 potential (solid line) with the similarly shallow He2 ground state potential, LM2M2 from Aziz and Slaman  (dotted line). The T$``$ dimer and the He dimer potentials have a well depth of $`D_e=4.6`$cm-1 and $`7.6`$cm-1, respectively. Note that on the scale shown in Fig. 1, the H$``$ and D$``$ dimer potentials would be indistinguishable from the T$``$ dimer potential. The minimum of the (T$``$)2 potential lies at a significantly larger interparticle distance ($`r_e=7.8`$a.u.) than for the He dimer ($`r_e=5.6`$a.u.). Given the lighter mass of the tritium atom ($`m=5496.9`$a.u.) compared to that of the 4He atom ($`m=7296.3`$a.u.), it is not surprising that the T$``$ dimer is not bound, even though the tritium van der Waals coefficient $`C_6=6.499`$a.u. is larger than that for He, $`C_6=1.367`$a.u. (recall that the 4He dimer binding energy is only $`9.1\times 10^4`$cm-1 ). For comparison, the inset of Fig. 1 shows the tritium triplet potential (solid line) together with the tritium singlet potential (dashed line). The singlet curve is almost four orders of magnitude deeper than the triplet curve and supports 27 vibrational $`s`$-wave bound states.

In agreement with values tabulated in the literature , we calculate the two-body $`s`$-wave triplet scattering lengths $`a_t`$ for (H$``$)2 to be $`a_t=1.33`$a.u. \[$`1.33`$a.u.\], and for (D$``$)2 to be $`6.89`$a.u. \[$`6.88`$a.u.\], using the reduced atomic \[reduced nuclear\] mass. From symmetry considerations, the $`s`$-wave scattering length for (D$``$)2 is not an observable, and is given here for diagnostic purposes only. For tritium, we predict a positive singlet scattering length, $`a_s=34.6`$a.u. \[$`35.8`$a.u.\], and a large negative $`s`$-wave scattering length, $`a_t=82.1`$a.u. \[$`81.9`$a.u.\]. The T$``$ dimer does not possess a bound state, but this large negative $`a_t`$ indicates that it is only “slightly short of binding”.

Despite the fact that the tritium triplet scattering length is negative — implying an unstable condensate — it may be possible to form a stable tritium condensate utilizing a Feshbach resonance. Coupled-channel scattering calculations that couple the singlet and triplet subspaces reveal such a Feshbach resonance, i.e. a diverging scattering length for two atoms characterized by quantum numbers $`F`$ and $`m_F`$, as a function of the magnetic field strength. The coupling arises through the atomic hyperfine interaction, and has to be accounted for by an effective two-atom Hamiltonian . Feshbach resonances have been observed experimentally for 23Na and 85Rb  among others. The latter paved the way for a series of important BEC experiments, some of which enter the large interaction regime  while others probe the collapse regime .

We find that collisions between two tritium atoms, each in their $`|F,m_F=|0,0`$ state, results in a scattering lenth of $`a(00,00)=57`$a.u. \[64a.u.\] for zero magnetic field strength $`B`$, using the atomic \[nuclear\] mass. As $`B`$ increases, the scattering length $`a(00+00)`$ rises and eventually goes through infinity across a broad Feshbach resonance centered at $`B=870`$G \[810G\] (see Fig. 2). We also looked for, but could not find, a similar resonance for collisions of hydrogen atoms. Note that our predictions are not sensitive to a replacement of the atomic mass by the nuclear mass.

Formation of a tritium condensate in the high-field-seeking $`|0,0`$ state would require some variety of nonmagnetic trap, such as the dipole CO2 laser trap that has already produced a 87Rb condensate . Since the electric dipole polarizability of tritium is only 4.5 a.u., the 12 W CO2 laser setup of  would only produce trap depths of the order of a few microkelvin. There appears to be no reason why much stronger CO2 lasers could not be utilized, however. Realistically, the formation of a tritium condensate will probably require laser intensities at least an order of magnitude more intense. The large magnitude of the zero-field scattering length would make evaporative cooling far more effective than is the case for spin-polarized hydrogen. Another possible way to cool the spinless substate of tritium would be to implement a recent proposal to cool an atomic gas through magnetic field ramps across a Feshbach resonance . Despite these technical difficulties to be overcome, recent improvements in trapping and cooling technology would appear to make the creation of a tritium condensate a viable possibility.

To determine the bound state properties of (T$``$)N clusters with $`N>2`$, we first have to investigate the importance of non-additive contributions to the many-body interaction potential. Although three-body contributions, i.e., the Axilrod-Teller term  and three-body exchange terms , are significant for the hydrogen trimer in its electronic ground state, they should be less important for the spin-polarized trimer since the classical atom-atom equilibrium distance of the quartet trimer is more than five times as large as that for the doublet ground state trimer. Our calculations show that inclusion of the damped Axilrod-Teller term  raises the ground state energy of the spin-polarized tritium trimer by roughly 1.6%, and that of the larger clusters slightly more, e.g., by about 6 % for $`N=40`$. To describe (T$``$)N clusters, we assume in the following a simple pairwise additive potential energy surface, $`V=_{i<j}^NV(r_{ij})`$, where $`V(r)`$ denotes the triplet $`b^3\mathrm{\Sigma }_u^+`$ two-body potential. Conceivably, a more sophisticated many-body potential energy surface, which includes effects beyond the two-body potential, could modify our quantitative results somewhat, but we do not expect qualitative changes.

For the T$``$ trimer, we use the adiabatic hyperspherical representation . Including only one adiabatic channel yields a single bound state with energy $`1.60\times 10^3`$ cm-1. Coupling 25 adiabatic channels results in an energy of $`3.19\times 10^3`$ cm-1 with an uncertainty of 10-5 cm-1 and still no excited states. In the limit that an infinite number of channels are coupled, the bound state energy becomes exact, so the uncertainty is the result of including a finite number of channels. Upon inclusion of the Axilrod-Teller three-body term, the ground state energy is raised to $`3.14\times 10^3`$ cm-1 — again with uncertainty in the last digit. There is thus a single $`L=0`$ bound state for quartet tritium ($`L`$ is the total orbital angular momentum); no $`L>0`$ bound states are expected since none exist for the 4He trimer .

Since the T$``$ dimer is unbound, the T$``$ trimer is a Borromean state . One may then ask: does the T$``$ trimer state have Efimov character ? To investigate this question, we apply a simple quantitative criterion , although others exist . If the bound state disappears when one makes the potential more attractive (here achieved by simply multiplying the two-body potential with an overall scaling factor greater than 1), then the state under investigation is an Efimov state; if the bound state does not disappear, it is not an Efimov state. Applying this criterion, our coupled-channel calculations indicate that the bound T$``$ trimer state is not an Efimov state. In short, the T$``$ trimer has exactly one $`L=0`$ bound state, a Borromean state that is highly diffuse spatially.

To calculate the energetics and structural properties of (T$``$)N clusters with up to $`N=40`$ atoms, we employ the diffusion quantum Monte Carlo (DMC) technique . This method solves the time-independent many-body Schrödinger equation essentially exactly, to within a statistical error. Here, we employ the DMC method with importance sampling , using a descendant weighting scheme  for the extrapolation of structural properties. Our guiding wave functions , which enter the DMC calculation, have the analytical form given in Eq. (5) of Ref. , and recover between 84 and 96% of the DMC energy when used in a variational quantum Monte Carlo calculation. The DMC technique as implemented here determines only the ground state of the system; accessing excited state properties is in general not possible without algorithmic modifications.

For comparison, we find a DMC binding energy of $`2.9(5)\times 10^3`$cm-1 for the T$``$ trimer, in agreement with our hyperspherical calculation. The number in brackets denotes the statistical uncertainty. Figure 3 compares the ground state energy per particle $`E_0/N`$ of (T$``$)N clusters (pluses) with those of 4HeN clusters (using the LM2M2 potential , diamonds). This shows that (T$``$)N clusters are even more weakly bound than 4HeN clusters with the same number of atoms.

The comparison between (T$``$)N and 4HeN clusters can be extended by considering their structural properties. For example, we find that the (T$``$)5 cluster has an average interparticle distance of $`r_{ij}`$ of $`22.2`$a.u.; the 4He5 cluster, on the other hand, is $`r_{ij}=13.6`$a.u. Even if one takes into account that the classical equilibrium distance of the tritium triplet potential is about 2.2a.u. larger than that of the He dimer potential, the difference between the expectation values of the interparticle distance for these $`N=5`$ clusters indicates that the spin-polarized tritium system is even more diffuse than the 4He5 cluster. To illustrate this aspect further, the inset of Fig. 3 compares the pair distribution of the (T$``$)5 cluster (solid line) with that of the 4He5 cluster (dotted line). Clearly, the pair distribution of the (T$``$)5 cluster is much broader than that of the 4He5 cluster. We find similar behavior for clusters with more particles.

The tritium trimer is the smallest spin-polarized cluster (the dimer is unbound). As discussed above, the two-body potential for spin-polarized hydrogen is almost identical to that of tritium. The hydrogen atom, however, is about a factor of three lighter. Consequently, there is no bound state for the H$``$ trimer. An interesting question to ask is the following: how many atoms are needed to form a bound system of spin-polarized hydrogen atoms? Initial exploratory studies show that more than 100 atoms are needed to form a bound spin-polarized hydrogen cluster. The smallest spin-polarized hydrogen cluster could then be thought of as a “super-Borromean” cluster for which all smaller subsystems are unbound. A detailed study will be published elsewhere.

In summary, this Letter proposes a number of intriguing possibilities for the physics of tritium systems. We point out the possibility for forming a tritium condensate with controllable parameters via a Feshbach resonance. Down the road, one can imagine trapping an atomic hydrogen gas together with an atomic tritium gas, or possibly including deuterium to study fermion systems. We further found that the spin-polarized trimer possesses a Borromean, or halo, state. In addition, we mapped out the properties of larger “exotic” spin-polarized tritium clusters. Studies of spin-polarized clusters are interesting by themselves , as they enter new many-body physics regimes. For instance, unexpected physics may emerge from manipulating the two-body interaction via the Feshbach resonance in a cluster. In short, tritium offers a wealth of interesting physics by virtue of its weak attraction.

Acknowledgments: This work was supported by the National Science Foundation. BDE acknowledges support from the Research Corporation. We thank V. Kokoouline for help with the two-body coupled-channel calculations. Fruitful discussions with L. W. Bruch at an early stage of this project are acknowledged.


text:
# Some problems when calculating the quantum corrections to the classical ’t Hooft-Polyakov monopole

## 1 Introduction

While quantum corrections to the classical ’t Hooft-Polyakov monopole  have not been explicitly found, they have been for another soliton — the kink . It seems reasonable that the calculations for the monopole and kink should be similar. Assuming such a similarity, a suitable framework would start with the classical monopole, then determine the quantum fluctuations and, from these, derive the quantum corrections.

In this paper we draw attention to a feature of the quantum fluctuations that significantly differs between the monopole and kink. While closer examination may reveal this difference to be merely a technical issue, we think it could reflect a subtlety in the way semi-classical solitons behave.

To introduce our arguments, we briefly describe the quantum fluctuations around a kink soliton. Then we cover the analogous situation for the ’t Hooft-Polyakov monopole and indicate some problems. Finally, we discuss the cause of these problems.

## 2 The kink and its quantum corrections

The kink soliton  is a classical solution to a (1+1)-dimensional scalar field theory with Lagrangian density

$$[\mathrm{\Phi }]=\frac{1}{2}\dot{\mathrm{\Phi }}^2\frac{1}{2}\mathrm{\Phi }^2V[\mathrm{\Phi }],V[\mathrm{\Phi }]=\frac{1}{4}\lambda \left(\mathrm{\Phi }^2v^2\right)^2,$$
(1)
and has a classical scalar field profile

$$\mathrm{\Phi }_{\mathrm{cl}}(x)=\frac{m}{\sqrt{\lambda }}\mathrm{tanh}\frac{mx}{\sqrt{2}},m^2=\lambda v^2.$$
(2)
This kink has a classical mass defined from the Hamiltonian

$$M_{\mathrm{cl}}=H[\mathrm{\Phi }_{\mathrm{cl}}]=d^3x[\mathrm{\Phi }_{\mathrm{cl}}]=\frac{2\sqrt{3}m^3}{3\lambda }.$$
(3)
Here the Hamiltonian density is defined in the usual way

$$[\mathrm{\Phi },\mathrm{\Pi }]=\mathrm{\Pi }\dot{\mathrm{\Phi }}=\frac{1}{2}\mathrm{\Pi }^2+\frac{1}{2}\mathrm{\Phi }^2+V[\mathrm{\Phi }],$$
(4)
where $`\mathrm{\Pi }=/\dot{\mathrm{\Phi }}=\dot{\mathrm{\Phi }}`$ is the canonical momentum.

Quantum corrections to the classical kink originate from fluctuations $`\varphi `$ around the classical background

$$\mathrm{\Phi }=\mathrm{\Phi }_{\mathrm{cl}}+\varphi .$$
(5)
Then the Hamiltonian for the fluctuations is found by substituting (5) into (4) and expanding to $`\mathrm{O}(\varphi ^3)`$

$`H[\varphi ,\pi ]=M_{\mathrm{cl}}+{\displaystyle d^3r\left[\frac{1}{2}\pi ^2+\frac{1}{2}\varphi \left(_x^2+𝒱\right)\varphi \right]}`$ (6)
with $`\pi =\dot{\varphi }`$ the canonical momentum and $`𝒱`$ a potential for the fluctuations

$$𝒱(x)=3m^2\mathrm{tanh}^2\frac{mx}{\sqrt{2}}m^2.$$
(7)
This potential is a well with asymptote $`𝒱(\mathrm{})=2m^2`$.

It is crucial that the following operator is Hermitian:

$$d^3x\psi (_x^2+𝒱)\varphi =d^3x\varphi (_x^2+𝒱)\psi .$$
(8)
This relation follows from partial integration, neglecting surface terms. There is, therefore, a complete set of orthonormal eigenmodes

$$\left(_x^2+𝒱\right)\varphi _s=\omega _s^2\varphi _s,d^3x\varphi _s\varphi _t=\delta _{st}.$$
(9)
Expanding the fluctuation field as a sum $`\varphi =_sq_s\varphi _s`$ over these eigenmodes, we find

$$H(p_s,q_s)=M_{\mathrm{cl}}+_s\left[p_s^2+\omega _s^2q_s^2\right]+\mathrm{},$$
(10)
where $`q_s=d^3x\varphi \varphi _s`$; $`p_s=\dot{q}_s`$; the classical mass $`M_{\mathrm{cl}}`$ is from (3) and the dots represent $`\mathrm{O}(\varphi ^3)`$ terms.

Therefore the Hamiltonian (6) becomes a sum over individual Hamiltonians $`H(p_s,q_s)`$ that describe Harmonic oscillators. These oscillators can be quantized with Heisenberg’s commutation relation $`[p_s,q_t]=i\mathrm{}\delta _{st}`$, giving

$$H(n_{\omega _s})=M_{\mathrm{cl}}+_s(n_{\omega _s}+\frac{1}{2})\mathrm{}\omega _s+\mathrm{}.$$
(11)
Here $`n_{\omega _s}`$ is the occupation number of the Harmonic oscillator with frequency $`\omega _s`$.

In conclusion, the fluctuations around the classical kink can be written as a sum over eigenmodes that are stationary states of a time-independent Schrödinger equation (9). When these states are quantized, the meson spectrum has energies $`(n_{\omega _s}+\frac{1}{2})\mathrm{}\omega _s`$. Then the quantum correction to the classical kink mass is a sum over the differences between the meson’s zero-point and vacuum energies

$$M=M_{\mathrm{cl}}+_s\frac{1}{2}\mathrm{}\delta \omega _s+\mathrm{}.$$
(12)
Regularization techniques are used to evaluate this sum, which gives a finite value of order $`\mathrm{}m`$.

In the above discussion, there are a couple of points that will be relevant to our later arguments:
(a) The operator in (8) is Hermitian — then the fluctuations can be expanded over a complete set of eigenmodes.
(b) The quantum fluctuations solve a time-independent Schrödinger equation (9) — giving a connection between quantum field theory and quantum mechanics.

Later on in this paper, we compare these points with the behaviour of the ’t Hooft-Polyakov monopole.

## 3 The ’t Hooft-Polyakov monopole

Now we attempt to repeat the above arguments for a ’t Hooft-Polyakov monopole. Our aim is to see whether there are any significant differences between the kink and monopole.

The ’t Hooft-Polyakov monopole  is a classical solution to a (3+1)-dimensional $`\mathrm{SU}(2)`$ scalar-gauge field theory with Lagrangian density

$`[\mathrm{\Phi }_a,A_a^i]=\frac{1}{4}F_a^{\mu \nu }F_{a\mu \nu }+\frac{1}{2}D^\mu \mathrm{\Phi }_aD_\mu \mathrm{\Phi }_aV[\mathrm{\Phi }_a],`$ (13)
$`F_a^{\mu \nu }=^\mu A_a^\nu ^\nu A_a^\mu +eϵ_{abc}A_b^\mu A_c^\nu ,`$ (14)
$`D^\mu \mathrm{\Phi }_a=^\mu \mathrm{\Phi }_a+eϵ_{abc}A_b^\mu \mathrm{\Phi }_c,`$ (15)
$`V[\mathrm{\Phi }_a]=\frac{1}{4}\lambda (\mathrm{\Phi }_a\mathrm{\Phi }_av^2)^2,`$ (16)
and has classical scalar and gauge fields

$$(\mathrm{\Phi }_a)_{\mathrm{cl}}=\frac{H(r)}{er}\widehat{r}_a,(A_a^i)_{\mathrm{cl}}=\frac{K(r)1}{er}ϵ_{iab}\widehat{r}_b$$
(17)
in the radial $`ϵ_{abc}\mathrm{\Phi }_b\widehat{𝒓}_c`$ and temporal $`A_a^0=0`$ gauge. Single-valuedness requires $`H(0)=K(0)1=0`$, while the asymptotic boundary conditions are $`H/rev`$ and $`K0`$ as $`r\mathrm{}`$.

This monopole has magnetic charge $`4\pi /e`$ and vanishing electric charge, both of which are defined in a gauge-invariant way

$`g={\displaystyle \frac{1}{v}}{\displaystyle _S_{\mathrm{}}}𝑑S^i\mathrm{\Phi }_aB_a^i={\displaystyle \frac{4\pi }{e}},q={\displaystyle \frac{1}{v}}{\displaystyle _S_{\mathrm{}}}𝑑S^i\mathrm{\Phi }_aE_a^i=0.`$ (18)
These charges are defined by the non-Abelian electric and magnetic fields, $`E_a^i=F_a^{0i}`$ and $`B_a^i=\frac{1}{2}ϵ_{ijk}F_a^{jk}`$.

The classical mass of the monopole is defined from the Hamiltonian

$`M_{\mathrm{cl}}=H[(\mathrm{\Phi }_a)_{\mathrm{cl}},(A_a^i)_{\mathrm{cl}}]={\displaystyle d^3x[(\mathrm{\Phi }_a)_{\mathrm{cl}},(A_a^i)_{\mathrm{cl}}]},`$ (19)
which is specified by the Hamiltonian density

$$[\mathrm{\Phi }_a,\mathrm{\Pi }_a;A_a^i,\mathrm{\Pi }_a^i]=\mathrm{\Pi }_a\dot{\mathrm{\Phi }}_a+\mathrm{\Pi }_a^i\dot{A}_a^i[\mathrm{\Phi }_a,A_a^i]$$
(20)
with canonical momenta $`\mathrm{\Pi }_a=/\dot{\mathrm{\Phi }}_a=\dot{\mathrm{\Phi }}_a`$ for the scalar field and $`\mathrm{\Pi }_a^i=/\dot{A}_a^i=E_a^i`$ for the gauge field.

## 4 Fluctuations around the monopole

Quantum corrections to the classical monopole mass originate from the fluctuations $`\varphi _a`$ and $`a_a^i`$ around the classical fields $`(\mathrm{\Phi }_a)_{\mathrm{cl}}`$ and $`(A_a^i)_{\mathrm{cl}}`$. To ease the following calculations, we consider these fluctuations in the radial $`ϵ_{abc}\mathrm{\Phi }_b\widehat{r}_c=0`$ and temporal $`A_a^0=0`$ gauge

$`\mathrm{\Phi }_a=(\mathrm{\Phi }_a)_{\mathrm{cl}}+\varphi \widehat{r}_a,A_a^i=(A_a^i)_{\mathrm{cl}}+a_a^i.`$ (21)
This gauge choice reduces the scalar field to only one degree-of-freedom.

To describe these fluctuations, we construct the Hamiltonian from the fluctuation fields and their canonical momenta \[like with the kink in eqs. (5-7)\]. These canonical momenta are related to the fluctuations (21) by

$`\mathrm{\Pi }_a=(\mathrm{\Pi }_a)_{\mathrm{cl}}+\pi \widehat{r}_a=\pi \widehat{r}_a,\mathrm{\Pi }_a^i=(\mathrm{\Pi }_a^i)_{\mathrm{cl}}+\pi _a^i=\pi _a^i,`$ (22)
where $`\pi =\dot{\varphi }`$ and $`\pi _a^i=\dot{a}_a^i`$ in the temporal gauge. Then the Hamiltonian separates into zeroth-, linear- and quadratic-order terms in the fluctuations $`\varphi `$, $`\pi `$, $`a_a^i`$, $`\pi _a^i`$:

$`H[\mathrm{\Phi }_a,\mathrm{\Pi }_a;A_a^i,\mathrm{\Pi }_a^i]=H_{\mathrm{cl}}+H_{\mathrm{lin}}+H_{\mathrm{qu}}+\mathrm{}.`$ (23)
Clearly,

$$H_{\mathrm{cl}}=H[(\mathrm{\Phi }_a)_{\mathrm{cl}},(A_a^i)_{\mathrm{cl}}]=M_{\mathrm{cl}},$$
(24)
while it can be shown that

$$H_{\mathrm{lin}}=0.$$
(25)
This linear Hamiltonian vanishes by $`_{\mathrm{lin}}=0`$, which is because the monopole (17) solves the equations of motion and is then a stationary point of the action.

Therefore we are left with only a quadratic-order Hamiltonian $`H_{\mathrm{qu}}`$ for the fluctuations, which has a density

$$_{\mathrm{qu}}=\pi _a\pi _a+\pi _a^i\pi _a^i_{\mathrm{qu}}.$$
(26)
Here $`_{\mathrm{qu}}`$ contains the quadratic-order fluctuation terms in the Lagrangian density

$$_{\mathrm{qu}}=\frac{1}{4}\left(F_a^{\mu \nu }F_{a\mu \nu }\right)_{\mathrm{qu}}+\frac{1}{2}\left(D^\mu \mathrm{\Phi }_aD_\mu \mathrm{\Phi }_a\right)_{\mathrm{qu}}V_{\mathrm{qu}}.$$
(27)
We now determine each of these terms.

First, we consider the components of $`\frac{1}{4}F_a^{\mu \nu }F_{a\mu \nu }`$ that are quadratic in the fluctuations. These components are more easily found from the non-Abelian electric and magnetic fields

$`\frac{1}{4}(F_a^{\mu \nu }F_{a\mu \nu })_{\mathrm{qu}}`$ $`=`$ $`\frac{1}{2}(E_a^i)_{\mathrm{lin}}(E_a^i)_{\mathrm{lin}}+(E_a^i)_{\mathrm{cl}}(E_a^i)_{\mathrm{qu}}`$ (29)
$`\frac{1}{2}(B_a^i)_{\mathrm{lin}}(B_a^i)_{\mathrm{lin}}(B_a^i)_{\mathrm{cl}}(B_a^i)_{\mathrm{qu}}.`$
Here the zeroth-, linear- and quadratic-order fields are found from $`E_a^i=F_a^{0i}`$ and $`B_a^i=\frac{1}{2}ϵ_{ijk}F_a^{jk}`$. For the non-Abelian electric field, we find

$$(E_a^i)_{\mathrm{lin}}=\dot{a}_a^i=\pi _a^i,(E_a^i)_{\mathrm{qu}}=0;$$
(30)
while for the non-Abelian magnetic field,

$`(B_a^i)_{\mathrm{lin}}`$ $`=`$ $`ϵ_{ijk}^ja_a^k+eϵ_{ijk}ϵ_{abc}(A_b^j)_{\mathrm{cl}}a_c^k,`$ (31)
$`(B_a^i)_{\mathrm{qu}}`$ $`=`$ $`\frac{1}{2}eϵ_{abc}ϵ_{ijk}a_b^ja_c^k.`$ (32)
To keep things simple, we shall write the classical fields as $`(A_a^i)_{\mathrm{cl}}`$ and $`(B_a^i)_{\mathrm{cl}}`$. Substituting (30-32) into (29), we then find

$`\frac{1}{4}(F_a^{\mu \nu }F_{a\mu \nu })_{\mathrm{qu}}=\pi _a^i\pi _a^i\frac{1}{2}(B_a^i)_{\mathrm{lin}}(B_a^i)_{\mathrm{lin}}`$ (33)
$`(B_a^i)_{\mathrm{cl}}(B_a^i)_{\mathrm{qu}},`$ (34)
where upon partial integration and then neglecting the surface terms

$`(B_a^i)_{\mathrm{lin}}(B_a^i)_{\mathrm{lin}}`$ $`=`$ $`a_a^i\left(^2\delta ^{ij}^j^i\right)a_a^j`$ (35)
$`+`$ $`a_a^i\left[e^2ϵ_{ikm}ϵ_{jlm}ϵ_{acd}ϵ_{bce}(A_d^k)_{\mathrm{cl}}(A_e^l)_{\mathrm{cl}}\right]a_b^j`$ (36)
$``$ $`a_a^i\left[2eϵ_{ikm}ϵ_{abc}ϵ_{jlm}(A_c^k)_{\mathrm{cl}}^l\right]a_b^j,`$ (37)
$`(B_a^i)_{\mathrm{cl}}(B_a^i)_{\mathrm{qu}}`$ $`=`$ $`\frac{1}{2}eϵ_{abc}ϵ_{ijk}(B_c^k)_{\mathrm{cl}}a_a^ia_b^j.`$ (38)
Note the linear-derivative term in (35) — this will be troublesome later.

Next, we consider the components of $`\frac{1}{2}D^\mu \mathrm{\Phi }D_\mu \mathrm{\Phi }`$ that are quadratic in the fluctuations. These components are found by expanding

$`\frac{1}{2}(D^\mu \mathrm{\Phi }_aD_\mu \mathrm{\Phi }_a)_{\mathrm{qu}}`$ $`=`$ $`\frac{1}{2}\pi _a\pi _a\frac{1}{2}(D^i\mathrm{\Phi }_a)_{\mathrm{lin}}(D^i\mathrm{\Phi }_a)_{\mathrm{lin}}`$ (40)
$`(D^i\mathrm{\Phi }_a)_{\mathrm{cl}}(D^i\mathrm{\Phi }_a)_{\mathrm{quad}}.`$
Here we write $`\varphi _a=\varphi \widehat{r}_a`$, which gives $`D^0\varphi _a=\dot{\varphi }_a=\pi _a`$ in a temporal gauge. Similarly, the linear- and quadratic-order components of the covariant derivative are

$`(D^i\mathrm{\Phi }_a)_{\mathrm{lin}}`$ $`=`$ $`^i\varphi _a+eϵ_{abc}a_b^i(\mathrm{\Phi }_c)_{\mathrm{cl}}+eϵ_{abc}(A_b^i)_{\mathrm{cl}}\varphi _c,`$ (41)
$`(D^i\mathrm{\Phi }_a)_{\mathrm{qu}}`$ $`=`$ $`eϵ_{abc}a_b^i\varphi _c.`$ (42)
To keep things simple, we shall write the classical component of the covariant derivative as $`(D^i\mathrm{\Phi }_a)_{\mathrm{cl}}`$. Then the expressions in (40) are

$`(D^i\mathrm{\Phi }_a)_{\mathrm{lin}}(D^i\mathrm{\Phi }_a)_{\mathrm{lin}}=\varphi _a^2\varphi _a`$ (43)
$`+e^2ϵ_{abc}ϵ_{ade}\left[(A_b^i)_{\mathrm{cl}}(A_d^i)_{\mathrm{cl}}\varphi _c\varphi _e+(\mathrm{\Phi }_c)_{\mathrm{cl}}(\mathrm{\Phi }_e)_{\mathrm{cl}}a_b^ia_d^i\right]`$ (44)
$`+2e^i\varphi _aϵ_{abc}\left[a_b^i(\mathrm{\Phi }_c)_{\mathrm{cl}}+(A_b^i)_{\mathrm{cl}}\varphi _c\right]`$ (45)
$`+2e^2ϵ_{abc}ϵ_{ade}(\mathrm{\Phi }_c)_{\mathrm{cl}}(A_d^i)_{\mathrm{cl}}a_b^i\varphi _e,`$ (46)
$`(D^i\mathrm{\Phi }_a)_{\mathrm{cl}}(D^i\mathrm{\Phi }_a)_{\mathrm{qu}}=eϵ_{abc}(D^i\mathrm{\Phi }_a)_{\mathrm{cl}}a_b^i\varphi _c.`$ (47)
Again we find several terms that are linear in the field derivatives. However, this time the radial gauge expression $`^i\varphi _a=r^i(\varphi /r)\widehat{r}_a+(\varphi /r)\delta _{ia}`$ implies

$`eϵ_{abc}(A_b^i)_{\mathrm{cl}}(^i\varphi _a)\varphi _c`$ $`=`$ $`eϵ_{ibc}(A_b^i)_{\mathrm{cl}}\widehat{r}_c\varphi ^2/r,`$ (48)
$`eϵ_{abc}(\mathrm{\Phi }_c)_{\mathrm{cl}}(^i\varphi _a)a_b^i`$ $`=`$ $`eϵ_{ibc}(\mathrm{\Phi }_c)_{\mathrm{cl}}a_b^i\varphi /r.`$ (49)
Therefore these terms are not problematic.

Finally, there are also quadratic fluctuation terms in the scalar potential,

$$V_{\mathrm{qu}}[\varphi ]=\frac{3}{2}\lambda (\mathrm{\Phi }_a)_{\mathrm{cl}}(\mathrm{\Phi }_a)_{\mathrm{cl}}\varphi ^2\frac{1}{2}\lambda v^2\varphi ^2.$$
(50)
These are found by substituting (17) into $`V[\mathrm{\Phi }]`$.

Hence the final expression for the Hamiltonian is found by putting everything between (19) and (50) together into

$$H=M_{\mathrm{cl}}+H_{\mathrm{qu}}[\varphi ,\pi ;a_a^i,\pi _a^i]+\mathrm{}.$$
(51)
This $`H_{\mathrm{qu}}`$ splits into scalar, gauge and mixed parts

$$H_{\mathrm{qu}}=H_{\mathrm{qu}}[\varphi ,\pi ]+H_{\mathrm{qu}}[a_a^i,\pi _a^i]+I_{\mathrm{qu}}[\varphi ,\pi ;a_a^i,\pi _a^i],$$
(52)
specified by a Hamiltonian density for scalar fluctuations

$`_{\mathrm{qu}}[\varphi ,\pi ]`$ $`=`$ $`\frac{1}{2}\pi _a\pi _a\frac{1}{2}\varphi _a^2\varphi _a+\frac{1}{2}𝒱_\varphi \varphi ^2,`$ (53)
$`𝒱_\varphi (r)`$ $`=`$ $`e^2ϵ_{abc}ϵ_{ade}(A_b^i)_{\mathrm{cl}}(A_d^i)_{\mathrm{cl}}+3\lambda (\mathrm{\Phi }_a)_{\mathrm{cl}}(\mathrm{\Phi }_a)_{\mathrm{cl}}`$ (55)
$`+\mathrm{\hspace{0.17em}2}eϵ_{ibc}(A_b^i)_{\mathrm{cl}}\widehat{r}_c/r\lambda v^2;`$
a Hamiltonian density for gauge fluctuations

$`_{\mathrm{qu}}[a_a^i,\pi _a^i]`$ $`=`$ $`\frac{1}{2}\pi _a^i\pi _a^i\frac{1}{2}a_a^i\left(^2\delta ^{ij}^j^i\right)a_a^j`$ (57)
$`+\frac{1}{2}a_a^i𝒱_{ab}^{ij}a_b^j+\frac{1}{2}a_a^i𝒟_{ab}^{ij}a_b^j,`$
$`𝒱_{ab}^{ij}(r)`$ $`=`$ $`e^2ϵ_{ikm}ϵ_{jlm}ϵ_{bcd}ϵ_{ace}(A_d^k)_{\mathrm{cl}}(A_e^l)_{\mathrm{cl}}`$ (59)
$`+eϵ_{abc}ϵ_{ijk}(B_c^k)_{\mathrm{cl}}+e^2ϵ_{bcd}ϵ_{aed}(\mathrm{\Phi }_c)_{\mathrm{cl}}(\mathrm{\Phi }_e)_{\mathrm{cl}}\delta ^{ij},`$
$`𝒟_{ab}^{ij}`$ $`=`$ $`2eϵ_{ikm}ϵ_{abc}ϵ_{jlm}(A_c^k)_{\mathrm{cl}}^l;`$ (60)
and a scalar-gauge mixing term

$`_{\mathrm{qu}}[\varphi ,\pi ;a_a^i,\pi _a^i]`$ $`=`$ $`𝒱_a^ia_a^i\varphi ,`$ (61)
$`𝒱_a^i`$ $`=`$ $`e^2ϵ_{abc}ϵ_{bde}(\mathrm{\Phi }_c)_{\mathrm{cl}}(A_d^i)_{\mathrm{cl}}+eϵ_{abc}(D^i\mathrm{\Phi }_b)_{\mathrm{cl}}`$ (63)
$`+eϵ_{aic}(\mathrm{\Phi }_c)_{\mathrm{cl}}/r.`$
Collectively, the above equations (51-63) are analogous to (6) and (7) for the kink.

## 5 Problems with quantization

We now comment on these equations and indicate some differences from the kink in section 2.

The Hamiltonian (51) has a similar layout to that for the kink: with one part the classical mass and the other a quadratic Hamiltonian for the fluctuations. Therefore one might expect the next step would be to expand

$$\varphi =\underset{s}{}\varphi _sq_s,a_a^i=\underset{s}{}(a_a^i)_sQ_s$$
(64)
with $`\varphi _s`$, $`(a_a^i)_s`$ the eigenfunctions of

$`\left(^2+𝒱_\varphi \right)\varphi _s`$ $`=`$ $`\omega _s^2\varphi _s,`$ (65)
$`\left(\delta _{ab}^j^i^2\delta _{ab}\delta ^{ij}+𝒱_{ab}^{ij}+𝒟_{ab}^{ij}\right)(a_a^i)_s`$ $`=`$ $`\mathrm{\Omega }_s^2(a_b^j)_s.`$ (66)
Afterwards we could write each quadratic-order Hamiltonian as a sum over Harmonic oscillators. These oscillators could then be quantized.

However, there are a couple of obstacles to this procedure.

First, there is a a scalar-gauge mixing term (61). Then the eigenmodes are not separate scalar and gauge modes, but are instead a mixture of the two. This problem is mainly technical and makes any calculation more difficult.

Secondly, the linear derivative term $`𝒟_{ab}^{ij}`$ in (60) means the operator

$$^2\delta _{ab}\delta ^{ij}+\delta _{ab}^j^i+𝒱_{ab}^{ij}+𝒟_{ab}^{ij}$$
(67)
is not Hermitian because

$$d^3xa_a^i𝒟_{ab}^{ij}b_b^jd^3xb_a^i𝒟_{ab}^{ij}a_b^j.$$
(68)
In fact, examination of (68) reveals Hermiticity only holds for gauge perturbations that satisfy

$$(A_a^i)_{\mathrm{cl}}^la_b^j=(^iA_a^l)_{\mathrm{cl}}a_b^j+(A_a^l)_{\mathrm{cl}}^ia_b^j.$$
(69)
One cannot expect such an arbitrary constraint to hold for general gauge perturbations.

Hence it seems there is a fundamental problem with quantizing the fluctuations. For quantization, the fluctuation fields should be expanded over the eigenmodes. Then the Hamiltonian can be expanded as a sum over Harmonic oscillators $`H=_sp_s^2+\omega _s^2q_s^2`$, which can then be quantized with $`[p_s,q_t]=i\mathrm{}\delta _{st}`$. But if the operator (67) is not Hermitian there are problems with the eigenmode expansion. In particular:
(a) The eigenmodes may no longer be complete — then the expansions in (64) cannot be justified since there may be gauge fluctuations without an eigenmode sum.
(b) The eigenmodes may no longer be orthogonal — then the coefficients $`Q_s`$ may not be written as $`Q_s=d^3xaa_s`$. This problem is less serious than (a).
(c) We lose the connection between quantum field theory and quantum mechanics.

Without fully examining the eigenmodes of (67), which is a very difficult problem, one cannot really say how the analysis should be modified. It is possible that if there are fewer quantum fluctuations then there could be a smaller quantum mass-correction. In addition, there could be other effects that are difficult to anticipate.

We should mention that a more detailed calculation for the ’t Hooft-Polyakov monopole will be presented elsewhere . In that paper, many expressions given here will be presented in terms of the profile functions, rather than the classical fields. We will also estimate the quantum correction to the monopole mass.

## 6 Discussion

To complete our discussion, we examine when this loss of Hermiticity can occur in other situations by taking two illustrative examples: purely scalar field theory and Maxwell theory.

For the first example, we consider a scalar field theory with potential $`V[\mathrm{\Phi }]`$

$$[\mathrm{\Phi }]=\frac{1}{2}_\mu \mathrm{\Phi }^\mu \mathrm{\Phi }V[\mathrm{\Phi }].$$
(70)
We take an arbitrary classical background $`\mathrm{\Phi }_{\mathrm{cl}}`$ and introduce a fluctuation field $`\varphi `$ with

$$\mathrm{\Phi }=\mathrm{\Phi }_{\mathrm{cl}}+\varphi .$$
(71)
While $`\varphi \mathrm{\Phi }_{\mathrm{cl}}`$, the potential is to $`\mathrm{O}(\varphi ^3)`$

$`V[\mathrm{\Phi }]=V[\mathrm{\Phi }_{\mathrm{cl}}]+{\displaystyle \frac{V}{\mathrm{\Phi }}}[\mathrm{\Phi }_{\mathrm{cl}}]\varphi +{\displaystyle \frac{1}{2!}}{\displaystyle \frac{^2V}{\mathrm{\Phi }^2}}[\mathrm{\Phi }_{\mathrm{cl}}]\varphi ^2.`$ (72)
Hence the quadratic terms are

$$_{\mathrm{qu}}[\varphi ]=\frac{1}{2}_\mu \varphi ^\mu \varphi \frac{1}{2!}\frac{^2V}{\mathrm{\Phi }^2}[\mathrm{\Phi }_{\mathrm{cl}}]\varphi ^2,$$
(73)
and the quadratic Hamiltonian is

$$_{\mathrm{qu}}[\varphi ,\pi ]=\frac{1}{2}\pi ^2+\frac{1}{2}\varphi \left(^2+\frac{^2V}{\mathrm{\Phi }^2}[\mathrm{\Phi }_{\mathrm{cl}}]\right)\varphi .$$
(74)
In this purely scalar field theory, the relevant, bracketed, operator is Hermitian. Therefore the loss of Hermiticity for the monopole appears to be a gauge sector effect.

In our second example, we consider a similar analysis to above, but for Maxwell theory

$$=\frac{1}{4}F^{\mu \nu }F_{\mu \nu },F^{\mu \nu }=^\mu A^\nu ^\nu A^\mu .$$
(75)
Again we take an arbitrary classical background $`A_{\mathrm{cl}}^\mu `$ and introduce a fluctuation field $`a^\mu `$ with

$$A^\mu =A_{\mathrm{cl}}^\mu +a^\mu .$$
(76)
Substituting this into (75), we find a quadratic Lagrangian density

$$_{\mathrm{qu}}[a^\mu ]=\frac{1}{4}f^{\mu \nu }f_{\mu \nu },f^{\mu \nu }=^\mu a^\nu ^\nu a^\mu .$$
(77)
Because of the linearity of Maxwell theory, the fluctuation field decouples from the classical background. In Maxwell theory there is no problem with quantizing the fluctuations.

Therefore the loss of Hermiticity in the fluctuation Hamiltonian appears to be caused by the non-Abelian nature of the gauge theory. In purely scalar or Abelian theories no such loss of Hermiticity occurs.

For a final comment, we note this effect could be relevant to other classical solutions that contain non-Abelian gauge fields — an example being Skyrmions in the gauged sigma model.

<sup>1</sup><sup>1</sup>footnotetext:  email: n$`\mathrm{\_}`$lepora@hotmail.com


text:
# Growing smooth interfaces with inhomogeneous, moving external fields: dynamical transitions, devil’s staircases and self-assembled ripples

## Abstract

We study the steady state structure and dynamics of an interface in a pure Ising system on a square lattice placed in an inhomogeneous external field with a profile designed to stabilize a flat interface, and translated with velocity $`v_e`$. For small $`v_e`$, the interface is stuck to the profile, is macroscopically smooth, and is rippled with a periodicity in general incommensurate with the lattice parameter. For arbitrary orientations of the profile, the local slope of the interface locks in to one of infinitely many rational values (devil’s staircase) which most closely approximates the profile. These “lock-in” structures and ripples dissappear as $`v_e`$ increases. For still larger $`v_e`$ the profile detaches from the interface.

The ability to grow flat solid surfacesmagnet  is often of major technological concern, for example, in the fabrication of magnetic materials for recording devices where surface roughnessdegrade  causes a sharp deterioration of magnetic properties. Most growing surfaces or interfaces on the other hand coarsenBarabasi  with a width, $`\sigma `$, which diverges with system size as $`L^\alpha `$ and time as $`t^\beta `$, where $`\alpha `$ and $`\beta `$ ($`\alpha ,\beta 0`$) are the roughness and dynamical exponents respectively. Is it possible to drive such an interface with a pre-determined velocity $`v_f`$ and, at the same time, keep it flat (i.e. $`\alpha =\beta =0`$) ? In this Letter, we explore this possiblity by studying growing interfaces in a non-uniform field with an appropriately shaped profile moving without change of shape at an externally controllable velocity $`v_e`$. We find that for $`v_e`$ less than a limiting value $`v_{\mathrm{}}`$, it is possible to produce a macroscopically flat, interface growing with average velocity $`v_f=v_e`$. Microscopically, however, the interface shows an infinity of dynamical ripple structures similar to self-assembled, commensurate -incommensurate(C-I) FK  domains produced by atomic mismatchmismatch . The ripples vanish with increasing $`v_e`$ through a fluctuation induced C-I transition.

Real solid-solid interfaces being complexsutton , an understanding of the dynamics of such interfaces in a general time -dependent, inhomogeneous field may only be achieved by beginning with a relatively simple, but nontrivial, system viz. an interface in the two -dimensional (2-d) Ising model. Our results, therefore, concern mainly this model system, though towards the end we discuss possible modifications, if any, for solid interfaces.

The dynamics of an Ising interface in a (square) lattice driven by uniform external fields is a rather well studiedBarabasi  subject. The velocity, $`v_{\mathrm{}}`$ of the interface depends on the applied field, $`h`$ and the orientation $`\theta `$ measured with respect to the underlying lattice. The interface is rough and coarsens with KPZBarabasi  exponents $`\alpha =1/2`$ and $`\beta =1/3`$.

Consider, an interface $`Y(x,t)`$ between phases with magnetization, $`\varphi (x,y,t)>0`$ and $`\varphi (x,y,t)<0`$, in a 2-d square latticemexpln  obeying single-spin flip Glauber dynamics Men  in the limit $`h/J,T/J0`$. Here $`J`$ is the Ising exchange coupling and $`T`$ the temperature. An external non-uniform field with a profile $`h(y,t)=h_{max}f(y,t)`$ where $`f(y,t)=\mathrm{tanh}((yv_et)/\chi )`$ and $`\chi `$ is the width of the profile (see Fig 1(a)) is applied such that $`h=h_{\mathrm{𝑚𝑎𝑥}}`$ in the $`\varphi >0`$ and $`h_{\mathrm{𝑚𝑎𝑥}}`$ in the $`\varphi <0`$ regions separated by a sharp edge. The driving force depends on the relative local position of $`Y(x,t)`$ and the edge. In the low temperature limit the interface moves solely by random corner flipsBarabasi  (Fig. 1(b)), the fluctuations necessary for nucleating islands of the minority phase in any region being absent. We study the behaviour of $`v_f`$ and the structure of the interface as a function of $`v_e`$ and orientation.

Naively, one expects fluctuations of the interfacial coordinate $`Y(x,t)`$ to be completely suppressed in the presence of $`h(y,t)`$. Indeed, as we show below, a mean field theory obtains the exact behaviour of $`v_f`$ as a function of $`v_e`$ (Fig. 1 (c)). Using model A dynamicsChaikin  for a coarse-grained Hamiltonian (which, for the moment, ignores the lattice) of an Ising system in a external non-uniform field together with the assumption that the magnetization $`\varphi `$ is uniform everywhere except near the interface, one can derivephysica  an equation of motion for the interface.

$$\frac{Y}{t}=\lambda _1\frac{^2Y}{x^2}\lambda _2\left(\frac{Y}{x}\right)^2f(Y,t)\lambda _3f(Y,t)+\zeta (x,t)$$
(1)
where $`\lambda _1`$,$`\lambda _2`$ and $`\lambda _3`$ are parameters and $`\zeta `$ is a Gaussian white noise with zero mean. Note that Eq. (1) lacks Galilean invariancemultif  $`Y^{}Y+ϵx,x^{}x\lambda _2ϵt,t^{}t`$. A mean field calculation amounts to taking $`YY(t)`$ i.e. neglecting spatial fluctuations of the interface and noise. For large times ($`t\mathrm{}`$), $`Yv_ft`$, where $`v_f`$ is obtained by solving the self-consistency equation; $`v_f=lim_t\mathrm{}\lambda _3\mathrm{tanh}\{(v_fv_e)t/\chi \}=\lambda _3\mathrm{sign}(v_fv_e)`$. For small $`v_e`$ the only solution (Fig. 1(c), inset) is $`v_f=v_e`$ and for $`v_e>v_{\mathrm{}}`$, where $`v_{\mathrm{}}=\lambda _3`$ we get $`v_f=\lambda _3=v_{\mathrm{}}`$. We thus have a sharp transition (Fig. 1(c)) from a region where the interface is stuck to the edge to one where it moves with a constant velocity. How is this result altered by including spatial fluctuations of $`Y`$ ? This question is best answered by mapping the interface problem to a 1-d cellular automatonBarabasi ; AEP .

The interface coordinate $`Y(x,t)`$ in a square lattice is representedAEP ; Maj1  by the set of integers $`\{y_i\},\mathrm{\hspace{0.17em}\hspace{0.17em}1}<i<N_p`$ denoting positions of $`N_p`$ hard-core ($`y_{i+1}y_i+1`$) particles in a 1-d lattice of $`N_s`$ sites. The particle density $`\rho =N_p/N_s`$ determines the mean slope of the interface $`\mathrm{tan}\theta _f=1/\rho `$ and the motion of the interface by corner flips corresponds to the hopping of particles with right and left jump probabilities $`p`$ and $`q`$ ($`p+q=1`$). Trial moves are attempted sequentially on randomly chosen particlesAEP  and $`N_p`$ attempted hops constitute a single time step. In our case, $`p`$ and $`q`$ are position dependent such that, $`\mathrm{\Delta }_i(v_et)=pq=\mathrm{\Delta }\mathrm{sign}(y_ii/\rho v_et)`$ with $`\mathrm{\Delta }=1`$. The bias $`\mathrm{\Delta }_i(v_et)`$ is appropriate for a step function ($`\chi =0`$) profile with the slope of the profile edge equal to the average slope of the interface. We track the instantaneous particle velocity $`v_f(t)`$ defined as the number of particles moving right per unit time, the average position $`<y(t)>=N_p^1_{i=1,N_p}y_i(t)`$, the width $`\sigma ^2(t)=N_p^1_{i=1,N_p}<(y_i(t)<y_i(t)>)^2>`$ and the local slope of the interface given by the local density of particles. Angular brackets denotes an average over the realizations of the random noise. The usual particle hole symmetry for an exclusion processBarabasi ; AEP ; Maj1  is violated since exchanging particles and holes alters the relative position of the interface compared to the edge.

We perform numerical simulationsAEP  of the above model for $`N_s`$ upto $`10^4`$ to obtain $`v_f`$ for the steady state interface as a function of $`v_e`$ as shown in Fig. 1(c). A sharp dynamical transition from an initially stuck interface with $`v_f=v_e`$ to a free, detached interface with $`v_f=v_{\mathrm{}}=\mathrm{\Delta }(1\rho )`$ is clearly evident as predicted by mean field theory. The detached interface coarsens with KPZ exponentsphysica . Note that, even though the mean field solution for $`v_f(v_e)`$ neglects the fluctuations present in our simulation, it is exact. The detailed nature of the stuck phase ($`v_f=v_e`$ and $`\sigma `$ bounded) is, on the other hand, considerably more complicated than the mean field assumption $`Y(x,t)=Y(t)`$.

The ground state of the interface in the presence of a stationary ($`v_e=0`$) field profile is obtained by minimizing $`E=1/N_p_i(y_ii/\rho c)^2`$ with respect to the set $`\{y_i\}`$ and the constant, $`c`$, subject to the constraint that $`y_i`$ are integers. The form of $`E`$ leads to an additional non-local, repulsive, interaction between particles. The minimized energy may be calculated exactly, $`E(\rho =m/n)=\frac{1}{6}(\frac{1}{2}\frac{1}{m})(1\frac{1}{m})+\frac{1}{4m}\frac{1}{4m^2}`$ for $`m`$ even and $`\frac{1}{12}(1\frac{1}{m^2})`$ for $`m`$ odd, where the density $`\rho =m/n`$, is a rational fraction. The energy satisfies the bounds $`E(1/n)=0<E(m/n)<lim_m\mathrm{}E(m/n)=1/12`$ where the upper bound is for an irrational density. For an arbitrary $`0<\rho <\mathrm{\hspace{0.17em}1}`$ the system ($`\{y_i\}`$) prefers to distort, conforming within local regions, to the nearest low-lying rational slope $`1/\stackrel{~}{\rho }`$ interspersed with an ordered array of “discommensurations” of density $`\rho _d=|\rho \stackrel{~}{\rho }|`$ which appear as long wavelength ripples (see Fig. 3 inset (c)). The $`\stackrel{~}{\rho }`$ as a function of $`\rho `$ shows a “devil’s staircase” structure (complete for $`v_e0`$) with a multifractal CDW  measure. We observe this in our simulation by analyzing the instantaneous distribution of the local density of particles to obtain weights for various simple rational fractions upto generation $`g=9`$ in the Farey tree of rationalsfarey . A time average of the density corresponding to the fraction with largest weight at any $`t`$, then give us the most probable density $`\stackrel{~}{\rho }`$ — distinct from the average $`\rho `$ which is constrained to be the inverse slope of the profile. For small $`v_e`$ the interface is more or less locked in to a single $`\stackrel{~}{\rho }`$, shown as white regions in the phase diagram (Fig. 2) in the $`v_e\rho `$ plane.

For low velocities and density where correlation effects due to the hard core constraint are negligible, the dynamics of the interface may be obtained exactlyphysica . Under these circumstances the $`N_p`$ particle probability distribution for the $`y_i`$’s, $`P(y_1,y_2,\mathrm{},y_{N_p})`$ factorizes into single particle terms $`P(y_i)`$. Knowing the time development of $`P(y_i)`$ and the ground state structure the motion of the interface at subsequent times may be trivially computed as a sum of single particle motions. A single particle (with say index $`i`$) moves with the bias $`\mathrm{\Delta }_i(v_et)`$ which, in general, may change sign at $`y<i/\rho +v_et<y+1`$. Solving the appropriate set of master equationsphysica  we obtain, for $`v_e<<1`$ the rather obvious steady state solution $`P(y_i)=1/2(\delta _{y_i,y}+\delta _{y_i,y+1})`$ and the particle oscillates between $`y`$ and $`y+1`$. Subsequently, when $`i/\rho +v_ety+1`$, the particle jumps to the next position and $`P(y_i)`$ relaxes exponentially with a time constant $`\tau =1`$ to it’s new value with $`yy+1`$. In general, for rational $`\rho =m/n`$, the motion of the interface is composed of the independent motions of $`m`$ particles each separated by a time lag of $`\tau _L=1/mv_e`$. The result of the analytic calculation for small $`v_e`$ and $`\rho `$ has been compared to those from simulations in Fig. 3 for $`\rho =1/5`$ and $`2/5`$. For a general irrational $`\rho <1/2`$, $`m\mathrm{}`$ consequently, $`\tau _L0`$. The $`y_i`$’s are distributed uniformly around the mean implying $`\sigma ^2=1/3`$ independent of system size and time. For $`\rho >1/2`$ the width $`\sigma ^2=(1\rho )/3\rho `$ since the number of mobile particles decreases by a factor of $`(1\rho )/\rho `$.

The forward motion of an irrational interface is accompanied by the motion of discommensurations along the interface with velocity $`v_e`$ which constitutes a kinematic waveBarabasi ; Maj1  parallel to the interface.

As the velocity $`v_e`$ is increased the system finds it increasingly difficult to maintain its ground state structure and for $`\tau \tau _L`$ the instantaneous value of $`\stackrel{~}{\rho }`$ begins to make excursions to other nearby low-lying fractions and eventually becomes free. Steps corresponding to $`\stackrel{~}{\rho }=m/n`$ dissappear (i.e. $`\stackrel{~}{\rho }\rho `$) sequentially in order of decreasing $`m`$ and the interface loses the ripples. The interface is disordered though $`\alpha `$ and $`\beta `$ continue to be zero (black region in Fig. 2). The locus of the discontinuities (within an accuracy of $`1/N_p`$) in the $`\stackrel{~}{\rho }(\rho )`$ curve for various velocities $`v_e`$ gives the limit of stability of the lock-in rippled phases.

While the stability of mismatch domainsmismatch  is decided, mainly, by competition between mechanical, long-ranged (elastic) and short-ranged (atomistic) interactionsFK , dynamical ripples vanish with increasing $`v_e`$ through increased fluctuations. We argue that it is sufficient to project the entire configuration space $`y_i`$ of the stuck interface onto the single variable $`\rho `$. In Fig. 4 (inset) we have plotted the energy $`E(\rho )`$ for the ground-state configuration with density $`\rho `$. As is obvious from the figure, $`E(\rho )`$ has a structure similar to the free energy surface of a 1-d “trap” modeltrap  often used to describe glassy dynamics. The distinction, of course, is the fact that the energies of the traps in this case are highly correlated. We may then describe the dynamics of the stuck interface as the Langevin dynamicsChaikin  of a single particle with coordinate $`\rho ^{}`$ diffusing on a energy surface given by, $`F(\rho ^{})=E(\rho ^{})+\kappa (\rho ^{}\rho )^2,`$ kicked by a Gaussian white noise of strength $`T`$. The second term, containing the modulus $`\kappa `$, ensures that $`\rho ^{}\rho `$ for $`t\mathrm{}`$. At intermediate times, however, the system may get trapped indefinitely in some nearby low-lying minimum with $`\rho ^{}=\stackrel{~}{\rho }`$ if $`T`$ ($`v_e^2`$ by symmetry) is not large enough. As $`T`$ increases, the time spent in jumping between minima may exceed the residence time in the minimum resulting in a noise induced C-I transition (Fig. 4) from $`\stackrel{~}{\rho }`$.

In this Letter, we have studied the static and dynamic properties of an Ising interface in 2-d subject to a non-uniform, time-dependent external magnetic field. The system has a rich dynamical phase diagram with infinitely many steady states. The nature of these steady states and their detailed dynamics depend on the orientation of the interface and the velocity of the external field profile. How are our results expected to change for real driven solid interfaces? Firstly, real field profiles would have a finite width $`\chi >0`$. However, as long as $`\chi `$ is comparable to atomic dimensions we find no appreciable change in the results. Secondly, elastic interaction between “particles” or steps in the interfacesutton  may smoothen the devil’s staircase though we expect that for typical solids this will be minimal. Thirdly, the structure of the underlying lattice may influence the stability of particular orientations. Finally, exciting new physics may come into play as new modes e.g. point and line defectsmart , as well as phonon degrees of freedom (leading to acoustic emmissionssutton ) are accessed as $`v_e`$ increases.
Acknowledgement:  The authors thank M. Barma, J. K. Bhattacharya, J. Krug, S. N. Majumdar, A. Mookerjee and M. Rao for useful discussions; A. C. thanks C.S.I.R., Govt. of India for a fellowship.