text:
# Reassessment of the GZK cutoff in the spectrum of UHE cosmic rays in a universe with low photon-baryon ratio

## 1 Introduction

Experiments to measure the flux of UHECR have detected a number of events with energies $`>10^{20}eV`$ (Watson, 2000). These observations present something of a puzzle. There are no obvious mechanisms within our galaxy that could accelerate UHECRs to these energies. It is reasonable to postulate potential extra-galactic sources of UHECRs, e.g. quasars, but prevailing theory suggests that these particles should never reach the Earth because of interactions with the photons that constitute the cosmic microwave background radiation (CMBR). This prediction was made independently by Greisen (1966) and by Zatsepin & Kuzmin (1966), based on an analysis of the process $`p+\gamma p+\pi `$ (assuming that the majority of UHECRs at these energy levels are protons). This showed that UHECRs with energies $`>5\times 10^{19}`$ would exceed the energy threshold for photopion production (the GZK cutoff). Taking into account the collision cross section for the photopion reaction, and the CMBR photon density, the mean free path for these UHECRs would be of the order of $`10^{23}m`$ (see Appendix A for analysis), implying that they could not have originated from sources outside our galaxy.

The present generation of UHECR experiments, such as AGASA and HiRes, are not able to generate a large enough sample from the relatively low flux of UHECRs at these high energies (De Marco, Blasi & Olinto, 2003). The low statistical significance of the few observations of UHECRs with energies $`>10^{20}eV`$ means that no firm conclusions can be drawn at this stage. Larger experiments are planned, such as the Pierre Auger Observatory (PAO) and EUSO, which will increase the number of detected events by 2 orders of magnitude. Only then will it be possible to to prove conclusively whether or not the predicted GZK cutoff feature is present in the UHECR spectrum. Nevertheless, it is reasonable to evaluate potential cosmological scenarios that might result in the detection of significant UHECR flux above the presumed GZK cutoff.

There have been suggestions that the observed UHECRs with energies $`>10^{20}eV`$ might be due to Lorentz symmetry breaking (Kifune, 1999; Magueijo & Smolin, 2002), or deformation of Lorentz symmetry (Amelino-Camelia, 2002). These predictions are based on a class of models that postulate deviations from standard special relativity for particles that have attained very high energies, such that $`m/R<E/E_P`$, where $`m`$ and $`E`$ are respectively the particle mass and total energy, and $`E_P`$ is the Planck energy. These are generally referred to as Double Special Relativity (DSR) models. The main criticism of this approach is that it requires either a preferred class of inertial observers, or an ad-hoc construction of the Lorentz deformation function to generate predicted effects that would result in a shift of the GZK cutoff.

An alternative solution to the GZK cutoff problem is proposed in the following sections, based on a significantly lower CMB photon number density that arises naturally from an alternative cosmological model (Booth, 2002).

## 2 The Quasi-Static Universe model

The Quasi-Static Universe (QSU)is used here as a shorthand to refer to the paradigm described in Booth (2002). The main feature of this model is that the Planck scale is decoupled from the atomic scale conventionally used as the as the basis for our measurement system. The changes in the dynamical behaviour of the universe that result from this modification can not only account for most of the problems associated with the Big Bang model, but are also able provide a very simple explanation for the apparent acceleration of the expansion rate of the universe that has been observed in various high redshift supernovae studies in recent years. One of the principal consequences of the QSU paradigm is that the atomic scale, as defined by the de Broglie wavelength of sub-atomic particles, is not an appropriate reference frame for measuring gravitational phenomena or the behaviour of photons. The QSU model is based on the postulate that the correct reference frame for these phenomena is in fact a cosmological frame based on the mass and size of universe as a whole. In such a frame, photons do not undergo any change in frequency, since as far as they are concerned, the universe is static. Hence, it is not meaningful to talk in terms of photons losing energy in this reference frame. In transforming from the cosmological reference frame to our conventional atomic frame, photons will be perceived to exhibit a redshift as the scale factor of the universe with respect to the atomic frame increases with time. The crucial difference between the QSU model and the conventional formalism is that the relationship $`E_\gamma =h\nu `$ no longer holds true for photons emitted at times $`t<t_0`$, where $`t_0`$ signifies the present time. The energy of such ’old’ photons will remain constant, at the same value they possessed when they were first emitted. However their power, i.e. the energy transferred per unit time, will be reduced in proportion to their redshift, such that $`W=W_0.\nu \mathrm{/}\nu _0`$ , where $`W_0`$ and $`\nu _0`$ are respectively the photon power and frequency at the time of emission. Clearly, such a modification to one the most fundamental equations in physics will have a very significant impact on any phenomena that involve redshifted photons, and in particular, the CMBR.

To understand the implications of the QSU model for the CMBR, we need to review the way in which the standard Planck black-body distribution law is applied. A summary of the standard derivation of this law is provided in Appendix B. Present day observations of the CMB give a value for the energy density of $`U4\times 10^{14}Jm^3`$, which from (B4) corresponds to a temperature of $`T=2.7^{}K`$. Conventionally, the next step is to take this temperature and use equation (B6) to calculate the photon number density, giving a result of $`N4\times 10^8m^3`$. However, these formulae are only valid for black-body radiation *that is in equilibrium with its surroundings*. It would be perfectly correct to use these equations if we wished to deduce the photon number density for a black-body with this temperature *today*. In erroneously applying them in the context of the relic CBMR generated by the Big Bang, we are perpetuating the assumption inherent in going from (B3) to (B5) - that photon energy is always equal to $`h\nu `$.

In the QSU model, red-shifted photons do not lose energy, so it is not possible to determine the energy of an observed photon merely by measuring its frequency. It is also necessary to know the thermal history of the photon, i.e. its frequency when it was originally emitted. It is this that determines the energy of the photon. In order to obtain the correct result for the photon number density today that corresponds to an observed energy density, we need to correct for the fact that the photon energy $`h\nu `$ applied in going from (B3) to (B5) should reflect the energy at the time of emission, or more precisely, its energy at the time it was in black-body equilibrium with its surroundings.

We therefore need to apply a factor of $`T_{obs}/T_{equi}`$ to the expression for photon number density in (B5), where of $`T_{obs}`$ is the observed absolute temperature of the CMBR today, and $`T_{equi}`$ is the temperature at the time of black-body equilibrium. Clearly we do not know this figure with any precision. However we can take an educated guess that is at least consistent with other observational data and with plausible models for primordial nucleosynthesis. One such model is the neutron decay variant of the Cold Big Bang (CBB), which predicts a maximum reaction temperature of $`10^{10}K`$, with an equilibrium temperature of the order of $`10^9K`$. Applying the correction factor to equation (B6) gives a calculated photon number density of $`0.3m^3`$ \- a factor of $`10^9`$ lower than for the standard Hot Big Bang (HBB) model. This is very close to the measured baryon number density of the universe, giving $`\eta _\gamma 1`$.

## 3 Implications for cosmological processes

Having established that the value of $`\eta _\gamma `$ in the QSU model is many orders of magnitude less than the conventionally accepted value, the natural question to ask is: how does this affect other cosmological processes?

### 3.1 Primordial nucleosynthesis

Arguably, one of the process most sensitive to changes in $`\eta _\gamma `$ is that of primordial nucleosynthesis. Applying the conventional value of $`\eta _\gamma 2\times 10^9`$ to the standard model for HBB nucleosynthesis results in predicted element abundances that are in good accord with observational data. Relatively small changes in $`\eta _\gamma `$ will result in large changes to the predicted abundances, and would therefore appear not to be compatible with observations. However, it has been pointed out by Aguirre (2001) that, provided that appropriate changes are made to a range of initial parameters, including $`\eta _\gamma `$, it is possible to construct alternative models for primordial nucleosynthesis that will produce predicted element abundances that are in accord with observational data. One such model is the CBB, which takes a value of $`\eta _\gamma 1`$ as one of its initial conditions. It is perhaps worth mentioning in passing that the neutron decay variant of the CBB model predicts a photon energy of $`0.78MeV`$, giving a baryon to photon energy ratio of $`E_\gamma /E_B1200`$. Under the QSU scenario, this ratio should persist from the nucleosynthesis epoch to the present day, and indeed, the observed value of $`E_\gamma /E_B`$ is very close to this value.

### 3.2 Hydrogen ionization

Another astrophysical measurement that is linked to $`\eta _\gamma `$ is the hydrogen ionization fraction. As the universe expands and cools, protons recombine with electrons to form neutral hydrogen when the energy of CMB photons falls below the hydrogen ionization energy threshold. The equilibrium ionization fraction $`\chi _e`$ as a function of temperature is give by the Saha equation

$$\frac{1\chi _e}{\chi _e}=\frac{4\sqrt{2}\varsigma (3)}{\sqrt{\pi }\eta _\gamma }\left(\frac{T}{m_e}\right)^{\frac{3}{2}}e^{B\mathrm{/}T}$$
(1)
from which it can be see that the ionization fraction is also dependent on $`\eta _\gamma `$. The ionization fraction can be expressed as a function of redshift using the relation $`T=2.73(1+z)K`$. This is plotted in Figure 1 for a standard HBB cosmology with $`\eta _\gamma =2\times 10^9`$, and a CBB cosmology with $`\eta _\gamma =1`$. This illustrates that a reduction in the photon to baryon ratio will cause recombination to take place at a much lower redshift, with the surface of last scattering for the CMB occurring at $`z2400`$, as compared to $`z1100`$ for the standard model.

At first sight it might be expected that such a change should inevitably result in very significant modifications to the observed CMB angular power spectrum. However, an analysis of the dynamics of the two cosmological models reveals that the linear expansion of a CBB universe, in comparison to the rapid initial expansion rate of a standard model HBB universe, results in the epoch of creation of the last scattering surface being similar, in terms of coordinate time, for both cosmologies. Further detailed numerical analysis is required, using a package such as CMBFAST, in order to accurately determine the nature of the impact on the predicted CMB angular power spectrum arising from a $`\eta _\gamma =1`$ cosmological model. Nevertheless, initial indications are that such a model may well be consistent with the observed CMB angular power spectrum.

### 3.3 Reassessment of the GZK cutoff

The interaction between CMB photons and UHECRs is dependent on three factors: the photon energy, the capture cross-section for the photon-proton collision, and the photon number density. We shall now examine how each of these factors varies in a CBB universe.

From (A4) it might be anticipated that if each CMB photon in a QSU is a factor of $`10^9`$ more energetic than its counterpart in a standard HBB universe, then the threshold energy for the photon-proton reaction will be reduced from $`10^{20}eV`$ to $`10^{11}eV`$. However, this is unlikely to be the case here since energy transfer from the CMBR photon to the proton will be dependent on the photon power, rather than its total energy. Since the photon power is, as in the standard case, proportional to $`h\nu `$, there will be no change in the GZK cutoff threshold under the QSU scenario.

The next factor to consider is the collision cross-section for the photon-proton interaction. As discussed in Appendix A, this varies by a factor of $`5`$ between energies corresponding to the main $`\mathrm{\Delta }^+`$ resonance peak and higher photon energy levels. It is not, therefore, a significant factor in differentiating between the two cosmological scenarios.

Finally, we must consider the impact of the CMB photon number density on the mean free path of the UHECRs. Under the QSU model, the CMBR photon number density is a factor of $`10^9`$ less than for the HBB model. From (A6), taking a value for $`\sigma =100\mu \text{barn}`$, and $`n=0.2m^3`$, gives a value for the mean free path of $`\lambda =10^{32}m`$. This is several orders of magnitude greater than the horizon distance of the universe. The universe is therefore effectively transparent to UHECRs, and the QSU model predicts that there will be no observable reduction in the flux of UHECRs as a result of the GZK cutoff.

## 4 Conclusions

The reduced CMB photon number density predicted by the QSU model provides a natural explanation for the observed UHECR flux with energies $`>10^{20}eV`$. However, the statistical significance of the measurements obtained from the existing AGASA and HiRes experiments is still to low to provide conclusive evidence for the presence or absence of the GZK cutoff feature. For this, we will have to wait for the results from the next generation of UHECR experiments, such as PAO and EUSO.

## Appendix A The photon-proton reaction

The principal reaction contributing to the attenuation of UHECRs is the interaction with CMB photons

$`p+\gamma `$ $``$ $`\mathrm{\Delta }^+p+\pi ^0`$ (A1)
$`p+\gamma `$ $``$ $`\mathrm{\Delta }^+n+\pi ^+`$ (A2)
The reaction will proceed when the combined centre-of-mass energy of the proton and photon is equal to or greater than the sum of the pion and proton (or neutron) mass, which can be expressed as

$$m_pm_\pi +\frac{m_\pi ^2}{2}q\left(\sqrt{p^2+m_p^2}p\mathrm{cos}\theta \right)$$
(A3)
where $`q`$ is the photon momentum along the $`x`$-axis and $`p`$ is the momentum of the proton hitting the photon at an angle of $`\theta `$ in the $`xy`$ plane. Since the pion mass is much smaller than the proton (or neutron) mass, this expression can be simplified to

$$E_pp\mathrm{cos}\theta \frac{m_pm_\pi }{q}$$
(A4)
For a thermal gas of relativistic bosons $`q2.7T`$, and with $`T_{CMBR}2.7K`$, corresponding to an energy of $`2.3\times 10^4eV`$, inserting the pion and proton rest masses gives a cut-off energy of

$$E_p10^{20}eV$$
(A5)
This defines the Greisen Zatsepin Kuzmin (GZK) limit.

The mean free path for this reaction is given by

$$\lambda =\frac{1}{n\sigma }$$
(A6)
where $`\sigma `$ is the capture cross-section for the $`p\gamma _{CMBR}`$ reaction, and $`n`$ is the number density of CMBR photons. $`\sigma `$ varies from a maximum of $`500\mu \text{barn}`$ at the main $`\mathrm{\Delta }^+`$ resonance to a level of $`100\mu \text{barn}`$ at higher energy levels (see for example (Mücke et al, 1999)).

## Appendix B Black Body Radiation

Planck’s black-body distribution law can be derived in two steps. First, by considering the number of radiation modes that can be supported in a black-body cavity of unit volume, an expression for the mode density can be obtained

$$dN(\nu )=\frac{8\pi \nu ^2}{c^3}d\nu $$
(B1)
The energy density as a function of frequency can then be calculated by multiplying the mode density by the average energy per mode. In classical terms this would simply be $`kT`$, where $`k`$ is the Boltzmann factor and $`T`$ is the absolute temperature. However, the quantum nature of photons results in the probability distribution being skewed, so that the correct mean energy per mode must be calculated from the Boltzmann distribution, giving

$$\overline{E}=\frac{h\nu }{e^{h\nu /kT}1}$$
(B2)
Combining (B1 and (B2) we obtain the final form of Planck’s law for the energy density of black-body radiation as a function of frequency

$$dU(\nu )=\frac{8\pi h\nu ^3}{c^3}.\frac{d\nu }{e^{h\nu /kT}1}$$
(B3)
The total energy density per unit volume is then simply obtained by integrating (B3) to give

$$U=\frac{8\pi ^5k^4T^4}{15c^3h^3}$$
(B4)
The final expression for the photon number density as a function of frequency is merely (B3) divided by the photon energy, $`h\nu `$

$$dN(\nu )=\frac{8\pi \nu ^2}{c^3}.\frac{d\nu }{e^{h\nu /kT}1}$$
(B5)
The photon number density per unit volume is therefore given by

$`N`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{8\pi \nu ^2}{c^3}}.{\displaystyle \frac{d\nu }{e^{h\nu /kT}1}}`$ (B6)
$`=`$ $`{\displaystyle \frac{16\pi k^3T^3\zeta (3)}{c^3h^3}}`$ (B7)


format: text

text:
# Effects of Bulk Viscosity in Non-linear Bubble Dynamics

## Abstract

The non-linear bubble dynamics equations in a compressible liquid have been modified considering the effects of compressibility of both the liquid and the gas at the bubble interface. A new bubble boundary equation has been derived, which includes a new term resulted from the liquid bulk viscosity effects. The influence of this term has been numerically investigated considering the effects of water vapor and chemical reactions on the bubble evolution. The results clearly indicate that the new term has an important damping role at the collapse, so that its consideration decreases the amplitude of the bubble rebounds after the collapse. This damping feature is more remarkable for higher deriving pressures.

When a small isolated gas bubble, immersed in a liquid, experiences a high amplitude spherical sound field, it grows and contracts non-linearly. Description of the dynamics of such non-linear motion is an old challenging problem. The radial dynamics of the bubble in an incompressible liquid is described by the well-known incompressible Rayleigh-Plesset equation Rayleigh:1917 ; Noltingk:1951 . The extension of this equation to the bubble motion in a compressible liquid has been studied by many previous authors Herring:1941 ; Keller:1980 . The most complete existing description was presented by Prosperetti and Lezzi Prosperetti:1987 . They used a singular-perturbation method of the bubble-wall Mach number and derived a one-parameter family of equations describing the bubble motion in the first order approximation of compressibility. This family of equations are written as:

$`(1(\eta +1){\displaystyle \frac{\dot{R}}{C}})R\ddot{R}+{\displaystyle \frac{3}{2}}(1{\displaystyle \frac{1}{3}}(3\eta +1){\displaystyle \frac{\dot{R}}{C}})\dot{R}^2=`$
$`{\displaystyle \frac{R}{\rho C}}{\displaystyle \frac{d}{dt}}\left(P_lP_a\right)+(1+(1\eta ){\displaystyle \frac{\dot{R}}{C}})\left({\displaystyle \frac{P_lP_aP_0}{\rho }}\right),`$
where, $`R`$, $`C`$, $`P_0`$, $`P_a`$, and $`\rho `$ are bubble radius, liquid sound speed, ambient pressure, driving pressure, and density of the liquid, respectively. Also, $`\eta `$ is an arbitrary parameter. Equation (Effects of Bulk Viscosity in Non-linear Bubble Dynamics) must be supplemented by a boundary condition equation at the bubble interface to relate the liquid pressure, $`P_l`$, to the gas pressure inside the bubble. Like all previous authors, Prosperetti and Lezzi Prosperetti:1987  used the following incompressible equation for this purpose:

$$P_l=P_g4\mu \frac{\dot{R}}{R}\frac{2\sigma }{R},$$
(2)
where, $`P_g`$, $`\mu `$, and $`\sigma `$ are gas pressure at the bubble interface, liquid viscosity coefficient, and surface tension, respectively. Most of the previously obtained equations belong to this single parameter family of equations, corresponding to different values of $`\eta `$. Moreover, $`\eta =0`$ yields results in closest agreement with the numerical simulation of full partial differential equations Prosperetti:1987 .

In all previous works Rayleigh:1917 ; Noltingk:1951 ; Herring:1941 ; Keller:1980 ; Prosperetti:1987 , an important approximation has been used in the derivation of the bubble dynamics equations. That is the incompressibility assumption of the liquid motion at the bubble interface, which has been used in the derivation of Eq’n. (2). Note that, all of the effects of the liquid compressibility in all previous papers have been resulted from the liquid motion around the bubble, but not from the bubble boundary condition equation. In fact, all previous authors, on one hand took into account the compressibility of the liquid motion around the bubble, but on the other hand neglected its consideration at the bubble interface.

In this paper, we have modified the bubble dynamics equations considering the effects of the liquid compressibility at the bubble interface. We have derived a new bubble boundary equation instead of Eq’n. (2). The new equation has new terms resulted from the effects of bulk viscosity of the liquid and the gas.

To derive the compressible bubble boundary equation, the continuity equation and the radial component of the stress tensor under the spherical symmetric condition can be written as:

$$\frac{1}{\rho }\left[\frac{\rho }{t}+u\frac{\rho }{r}\right]=\frac{u}{r}\frac{2u}{r}=\mathrm{\Delta },$$
(3)
$$T_{rr}=p+(\mu _b\frac{2\mu }{3})\mathrm{\Delta }+2\mu \left(\frac{u}{r}\right).$$
(4)
where, $`\rho `$, $`u`$, $`p`$, and $`\mathrm{\Delta }`$ are density, velocity, pressure, and divergence of the velocity, respectively. Also, $`\mu _b`$ is the bulk viscosity coefficient and is defined by $`\mu _b=\lambda +2\mu /3`$, where $`\lambda `$ is second coefficient of viscosity White:1991 . Inserting $`u/r`$ from Eq’n. (3), into Eq’n. (4) yields:

$`T_{rr}=p+(\mu _b+{\displaystyle \frac{4\mu }{3}})\mathrm{}4{\displaystyle \frac{\mu u}{r}}.`$ (5)
The velocity divergence, $`\mathrm{}`$, can be written as:

$$\mathrm{}=\frac{1}{\rho }\frac{d\rho }{dt}=\frac{1}{\rho c^2}\frac{dp}{dt},$$
(6)
where, the sound speed, $`c`$, is defined as $`c^2=dp/d\rho `$. The boundary continuity requirement at the bubble interface is:

$`T_{rr}(liquid)_R=T_{rr}\left(gas\right)_R+2{\displaystyle \frac{\sigma }{R}}.`$ (7)
Applying Eq’n. (5) for the gas and the liquid parts of Eq’n. (7) leads to:

$`P_l+4{\displaystyle \frac{\mu \dot{R}}{R}}`$ $``$ $`\left(\mu _b+{\displaystyle \frac{4\mu }{3}}\right)\mathrm{}_l=P_g+4{\displaystyle \frac{\mu _g\dot{R}}{R}}`$ (8)
$``$ $`\left(\mu _{bg}+{\displaystyle \frac{4\mu _g}{3}}\right)\mathrm{}_g2{\displaystyle \frac{\sigma }{R}},`$
where, $`\mu _g`$ and $`\mu _{bg}`$ are the viscosity and the bulk viscosity coefficients of the gas at the bubble interface, respectively. Also, $`\mathrm{\Delta }_l`$ and $`\mathrm{\Delta }_g`$ are the divergence of velocity of the liquid and the gas, respectively. Substituting the divergence of velocity for the liquid and the gas from Eq’n (6) into Eq’n. (8) yields:

$`P_l+4{\displaystyle \frac{\mu \dot{R}}{R}}`$ $`+`$ $`\left({\displaystyle \frac{\mu _b}{\rho C^2}}+{\displaystyle \frac{4\mu }{3\rho C^2}}\right){\displaystyle \frac{dP_l}{dt}}=P_g+4{\displaystyle \frac{\mu _g\dot{R}}{R}}`$ (9)
$`+`$ $`\left({\displaystyle \frac{\mu _{bg}}{\rho _g}}+{\displaystyle \frac{4\mu _g}{3\rho _g}}\right){\displaystyle \frac{d\rho _g}{dt}}2{\displaystyle \frac{\sigma }{R}},`$
where, $`\rho _g`$ is the gas density at the bubble interface. Equation (9) represents the bubble boundary equation containing all effects of the compressibility and the viscosity of both the liquid and the gas. Comparison of Eq’ns. (2) and (9) indicates the existence of three new terms in Eq’n. (9) due to the liquid and the gas compressibility and viscosity effects. Here, we concentrate on the effects of the new term arising from the liquid compressibility. Therefore, we neglect the gas viscosity because of its smallness relative to the liquid viscosity as in previous works Rayleigh:1917 ; Noltingk:1951 ; Herring:1941 ; Keller:1980 ; Prosperetti:1987 . Under this circumstance, Eq’n. (9) becomes:

$$P_l+\left(\frac{\mu _b}{\rho C^2}+\frac{4\mu }{3\rho C^2}\right)\frac{dP_l}{dt}=P_g4\frac{\mu \dot{R}}{R}2\frac{\sigma }{R}.$$
(10)
It should be mentioned that, although the effects of compressibility consideration in Eq’n. (Effects of Bulk Viscosity in Non-linear Bubble Dynamics) are in the first order approximation, but these effects have been introduced completely in Eq’n. (10).

To close the mathematical analysis, the gas pressure evolution at the bubble interface, $`P_g`$, must be specified. In the most complete approach, it can be determined from the simultaneous solution of the conservation equations for the bubble interior and the bubble radius equations C.C.WU ; Moss:1994 ; Kondict:1995 ; Voung:1996 ; Yuan:1998 ; Xu:2003 . Also, heat conduction and mass exchange between the bubble and the surrounding liquid affect the bubble evolution. In addition, chemical reactions occurring in the high temperature conditions at the end of the collapse, change the bubble content Kamath:1993 ; Yasui:1997 . All these complexities have been considered in a complete gas dynamics model by Storey and Szeri Storey:2000 .

On the other hand, strong spatial inhomogeneities inside the bubble are not remarkably revealed, unless at the end of an intense collapse Voung:1996 ; Yuan:1998 . Therefore, the uniformity assumption for the bubble interior seems to be useful and provides many features of the bubble motion Barber:1997 ; Brenner:2002 . Using this assumption, recently, Lohse and his coworkers presented an ODE model Toegel:2002 ; Toegel:2003 ; Lu:2003 , in which all effects of heat transfer at the bubble interface, phase change of water vapor, chemical reactions, and diffusion of reaction products have been considered. This model accurately describes various experimental phase diagrams Toegel:2003  and provides a good agreement with the complete direct numerical simulation of Storey and Szeri Toegel:2002 ; Storey:2000 .

Here, for describing the bubble interior evolution, we have used the Lohse’s group model (the same as what has been presented in Ref. Lu:2003 ). We do not repeat this model here and for more details we refer to Refs. Toegel:2003 ; Lu:2003 . The calculations were carried out under the framework of Eq’n (Effects of Bulk Viscosity in Non-linear Bubble Dynamics), ($`\eta =0`$), for both the new compressible (Eq’n. 10) and the old incompressible (Eq’n. 2) boundary conditions. We describe an argon bubble in water at room temperature, $`T_0=293.0K`$, and atmosphere pressure, $`P_0=1.0atm`$, under the conditions of Single Bubble Sonoluminescence Brenner:2002 ; Barber:1997 . The driving pressure was $`P_a(t)=P_a\mathrm{sin}\left(\omega t\right)`$, where $`\omega =2\pi \times 26.5kHz`$. The constants and the parameters were set accordingly CRC:1991 ; $`\rho =998.0kg/m^3`$, $`C=1483.0m/s`$, $`\mu =1.01\times 10^3kg/ms`$, $`\sigma =0.0707kgs^2`$. The bulk viscosity of water at room temperature was set to be $`\mu _b=4.1\times 10^3kg/ms`$ Karim:1952 . The constants and parameters of the gas evolution model were set the same as what has been presented in Ref. Lu:2003 .

Figures (1) and (2) illustrate the variations of the bubble characteristics (radius, total number of particle species, and temperature), for the two boundary condition cases. It is observed that the addition of the new viscous term in Eq’n. (10) considerably changes the bubble evolution after the collapse. The bubble motion is remarkably compressible during the collapse. Therefore, the new viscous term, which has been arisen from the liquid compressibility, is important in this time interval. This term exhibits a damping role and its consideration reduces the amplitude of the bubble rebounds. Also, the period of the rebounds decreases with the addition of the new term. Details of our calculations show that the minimum radius for the new case is about 10% greater than that of the old one. The difference between the two cases also appears on the variations of the total number of particle species (Ar and H2O plus reactions products) after the collapse (Fig. 1c). Note that, the difference gradually disappears as the bubble rebounds weaken.

In Fig. (2), details of the gas temperature evolution around the minimum radius have been demonstrated. Damping feature of the new term is clearly observed by a considerable decrease of the peak temperature (about 50%) and an increase of the temperature pulse width, at the collapse time. Also, the time of the peak temperature about 5 ns changes with the addition of the new term.

Figure (3) represents the effects of variations of $`P_a`$ on the bubble characteristics at the end of the collapse (peak temperature, mole fraction of H2O and reaction products, and minimum radius), for the two boundary condition cases. The ambient radius was fixed ($`R_0=5.0\mu m`$). Figure 3(a) shows that the peak temperature in the both cases grows as the driving pressure is increased. However, the rate of increase of the peak temperature for the new case is considerably smaller than that of the old one. This causes that the difference between the two cases becomes remarkable for the higher driving pressures (about 50% for $`P_a=1.5atm`$).

The bubble temperature at the end of the collapse is high enough to destroy the chemical bonds of water vapor molecules inside the bubble. The products of the dissociation of water vapor molecules are mainly H2, OH, H, O, and O2. The chemical reactions between the particle species existing inside the bubble affect the bubble content at the collapse and its peak temperature Kamath:1993 ; Yasui:1997 ; Storey:2000 ; Brenner:2002 ; Toegel:2002 ; Toegel:2003 ; Lu:2003 . Here, we have considered the effects of the reactions No. 1-8 of Refs. Yasui:1997 ; Toegel:2003 . The dependence of the mole fraction of H2O plus reactions products, which is defined as ($`N_{tot}N_{Ar}`$)/$`N_{tot}`$, to the driving pressure for the two boundary condition cases has been illustrated in Fig. 3(b). It shows that the mole fraction of H2O plus reactions products is similar for the two cases in low amplitudes. The difference appears for the higher deriving pressures. It is seen that, for the higher deriving pressures, the effect of H2O and reactions products is more important in the new case relative to the old one.

Figure 3(c) shows the variations of the minimum radius as a function of $`P_a`$, for the two cases. The trend of variations is similar for the two cases. But, the minimum radius for the new equation is more than that of the old equation because of the reduction of the collapse intensity.

A major deficiency of the old bubble dynamics equations is that for strongly driven bubbles, such as sonoluminescence bubbles, large amplitude rebounds are produced after the collapse, so that they often last until the next acoustic cycle of the periodic driving pressure. This is in contrast with the experimental results, which show rapidly damped rebounds Moss:2000 . By introducing a damping term arisen from the gas compressibility, Moss et. al Moss:2000  provided a typical solution for this problem. The effects of the suggested term by Moss et. al is very similar to the damping effects of the new term in this paper, (compare Fig. 1(b) with Figs. (3) and (4) of Ref. Moss:2000 ). It seems that the damping feature of the bulk viscosity is a better way for solving the mentioned problem. The reason is that Eq’n. (10) has been derived directly from the basic equations of fluid mechanics, on the contrary to Eq’n. (3.2) of Ref. Moss:2000 , which was derived by an approximate method.

According to the results of this paper, it is expected that the theoretical predictions of the bubble stability limits are affected by the addition of the new term to the bubble dynamics equations.

This work was supported by Sharif University of Technology and Bonab Research Center. Partial support of this work by Institute for Studies in Theoretical Physics and Mathematics is appreciated. The authors thank Andrea Prosperetti for his helpful comments.


text:
# Nonlinear optical beams carrying phase dislocations

## 1 Introduction

Phase dislocations carried by the wavefront of a light beam are associated with the zero-intensity points where the light intensity vanishes. The phase of the wave is twisted around such points creating a structure associated with an optical vortex. Optical beams with phase dislocations play an important role in linear singular optics .

In self-focusing nonlinear media, an intense finite-extent laser beam becomes localized due to the self-trapping mechanism which can compensate for the beam diffraction. The nonlinear self-action of light may result in the formation of stationary structures with both intensity and phase remaining unchanged along the propagation direction. Such self-trapped stationary structures of light beams are termed spatial optical solitons . When such solitons have phase singularities, they determine the internal structure of the beam; they can be stabilized by the light self-action generating nonlinear self-trapped optical beams carrying phase dislocations. Examples of such beams include vortex solitons  with point screw dislocations, multipole vector solitons  with $`\pi `$-edge dislocations, and more complicated ‘necklace’-type beams  and soliton clusters  with a combination of a screw dislocation at the beam center and, generally, $`\vartheta `$-edge dislocations, where $`\vartheta `$ is the phase jump between neighboring peaks in the intensity distribution .

The fundamental optical solitons show a fascinating combination of the properties of classical wave-packets together with a number of particle-like properties demonstrated in their elastic and inelastic interactions and mutual scattering, when each of the solitons preserves its identity. Moreover, the coherent interaction between the solitons depends strongly on a relative phase what provides an additional degree of freedom to control the interaction. We may draw an analogy between spatial soliton and the ‘atom of light’, and then the soliton collision and interaction can be treated in terms of the effective forces acting between these effective ‘atoms’. Following this concept, the higher-order multi-hump optical beams can be regarded as bound states of ‘atoms’ trapped by a common potential induced in a nonlinear medium. A balance of the interaction forces acting between the solitons is the necessary condition for the formation of the soliton ‘clusters’ or ‘molecules of light’.

In this paper we investigate the excitation of higher-order beams, including optical vortices and soliton clusters, through the inelastic soliton scattering and mutual trapping of initially well-separated fundamental solitons, the effect resembling a synthesis of ‘soliton molecules’. In addition, we propose the application of this effect in the context of ‘soliton algebra’  regarding the fundamental spatial solitons as the information carriers, and the transformation of an optical pattern induced by the soliton interaction as all-optical soliton switching.

## 2 Optical vortex solitons and soliton clusters

Optical vortices were introduced as the first example of a stationary light beam with the phase twisted around its core; the twist is proportional to $`2\pi `$ with integer $`m`$, the so-called topological charge of the phase dislocation . The physical model analyzing the evolution of the slowly varying field envelope $`E`$ is described by the nonlinear Schrödinger equation,

$$i\frac{E}{z}+\mathrm{\Delta }_{}E+f(I)E=0,$$
(1)
where $`\mathrm{\Delta }_{}`$ is the transverse Laplacian, and $`z`$ is the propagation distance measured in the units of the diffraction length. Function $`f(I)`$ describes the nonlinear properties of an optical medium, and it is assumed to depend on the total beam intensity, $`I=\left|E\right|^2`$. The simplest spatially localized solution of equation (1) carrying a phase dislocation, i.e. the vortex soliton, can be written in the form: $`E(r,\theta ,z)=A(r)\mathrm{}^{im\theta +i\beta z}`$, where $`A(r)`$ and $`\beta `$ are the beam amplitude and propagation constant, respectively, while $`r`$ and $`\theta `$ are the polar coordinates in the transverse plane.

In self-focusing nonlinear media, such vortex beams are the subject of the azimuthal modulational instability which result in splitting of the doughnut ring-like structure into a certain number of the fundamental solitons. The number of splitters and their dynamics are determined by the topological charge of the phase dislocation corresponding to the beam angular momentum (see, e.g.,  and references therein).

The simplest higher-order scalar stationary solutions found for the model (1) are a family of the radially symmetric solitons which includes radial modes with nodes in form of concentric rings. The important characteristic of these states is the soliton spin determined as a ratio of two conserved quantities, the beam angular momentum and the beam power. For the vortex soliton, the spin coincides with the topological charge $`m`$ of the phase dislocation carried by the vortex. Because of the condition of the field univocally, the topological charge $`m`$ is quantized and has integer value. The fundamental spatial solitons and their higher-order radial states have zero spin.

Novel types of the higher-order self-trapped optical beams can be introduced as azimuthally modulated self-trapped structures in the form of the so-called ‘necklace’ beams . However, it was found that a combination of the edge-phase dislocation with $`\pi `$-out-of-phase neighbor peaks cannot produce a stationary state, and the structure becomes slowly expanding . Such a stabilization is indeed possible for a more complicated system including the attraction between several incoherent beams . Another approach to this problem is to combine the screw dislocation in the origin of a ring-shaped beam with the edge dislocation within the necklace . The screw dislocation introduces a centrifugal force to the ring, being also responsible for spiralling and mutual repulsion of the solitons in the case of vortex break-up , and the edge dislocations prevents noise-induced instability break-up of the ring. Because of a nonzero angular momentum, the whole structure rotates with its propagation. As a result, the stabilization of the ring-shaped multihump beams requires complex phase distribution characterized by fractional value of the soliton spin .

A phase distribution required for the formation of quasi-stationary higher-order self-trapped optical beams was found in references  where the concept of soliton cluster was introduced. In this approach, the azimuthally modulated beam is regarded as a bound state of the interacting fundamental solitons. Because of phase-sensitive interaction, the requirement of the balance of the interaction forces between the solitons determines the beam phase in the form of a staircase-like screw dislocation. Figure 1 compares the vortex phase dislocation (left column) with the phase of a four-soliton cluster, having well defined $`\pi /2`$ steps between the soliton positions (right column). It was found  that a radially stable dynamical bound state is formed if these phase-jumps satisfy the condition $`\vartheta =2\pi m/N`$, with $`N4m`$ being the number of solitons in the ring.

Stability of the soliton clusters has been tested numerically for different nonlinear media, including cubic saturating, competing cubic self-focusing and quintic self-defocusing, and competing quadratic and cubic self-defocusing nonlinearities . The idea has been also extended to higher dimensions, covering the case of the spatio-temporal vortex solitons and light bullets. The common outcome of these studies is the confirmed robustness of the soliton clusters to random noise and strong radial perturbations. In the latter case, the pulsating states viewed as the radial excitations of ‘soliton molecule’ have been observed. Nevertheless, soliton clusters are the subject of modulational instability, and they are unstable with respect to azimuthal perturbations. The remarkable feature of this instability is that the number of the fundamental solitons flying off the main ring after the splitting is determined mainly by the topological charge $`m`$ instead of the initial number $`N`$ of solitons, similar to the vortex solitons. For what follows, we stress the fact that the conservation of the angular momentum of an optical beam in an isotropic medium determines the dynamics of splitters after break-up, so that the initial ‘spin’ angular momentum of vortex or cluster can be viewed as being transformed into the orbital momentum of the spiralling splitters .

## 3 Soliton molecules

Recent progress in both theoretical and experimental studies of the higher-order optical spatial solitons brought the soliton community to the gates of the direct search for all-optical soliton-based switching schemes, when the initial data carried by the light distribution on the front-face of the nonlinear medium can be processed, in a predictable and controllable way, by employing the light self-action effects. One of the examples of such an approach is the recently proposed concept of the ‘soliton algebra’ , based on the instability-induced break-up of optical vortices to the controllable number of the fundamental solitons. This idea also represents the example of a nontrivial approach to the soliton instability, when the symmetry-breaking instability, usually regarded as a serious disadvantage in using spatial solitons, is employed as a key physical mechanism for all-optical soliton switching from a given initial state (optical vortex) to the known final state defined by a certain number of fundamental solitons. This approach can be generalized to a broad variety of scalar and vector higher-order metastable solitons.

The symmetry-breaking soliton instability may serve as a physical mechanism for all-optical switching with only one disadvantage – it is a one-way process describing the transition from an initial higher-order state (a soliton molecule or cluster) to a number of simple stable states (dipole-mode and fundamental solitons). Below, we propose the opposite process, viewed as the excitation or ‘synthesis’ of higher-order states from predefined number of initially separated solitons, or ‘atoms of light’, in a nonlinear bulk medium. Indeed, introducing molecules of light would not be self-consistent without the possibility of mutual trapping of the free atoms or molecule synthesis.

To demonstrate this phenomenon, we propagate numerically the ring-shaped arrays of initially well separated coherently interacting fundamental solitons in a saturable nonlinear medium. Figure 2 shows a characteristic example of a set of six solitons which have their relative phases growing in steps of $`\vartheta =\pi /3`$ along the ring, being initially directed to collide with each other. This initial condition corresponds to the inversion of the instability-induced ring break-up, so that the solitons move towards the ring instead of flying away. We observe highly inelastic collision of the solitons when they strongly interact overlapping and loosing their identity. Nevertheless, the initial phases of the solitons are tilted in such a way that the total phase of the beams forms a screw dislocation in the ring origin which prevents a simple fusion of all solitons. Instead, the ring-shaped structure is formed. Similar situation is observed in figure 3 with an array of eight solitons and the formation of a metastable vortex ring. Due to large amplitude modulations, these intermediate (or metastable) structures never form ideal stationary states, and they quickly split off to a set of new isolated solitons, with the total number of splitters predefined by the initial conditions. In this way, we were able to produce, as a final state, the patterns with different number of solitons by changing initial parameters, including the number and phase tilt of solitons. Generally, the final number of solitons is determined by the ring instability mode with a largest growth rate, and in a saturable medium it is usually twice the topological charge . At the same time, we were able to force single-charged ‘meta-vortex’ in figures 2 and 3 to split to three and four fundamental solitons, respectively. The intermediate meta-state shows complex dynamics of the instability development which may continue for several tenths of diffraction lengths. Thus, the whole picture of the soliton collision, ring formation, and the ring splitting is somewhat similar to the ‘delayed-action interaction’ recently reported for interacting composite solitons carrying nonzero angular momentum . We note also that, in addition to the known transformation of ‘spin-to-orbital’ angular momentum , in figures 2 and 3 we observe a kind of ‘orbital-spin-orbital’ transformation. The change of the field momentum shown in figure 2 is about $`10\%`$ of the initial value, and it occurs only at the break-up stage.

Highly unstable ring-shaped beams are formed as a result of mutual trapping and inelastic interaction of coherent solitons, as shown in figures 2 and 3. They represent the intermediate steps in the process of the nonlinear switching between the states with a different number of solitons. In the context of the soliton algebra and all-optical switching, it might be of interest to study in more details and to determine the quantitative parameters necessary to obtain the final state with a given number of solitons. At the same time, the mutual trapping of coherent solitons is found to be too sensitive to the initial perturbations in order to produce the metastable clusters.

Nevertheless, it is indeed possible to generate higher-order optical beams, i.e. to synthesize the soliton molecules, from interacting fundamental coherent and incoherent solitons. In figure 4 we show the excitation of a long-lived four-lobe vector cluster. We note that the four solitons in the figure 4(b) at $`z=0`$ are directed exactly to the center so that there is no azimuthal tilt of their phase. As within the cluster in figure 1, each soliton has a phase growing by a jump $`\vartheta =\pi /2`$ along the ring, and this initial phase distribution guarantees the appearance of the screw phase dislocation, even the corresponding value of the total angular momentum is very small (a solid line in figure 4). Because only the total angular momentum is conserved, but not the partial momenta of the components, there is a freedom for components to symmetrically exchange the angular momentum during the beam propagation. After mutual trapping of all solitons at $`z25`$ the new-born vector cluster experiences strong radial oscillations between the states with maximal (e.g. at $`z=30`$ and $`z=60`$) and minimal partial angular momenta in the components, see the diagram in figure 4. At the same time, the $`\pi /2`$-phase jumps introduced initially between the solitons (edge dislocations) survive these strong oscillations, and the whole cluster preserves its structure for a distance exceeding 100 diffraction lengths. Therefore, we observe the formation of a composite state, the soliton molecule, by colliding simple solitons with a nontrivial phase pattern.

## 4 Conclusions

We have studied the scattering and mutual trapping of several fundamental solitons and the generation of soliton clusters and vortex solitons – the complex self-trapped states of light carrying phase dislocations in the wave front. Inelastic collision of solitons has been shown to result in the formation of ring-shaped beams or metastable vortices which subsequently break-up creating different number of fundamental solitons defined by the initial conditions. This kind of the soliton delayed-action interaction and the nonlinear transformation of the number of fundamental solitons has been analyzed in the context of the soliton algebra and all-optical soliton switching. We have also shown that the vectorial interaction between the field components provides an additional mechanism of the soliton cluster stabilization. In particular, we have demonstrated the formation of vector soliton clusters from colliding solitons, the process which can be regarded as a synthesis of ‘soliton molecules’.

A support from the Alexander-von-Humboldt Foundation and the Australian Research Council is acknowledged.

## References


text:
# VLBA Calibrator Survey: Astrometric and Image Results

## 1. Introduction

The purpose of the VLBA Calibrator Survey (Beasley et al. 2002), (Fomalont et al. 2003) was to extend the list of compact sources with coordinates known at the milliarcseconds level for use as calibrators for mapping and as target sources for geodetic applications. In total, the positions of 1608 new sources were obtained. As a result of processing VCS sessions the total number of sources with precise coordinates increased by a factor of three. The goals of the survey were: a) to increase the surface density of known geodetic-grade calibrators with mas-accurate positions; b) to facilitate routine phase-referencing allowing high-resolution radio imaging of weak target sources; c) to provide a uniform image database at 2.3 and 8.4 GHz for use in scientific applications.

## 2. Observations

Observations were carried out in twelve 24 hour sessions. Observations used the VLBA dual-frequency geodetic mode observing simultaneously at 2.3 and 8.4 GHz. Signal was recorded at four 8 MHz wide baseband channels over a large spanned bandwidth (100 MHz at S-band and 400 MHz at X-band). Each source was observed 2–3 times for 60–90 sec.

## 3. Astrometric catalogue

Before analysis of the VCS sessions the number of sources observed under geodesy/astrometry programs with positions known better than 1 nrad was 615. After analysis of the VCS campaign this number increased to 2011. Now there is at least one calibrator within $`4^{}`$ of any target direction at 90.1% of the sky at declinations $`\delta >40^{}`$, as shown in figure 1.

Some 15% of sources from the VCS lists were observed previously in astrometric/geodetic programs. Analysis of the differences between the VCS catalogues and Goddard gsf2003c catalogue allowed us to assess the error model. As a result of this analysis the formal uncertainties were inflated in this way:

$$\sigma (\alpha )_r=\sqrt{(1.5\sigma _\alpha )^2+E(\alpha )^2/cos^2\delta }\sigma (\delta )_r=\sqrt{(1.5\sigma _\delta )^2+E(\delta )^2}$$
(1)
$`E(\alpha )`$ and $`E(\delta )`$ are declination dependent:

## 4. Image results

Amplitude and initial phase calibration was done by AIPS software. Sources were imaged by the Caltech DIFMAP package. Some 70% of sources were imaged automatically, the others were imaged manually. The typical image rms is 2–3 mJy, with dynamic range 30:1 or better.

The image database and a source search engine are available on the Web at http://magnolia.vlba.nrao.edu/vlba\_calib . The astrometric catalogue is in http://gemini.gsfc.nasa.gov/solutions/astro .

## References

Beasley, A.J., Gordon, D., Peck, A.B., Petrov, L., MacMillan, D.S., Fomalont E.B., Ma, C., 2002, ApJS, 141, 13.

Fomalont E., Petrov L., McMillan D.S., Gordon D., Ma C., 2003, submitted to AJ.