text:
# Fundamental thermal fluctuations in microspheres

## I  Introduction

The development of modern technology in many fields leads to further miniaturization of components. This makes it necessary to take into account certain fundamental physical limitations. An example of such limitations are thermodynamical fluctuations of temperature in a small volume. These fluctuations are transformed to wideband noise in the output channel due to the temperature dependence of device parameters. The same limitations also frequently appear in experimental physics in macroscopic high–precision measurements. As it was recently shown, fundamental fluctuations with the same origin limit the sensitivity of gravitational wave antenna (international project LIGOA.Abramovici et al. (1992)), where thermal expansion, thermal dependence of refractive index and Young’s modulus give rise to different types of noise V.B.Braginsky et al. (1999, 2000); G.Cagnoli and P.A.Willems (2002). Thermorefractive fluctuations lead to phase noise in long fibersW.H.Glenn (1989); K.H.Wanser (1992), which were used to observe the effect for the first time K.H.Wanser et al. (1993).

Microspheres V.B.Braginsky et al. (1989) are relatively novel type of optical resonators uniquely combining small size (from tens to thousands of micrometers) and high quality–factor up to $`Q10^{10}`$ for the so–called “whispering–gallery modes”M.L.Gorodetsky et al. (1996) (WGMs). The small size of the effective volume, occupied by the e.m. field of the mode, makes for low thresholds of nonlinear effects such as bistability and oscillatory instability M.L.Gorodetsky and V.S.Ilchenko (1992), which are preconditioned by Kerr’s and thermal effects. The small volume makes it possible for such resonators to be used as a tool for detection and measurement of the thermorefractive noise M.L.Gorodetsky (2001). Apart from its importance in the LIGO project, the measurement of thermorefractive noise can serve as an innovative experimental examination of the theory of microscopic fluctuations of temperature. Thermorefractive noise should also be taken into account in possible applications of whispering–gallery modes resonators such as diode laser stabilization.

## II Thermorefractive noise in microresonators

One can understand the effect from the well known thermodynamical equation for the variance of temperature fluctuations $`u`$ in the volume $`V`$:

$$u^2=\frac{\kappa T^2}{\rho CV},$$
(1)
where $`T`$ is the temperature of the heat–bath, $`\kappa `$ is the Boltzmann constant, $`\rho `$ is density, and $`C`$ is specific heat capacity. By substituting in parameters for fused silica: $`\rho =2.2\text{g}/\text{cm}^3`$, $`C=6.7\times 10^6\text{erg}/(\text{g}\text{K})`$ and effective volume of field of the most localized WGM in the microsphere of radius $`R50\mu \text{m}`$ with $`V_{eff}10^9\text{cm}^3`$, we obtain the value of the standard deviation of temperature $`\sqrt{u^2}30\mu \text{K}`$. These temperature fluctuations, combined with the coefficient of thermal refraction $`dn/dT=1.45\times 10^5K^1`$ lead to the effect of relative eigen frequency fluctuations $`\delta \omega /\omega 3\times 10^{10}`$. The estimate is comparable with the bandwidth of resonances achievable in microspheres. We do not consider thermoelastic noise in microspheres, since the coefficient of thermal expansion $`\alpha =5.5\times 10^7\text{K}^1`$ is sufficiently smaller in fused silica than $`dn/dT`$. In subsequent sections we will present more rigorous analysis aimed at finding spectral properties of this noise, taking into account peculiar field distribution of WGM.

Variations of the refractive index $`n`$ in the dielectric cavity perturb its resonant frequencies. The perturbed wave–equation has the form:

$$\mathrm{\Delta }\stackrel{}{E}+(ϵ^0+2n\delta n)\frac{\omega ^2}{c^2}\stackrel{}{E}=0,$$
(2)
where $`\stackrel{}{E}`$ is the electric field strength in the cavity, $`ϵ^0=n^2`$ is permittivity, and $`\delta n=\frac{dn}{dT}u`$ is the variation of refractive index due to fluctuations of temperature $`u`$. If $`\stackrel{}{E}_0`$ is the ortho–normalized field distribution of an eigenmode of the unperturbed cavity ($`\stackrel{}{E}_i\stackrel{}{E}_j^{}𝑑\stackrel{}{r}=\delta _{ij}`$) and $`\omega =\omega _0+\delta \omega `$ is the frequency shift then after multiplication of this vector equation on complex conjugated vector $`\stackrel{}{E}_0^{}`$ and integration over the whole volume, neglecting the second order terms we obtain:

$$\frac{\delta \omega }{\omega _0}=\frac{1}{n}\underset{V}{}|\stackrel{}{E}_0^2|\delta n𝑑\stackrel{}{r}=\frac{1}{n}\frac{dn}{dT}\overline{u},$$
(3)
where $`\overline{u}`$ is the temperature deviation, averaged over the mode volume.

## III Power spectral density of thermal fluctuations

To calculate the effect of fluctuations of temperature the method of fluctuational thermal sources $`F(\stackrel{}{r},t)`$ may be used V.B.Braginsky et al. (1999, 2000):

$`{\displaystyle \frac{u}{t}}D\mathrm{\Delta }u=F(\stackrel{}{r},t),`$ (4)
where $`D=\lambda ^{}/(\rho C)`$ is thermal diffusivity and $`\lambda ^{}`$ is thermal conductivity ($`\lambda ^{}=1.4\times 10^5\text{erg}/(\text{cm}\text{s}\text{K})`$, $`D=9.5\times 10^3\text{cm}^2/\text{s}`$ for fused silica). This approach is analogous to the Langevin approach, which uses fluctuational forces in the equations of dynamics. It was shown before that if proper normalization of the sources is used:

$`B_{rt}^F`$ $`=`$ $`F(\stackrel{}{r},t)F(\stackrel{}{r}^{},t^{})`$ (5)
$`=`$ $`{\displaystyle \frac{2\kappa T^2D}{\rho C}}^2\delta (\stackrel{}{r}\stackrel{}{r}^{})\delta (tt^{}),`$
this approach leads to the correct results, which satisfy the Fluctuation–Dissipation theorem (FDT). In particular, it was shown that thermoelastic noise is associated through FDT with thermoelastic damping. It is also possible to show that thermorefractive noise is connected through FDT with electrocaloric losses M.L.Gorodetsky et al. .

Thermodynamical fluctuations of temperature, averaged over the mode volume, may be calculated as:

$$\overline{u}=u(\stackrel{}{r},t)|E_0(\stackrel{}{r})|^2𝑑\stackrel{}{r},$$
(6)
where $`|\stackrel{}{E}_0|^2𝑑\stackrel{}{r}=1`$ and

$$u(\stackrel{}{r},t)=\frac{F(\stackrel{}{\beta },\mathrm{\Omega })}{D\beta ^2+i\mathrm{\Omega }}e^{i\mathrm{\Omega }t+i\stackrel{}{\beta }\stackrel{}{r}}\frac{d\mathrm{\Omega }d\stackrel{}{\beta }}{(2\pi )^4},$$
(7)
is the general solution of (4). Spectral correlations of fluctuational forces satisfy the following condition:

$`B_{\beta \mathrm{\Omega }}^F`$ $`=`$ $`F(\stackrel{}{\beta }^{},\mathrm{\Omega }^{})F^{}(\stackrel{}{\beta },\mathrm{\Omega })`$ (8)
$`=`$ $`(2\pi )^4{\displaystyle \frac{2kT^2D}{\rho C}}\beta ^2\delta (\stackrel{}{\beta }\stackrel{}{\beta }^{})\delta (\mathrm{\Omega }\mathrm{\Omega }^{}).`$
We may now calculate the following averaged value $`B_\tau ^u=\overline{u}(t)\overline{u}(t+\tau )`$ (correlation function of temperature fluctuations averaged over the mode volume) and so, from Wiener–Hinchin, theorem the one–sided (hence additional factor 2) power spectral density $`S_{\overline{u}}(\mathrm{\Omega })`$ of fluctuations of temperature:

$`S_{\overline{u}}(\mathrm{\Omega })={\displaystyle \frac{4kT^2D}{\rho C}}`$ (9)
$`\times `$ $`{\displaystyle \frac{\beta ^2|E_0(\stackrel{}{r})|^2|E_0(\stackrel{}{r}^{})|^2}{D^2\beta ^4+\mathrm{\Omega }^2}e^{i\stackrel{}{\beta }(\stackrel{}{r}\stackrel{}{r}^{})}𝑑\stackrel{}{r}𝑑\stackrel{}{r}^{}\frac{d\stackrel{}{\beta }}{(2\pi )^3}}`$
$`=`$ $`{\displaystyle \frac{4kT^2D}{\rho C}}{\displaystyle \frac{\beta ^2|G(\stackrel{}{\beta })|^2}{D^2\beta ^4+\mathrm{\Omega }^2}\frac{d\stackrel{}{\beta }}{(2\pi )^3}},`$
where

$$G(\stackrel{}{\beta })=|\stackrel{}{E}_0|^2e^{i\stackrel{}{\beta }\stackrel{}{r}}𝑑\stackrel{}{r}$$
(10)
is normalized spatial spectrum of the energy distribution in the resonator. To verify this useful expression we may integrate it over all frequencies:

$`u^2={\displaystyle \underset{0}{\overset{\mathrm{}}{}}}S_{\overline{u}}(\mathrm{\Omega }){\displaystyle \frac{d\mathrm{\Omega }}{2\pi }}={\displaystyle \frac{kT^2}{\rho C}}{\displaystyle |G(\stackrel{}{\beta })|^2\frac{d\stackrel{}{\beta }}{(2\pi )^3}}.`$ (11)
Comparing the final expression with (1) we deduce that

$$V_{eff}^1=|G(\stackrel{}{\beta })|^2\frac{d\stackrel{}{\beta }}{(2\pi )^3}=|\stackrel{}{E}(\stackrel{}{r})|^4𝑑\stackrel{}{r}.$$
(12)
The same general expression for the effective volume of the mode in microresonator appears in the analysis of nonlinearity M.L.Gorodetsky and V.S.Ilchenko (1992) and scattering M.L.Gorodetsky et al. (2000) in microspheres.

It is important to note that the above expressions were obtained ignoring boundary conditions. If the field, is concentrated near thermally isolated surface, as in the case of microsphere, (8) should be modified, substituting in (8) $`[\delta (\beta _{}\beta _{}^{})+\delta (\beta _{}+\beta _{}^{})]\delta (\stackrel{}{\beta }_{}\stackrel{}{\beta }_{}^{})`$ instead of $`\delta (\stackrel{}{\beta }\stackrel{}{\beta }^{})`$, where $`\stackrel{}{\beta }_{}`$ is the component of wave–vector of fluctuations normal to the surface and $`\stackrel{}{\beta }_{}`$ are components, parallel to it.

In the analysis above, the medium was considered infinite (the volume of field localization is significantly smaller than the size of a device). Discrete spectrum of thermal waves should be considered for more accurate calculations, especially at very low frequencies:

$$F(\stackrel{}{r},t)=F_\nu (\mathrm{\Omega })\mathrm{\Phi }_\nu (\stackrel{}{r})e^{i\mathrm{\Omega }t}\frac{d\mathrm{\Omega }}{2\pi }.$$
(13)
$`B_{\nu \mathrm{\Omega }}^F`$ $`=`$ $`F_\nu (\mathrm{\Omega })F_\nu ^{}(\mathrm{\Omega }^{})`$ (14)
$`=`$ $`{\displaystyle B_{rt}^F\mathrm{\Phi }_\nu (\stackrel{}{r}^{})\mathrm{\Phi }_\nu ^{}(\stackrel{}{r})e^{i(\mathrm{\Omega }^{}t^{}\mathrm{\Omega }t)}𝑑\stackrel{}{r}𝑑t𝑑\stackrel{}{r}^{}𝑑t^{}}`$
$`=`$ $`2\pi {\displaystyle \frac{2\kappa T^2D}{\rho C}}\beta _\nu ^2\delta (\mathrm{\Omega }\mathrm{\Omega }^{})\delta (\nu ,\nu ^{}).`$
Analogously to (7)

$$u(\stackrel{}{r},t)=\underset{\nu }{}\frac{F_\nu (\mathrm{\Omega })}{D\beta _\nu ^2+i\mathrm{\Omega }}\mathrm{\Phi }_\nu (\stackrel{}{r})e^{i\mathrm{\Omega }t}\frac{d\mathrm{\Omega }}{2\pi },$$
(15)
and

$`\overline{u}(t)={\displaystyle \underset{\nu }{}\frac{F_\nu (\mathrm{\Omega })G_\nu }{D\beta _\nu ^2+i\mathrm{\Omega }}e^{i\mathrm{\Omega }t}\frac{d\mathrm{\Omega }}{2\pi }},`$ (16)
where $`G_\nu `$ are the coefficients of field’s intensity decomposition into normal thermal waves $`\mathrm{\Phi }_\nu (\stackrel{}{r})`$ of the finite medium.

Now, as before, we may calculate the correlation function of relative frequency fluctuations and power spectral density:

$`S_{\overline{u}}(\mathrm{\Omega })={\displaystyle \frac{4kT^2D}{\rho C}}{\displaystyle \underset{\nu }{}}{\displaystyle \frac{\beta _\nu ^2|G_\nu |^2}{D^2\beta _\nu ^4+\mathrm{\Omega }^2}}.`$ (17)
Applying the obtained expressions to the microsphere leads to the following expression for the power spectral density of relative frequency fluctuations:

$`S_{\delta \omega /\omega }(\mathrm{\Omega })`$ $`=`$ $`S_{\overline{u}}\left({\displaystyle \frac{dn}{dT}}\right)^2={\displaystyle \frac{\kappa T^2\sqrt{\mathrm{}}}{\pi ^{3/2}n^2R^2\sqrt{\lambda ^{}\rho C\mathrm{\Omega }}}}\left({\displaystyle \frac{dn}{dT}}\right)^2`$ (18)
$`\times `$ $`{\displaystyle \frac{1}{\sqrt{1b^2/d^2}}}{\displaystyle \frac{1}{(1+(\mathrm{\Omega }\tau _b)^{3/4})^2}},`$
where $`\mathrm{},m,b,d,\tau _b`$ are parameters, determined by the mode in microsphere (see Appendix A). This expression was experimentally verified in this work.

It is essential to note that in a finite body $`S_{\overline{u}}(\mathrm{\Omega })`$, contrary to (18), remains limited at zero frequency. However, instead of the sum (A), we used the former expression (18), obtained from the continuous spectrum, to evaluate the results of the measurements. There were important reasons for doing this. First of all, too many approximations were used when obtaining the sum. These are crude approximations especially for lower frequencies. However, it is shown in the Appendix and confirmed by numerical calculations that the sum (A) and integral (27) lead to the same asymptotic dependence at high frequencies. Secondly, microspheres during fabrication are formed on short fused silica stem which conducts heat and couples microspheres to heat–bath. Thermal exchange with the atmosphere due to convection is also not negligible. Therefore, microspheres are not thermally isolated. Thirdly, other parasitic noise effects dominate at low frequencies. Throughout the experiment, we only considered this phenomenon at frequencies of $`\mathrm{\Omega }>1/\tau _R=D/R^2400c^1`$ for ($`R=50\mu m`$).

## IV Experimental Setup

The idea of measurement of thermorefractive noise in optical microspheres is quite simple. If one tunes the measuring laser’s frequency on the slope of resonance curve of WGM, that is on the range, where amplitude of the signal being measured sharply depends on frequency, then the trembling of resonator’s eigen frequency transform into fluctuations of intensity of output radiation. These intensity fluctuations are recorded and further processed.

To measure thermorefractive noise we used small $`2R80÷260\mu m`$ fused silica microspheres with quality factors of $`Q10^9`$, manufactured with the hydrogen miniburner. In order to avoid degradation of Q–factor due to adsorption of atmospheric water M.L.Gorodetsky et al. (1996) microspheres were placed upon fabrication to the atmosphere of dry clean nitrogen into the special chamber where all measurements were conducted.

WGM in microspheres were excited with the prism coupler using He–Ne laser ($`\lambda =0.63\mu m`$) with piezo–driven front mirror as in M.L.Gorodetsky and V.S.Ilchenko (1994). Laser frequency could be tuned within the range of approximately $`0.8GHz`$. The output power of laser ($`540\mu W`$) was reduced to $`45\mu W`$ passing into the chamber. The reduction of power was necessary to weaken thermal and Kerr nonlinearities M.L.Gorodetsky and V.S.Ilchenko (1992) and was carried out with the neutral glass filters. Frequency of the laser was locked at the slope of resonance curve. This was done with the use of mini–spiral heater, placed nearby the resonator, which actively stabilized the temperature of the microsphere, hence stabilizing its eigen frequency. To provide the feedback loop a weak signal at $`91.2kHz`$ was used to modulate laser frequency. This frequency was chosen to coincide with one of acoustic resonances of the front mirror to facilitate laser frequency modulation. This signal is demodulated by the WGM and error signal was used to correct the current which feeds the heater. When the WGM is locked to the laser, the output of the microsphere comprises the thermorefractive amplitude noise. The output from the resonator was registered with Thorlabs PDA500 amplified GaAsP photodetector, with $`40dB`$ gain in amplitude and $`f_d=45kHz`$ bandwidth at this gain setting.

The signal from detector was digitized using ADC (ICP DAS PCI–$`1802L`$ computer board) with digitization speed of 333 kilosamples per second and 12–bits per sample resolution. Continuous sets of $`6553600`$ ($`19.66s`$) points were recorded in each session of measurement for further processing in each session of measurement. To find the estimate of noise frequency spectrum, the record was subdivided into many equal time intervals. Power spectral densities were calculated using Fast Fourier Transform (FFT) algorithm for each interval and then averaged over the entire set of intervals within each session. 100 intervals with $`65536=2^{16}`$ points in each were used to produce estimate of spectra with fine resolution of $`5Hz`$. To avoid the effect of frequency aliasing, the frequency band of the signal incoming from the detector was limited by the Butterworth filter of the eighth order. At stop–band frequency of $`166kHz`$ (the Nyquist frequency in our case) this filter attenuates the signal to -48dB level. First order RC–filter with time constant of $`\tau 2s`$ was used to get rid of the constant offset.

### IV.1 Calibration of spectra

In order to obtain absolute values of spectral densities for measured spectra, we have developed and used the method of calibration based on laser frequency modulation. In this method, weak sinusoidal voltage of known amplitude and frequency is admixed to piezo–actuator of laser’s mirror, thereby producing weak laser’s frequency modulation. Since the measured spectra represent relative frequency fluctuation, this additional sinusoidal modulation of laser frequency results in a narrow peak in frequency spectrum. If the amplitude of frequency modulation is known, it is possible to find the value of spectral density corresponding to the spike.

To find the coefficient of transformation of voltage to laser frequency modulation and to measure the quality factors of the resonators, we performed calibration of piezo–actuator using following technique. If laser frequency is swept by the saw–type voltage with known amplitude in the vicinity of eigen frequency of the microsphere and additional phase modulation with frequency higher than the bandwidth of resonant curves is applied to the pump, one can see on the detector two additional resonances apart from the central frequency. These sidebands produce frequency scale for the calibration of piezo–actuator and for the measurement of Q–factor from the width of resonant curves (see fig.2).

The coefficient of piezo–actuator has been measured and is equal to $`\frac{f}{v}=(3,63\pm 0.1)\frac{MHz}{V}`$.

If the amplitude of frequency modulation produced by piezo–actuator is known, the amplitude of the peak it produces in calculated power spectral densities can be used for their absolute calibration (see Appendix B).

To produce calibration peak we applied $`5÷15mV`$ modulation at $`172.9Hz`$ or $`1002Hz`$ (preferred frequency at later stages of experiment) to piezo–actuator of laser mirror, which resulted in relative frequency modulation of $`P=\delta \omega /\omega =(4÷11)\times 10^{11}`$. These frequencies were chosen to match that of bins of discrete Fourier transform of 65536 points (namely $`1001.99=197F_s/65536`$ and $`172.93=34F_s/65536`$). The higher frequency of $`1kHz`$ was chosen so that the calibration peak to be separated from low frequency setup noises in the frequency spectrum. It is important to note that the described calibration method does not require the quality factor of the microsphere as well as the transformation coefficients and gains in electrical tracts to be known.

After amplitude calibration is performed, additional correction is required to compensate for frequency filters in the detector by dividing the resulting spectrum by $`1/\sqrt{1+(f/f_d)^2}`$ — the frequency response of the detector. Frequency response of anti–aliasing filter was also compensated. To do this, it was measured and approximated with the polynomial of 7–th order so that the compensation procedure could be carried out by point–to–point multiplication.

### IV.2 Excess noises in the measurements

In this research, the superfluous noises were represented by the electrical interference (generally $`50Hz`$ and harmonics), acoustical and seismical noises, electromagnetic foil of devices and also by inherited amplitude and frequency noises of the laser. The frequency noise of the laser has turned out to be the most essential source of influence on the quantity being investigated. We have measured this noise together with all other technical noises of the setup. The same technique as for thermorefractive noise measurements was used with only difference of using relatively large microspheres with diameters of $`481,570,508,588,894,619\mu m`$ and modes with large $`lm`$ numbers of the order of 100 and above. In these modes, thermorefractive noise appears to be small enough for laser noise to dominate. Measurements of the laser noise for all aforementioned microspheres coincided well enough. The microspheres operated as the frequency discriminators in this case.

Laser noises calibration was carried out in the same way as it was done for thermorefractive noise calibration — with the use of calibrating peak. In measurements of the technical noise spectra calibrating spike at frequency of $`172Hz`$ was used.

Special measures were undertaken to weaken electromagnetic hindrances in the setup. We have applied proper signal shielding, moreover, we have shielded ADC board inside the computer, which resulted in the setup noise decrease of more than order of magnitude. The dynamic range of digitization system exceeded $`90dB`$ (taking into account averaging over many spectra).

### IV.3 Identification of WGMs and computer processing

Upon completion of recording process, the identification of WGM indices was carried out. The knowledge of the modes’ indices allows one to calculate the effective volume, occupied by the electromagnetic field of the mode, as well as theoretical energy distribution within the volume and therefore the theoretical spectral density. To perform the identification, filming with the digital camera of the speckled image of the mode was carried out. The speckled image is called forth by scattering of mode’s e.m. field at the residual molecular inhomogeneities on the surface of the microsphere M.L.Gorodetsky et al. (1996).

The approximation formulas were derived, linking the microsphere’s radius with $`\mathrm{}`$ index and the angular half–width of the WGM’s belt with $`\mathrm{}`$ and $`\mathrm{}m`$ numbers. (See Appendix C).

Upon recording, images were magnified and processed with the made–on–purpose computer program of graphical identification of the WGMs. The belt’s width ($`\mathrm{}m`$ value), simulating sphere size, slope and tilt of the belt could be controlled interactively, so that the theoretical framework projection of the mode could be made coinciding with the photographic projection of the real mode. This allowed one to obtain the parameters of the mode. The error of estimate of the $`\mathrm{}m`$ value is proportional to $`\sqrt{\mathrm{}m+1/2}`$ i.e. the width of the mode belt, but is strongly determined by the quality of the pictures obtained. Modes with $`\mathrm{}m3`$ can usually be precisely recognized while for $`3<\mathrm{}m<6`$ the error in obtaining the indices approaches 1. Index $`q`$ of the mode can be in principle found from the optimal coupling angle of incidence of the pump beam M.L.Gorodetsky et al. (1996). However, this possibility has not been realized in the present work. We just tried to excite and analyze the modes with smaller $`q`$ index which are usually characterized by higher Q and better coupling.

## V  Thermorefractive noise spectra

Figure 4 presents typical high resolution spectrum of measured noise in a microsphere 138$`\mu `$m in diameter and the spectrum of technical noises in large 894$`\mu m`$ sphere measured as described in IV.2. Dashed line shows theoretical curve for the spectrum calculated from recognized value of $`\mathrm{}m=4\pm 1`$. Good agreement of experimental data with theory for the frequencies $`>100Hz`$ is clearly observable. It can also be seen that the level of combined technical noises, which include frequency and power fluctuations of the laser and electrical noises in circuits is nearly $`20dB`$ lower than the observed effect, except for isolated spikes. Therefore, technical noises do not prevent the experimental test of thermorefractive model.

Calibrating spikes at $`1002Hz`$ for the smaller sphere and at $`172.93Hz`$ for the larger one can be seen on the graphics. At low frequencies, peaks of technical noise at AC $`50Hz`$–harmonics dominate. Additional sharp spikes could be seen at frequencies in the order of $`20kHz`$ which were identified as acoustical resonances of laser mirrors, responding to acoustical noise produced by laboratory equipment. One of these resonances at $`92kHz`$ was used in feedback loop for active mode frequency stabilization (see IV.1). These spikes, except the specially excited one at $`92kHz`$, are relatively small and narrow ($`20`$Hz) and could be seen only due to high resolution of our measurements ($`5Hz`$ in the range of $`10÷10^5Hz`$). They do not affect the result especially if additional averaging at high frequencies is applied.

Figure 5 depicts results of measurement of calibrated spectra of relative frequency fluctuations in four different microspheres for six different modes. The sizes of microspheres and estimated mode parameters are given on the figure.

Data for these graphics was obtained from the data presented on the previous figure using uniform in log–scale averaging according to the following algorithm: if $`f_{min}(f_{max}/f_{min})^{j/N}<f_i<f_{min}(f_{max}/f_{min})^{(j+1)/N}`$ then if $`n_j>0`$, $`f_j=\frac{1}{n_j}f_i`$; $`S_j=\frac{1}{n_j}S_i`$, where $`N`$ is the number of equidistant bins ($`N=200`$ for the presented figure) in log scale, $`n_j`$ is the number of points found in bin $`j=0..N1`$, $`f_{min}=10Hz`$ and $`f_{max}=10^5Hz`$. This approach allows one to estimate better the real spectra on higher frequencies while preserving frequency resolution at low frequencies. To make the figure less bulky, calibration peaks were digitally filtered out while processing.

## VI Discussion

We conclude that the results of measurements of frequency noises in microspheres confirm the theoretical model of thermorefractive fluctuations. This concerns both frequency dependence and dependence on the mode parameters $`\mathrm{}`$ and $`\mathrm{}m`$. The discrepancy is less than 30% in the frequency range of $`3\times 10^2÷10^5Hz`$ ($`\mathrm{\Omega }>1/\tau _R`$), where the model and approximations used are valid and which is free of installation noises. For frequencies above $`1kHz`$, differences are within the limits of calibration and WGM identification errors. Such a close agreement was obtained in many resonators of average size of $`120÷250\mu m`$ in diameter. In very small resonators at frequencies lower than $`1kHz`$, the noise is up to two times higher, and the dependence on frequency is more pronounced. In large resonators at lower frequencies the noise seems to be smaller. However, as it was noted, the theoretical model described below may be not valid when the frequency of noise is comparable with the inverted time of temperature relaxation of the whole microsphere via the microsphere stem ($`R_s20\mu `$ m). The influence of the stem for smaller spheres due to the ratio $`R_s/R`$ is larger. This looks compatible to the observed systematical deviation at low frequencies of experimental curves from the theoretical ones obtained for thermally non-isolated spheres. Other mechanisms of thermal noise, such as thermoelastic V.B.Braginsky et al. (1999) noise and thermal fluctuations of effective eccentricity of microspheres leading to additional frequency fluctuations proportional to $`\mathrm{}m`$ can also contribute to the effect. Additionally, there are other mechanisms of influence, which can modify spectra, such as convective thermal relaxation and influence of WGM locking scheme.

It is worth noting, that from the displayed spectra it is possible to calculate the spectra of microscopic fluctuations of temperature in the volume, occupied by the e.m. field of WGM. Using the known relation

$$\frac{\delta \omega }{\omega }=\frac{1}{n}\frac{dn}{dT}\delta T,$$
(19)
one can rescale the spectra with the factor of $`1/\left(\frac{1}{n}\frac{n}{T}\right)=10^5K`$ for fused silica. For instance, the mode of the $`138\mu m`$ resonator has an effective volume in the order of $`4,3\times 10^9`$cm<sup>3</sup>, which is 4 orders less than in former measurements of thermally induced phase noise in long fibers K.H.Wanser et al. (1993).

## Appendix A Calculation of thermorefractive noise in microspheres

To simplify the analysis we perform calculations below only for fundamental whispering–gallery mode $`TE_\mathrm{}\mathrm{}1`$, which has the smallest volume of localization. The result in principle may be extended to account for $`TE_{\mathrm{}mq}`$ modes with $`q,\mathrm{}m1`$. The field distribution of this mode may be approximated as follows:

$`\stackrel{}{𝐄}(r,\theta ,\varphi )\stackrel{}{𝐄}_\theta (r,\theta ,\varphi ){\displaystyle \frac{n\mathrm{}^{1/4}}{\sqrt{n^21}R^{3/2}\pi ^{3/4}}}`$ (22)
$``$ $`e^{\mathrm{}\mathrm{cos}^2\theta /2+i\mathrm{}\varphi }\{\begin{array}{cc}\hfill \mathrm{j}_{\mathrm{}}(knr)/\mathrm{j}_{\mathrm{}}(knR)\text{ for }rR& \\ \hfill e^{\gamma (rR)}\text{ for }r>R& \end{array},`$
$`knR`$ $``$ $`\mathrm{}+1/2+1.8558(\mathrm{}+1/2)^{1/3}{\displaystyle \frac{n}{\sqrt{n^21}}}.`$ (23)
However, even this approximation is too complex for analytical evaluation, so below we shall use the following Gaussian approximation of radial dependence:

$`\widehat{𝐞}_\theta (r,\theta ,\varphi ){\displaystyle \frac{1}{\pi \sqrt{2bdR_0}}}e^{\frac{(rR_0)^2}{2b^2}\frac{r^2\mathrm{cos}^2\theta }{2d^2}+i\mathrm{}\varphi }`$
$`knR_0\mathrm{}+1/2+0.71(\mathrm{}+1/2)^{1/3}`$
$`dR_0\mathrm{}^{1/2}`$
$`b{\displaystyle \frac{1}{kn}}\sqrt{{\displaystyle \frac{\mathrm{j}_{\mathrm{}}(knR_0)}{\mathrm{j}_{\mathrm{}}^{\prime \prime }(knR_0)}}}0.84R_0\mathrm{}^{2/3}.`$ (24)
This approximation describes rather adequately the distribution of optical energy inside the resonator and allows one to calculate Fourier integrals. Moreover, as it will be found below, because of the small depth of the field (parameter $`b`$), radial distribution practically does not influence the frequency fluctuations at frequencies of interest.

$`G(\stackrel{}{\beta })={\displaystyle \frac{1}{2\pi ^2bdR_0}}{\displaystyle \underset{0}{\overset{R}{}}}{\displaystyle \underset{0}{\overset{2\pi }{}}}{\displaystyle \underset{0}{\overset{\pi }{}}}e^{\frac{r^2\mathrm{cos}^2\theta }{d^2}\frac{(rR_0)^2}{b^2}}`$ (25)
$`\times e^{i\beta r(\mathrm{cos}\theta \mathrm{cos}\vartheta +\mathrm{sin}\theta \mathrm{sin}\vartheta \mathrm{cos}(\varphi \phi ))}\mathrm{sin}\theta d\theta d\varphi r^2dr`$
$``$ $`e^{\frac{(\beta d\mathrm{cos}\vartheta )^2}{4}}e^{\frac{(\beta b\mathrm{sin}\theta )^2}{4}}\sqrt{{\displaystyle \frac{1}{\pi \beta R_0\mathrm{sin}\vartheta }}}.`$
To obtain this result, we used several approximations while calculating the integrals. First integration was done over angle $`\theta `$, taking into account that the mode distribution is narrow and hence $`\mathrm{sin}(\theta )1`$. Second integral over $`\varphi `$ leads to Bessel function of the first order, which was approximated as spherical wave. Finally, integral over $`r`$ was found considering only fast varying functions. Other physical conditions were also used to neglect small terms: $`(RR_0)R`$, $`j(knR)j(knR_0)`$ and $`\mathrm{cos}(\beta R)1`$ — due to boundary conditions for thermal waves. In case of $`\beta 0`$, the last term in (25) should be equal to 1 and the final expression is incorrect due to used approximations for $`\mathrm{J}_0`$. However, this case is not very interesting for us. To correct the situation and obtain better approximation one may formally add 1 to the denominator.

Using (9) and keeping in mind additional factor 2 due to boundary conditions we obtain:

$`S_{\overline{u}}`$ $`=`$ $`{\displaystyle \frac{8kT^2D}{\rho C}}{\displaystyle \frac{\beta ^2|G(\stackrel{}{\beta })|^2}{a^4\beta ^4+\mathrm{\Omega }^2}\frac{d\stackrel{}{\beta }}{(2\pi )^3}}{\displaystyle \frac{2kT^2D}{\rho C\pi ^3R}}`$ (26)
$`\times `$ $`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \underset{0}{\overset{\pi }{}}}e^{\frac{(\beta d\mathrm{cos}\vartheta )^2}{2}}e^{\frac{(\beta b\mathrm{sin}\theta )^2}{2}}{\displaystyle \frac{\beta ^3d\vartheta d\beta }{D^2\beta ^4+\mathrm{\Omega }^2}},`$
$`\left[\underset{0}{\overset{\pi }{}}e^{x^2\mathrm{cos}^2(\vartheta )}𝑑\vartheta \frac{\sqrt{\pi }}{x}\right]`$,

$`S_{\delta \omega /\omega }(\mathrm{\Omega })`$ $``$ $`{\displaystyle \frac{\kappa T^2D}{\pi ^{5/2}n^2\rho CR}}{\displaystyle \frac{2}{\sqrt{d^2b^2}}}\left({\displaystyle \frac{dn}{dT}}\right)^2`$ (27)
$`\times `$ $`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\beta ^2e^{\beta ^2b^2/2}}{D^2\beta ^4+\mathrm{\Omega }^2}}{\displaystyle \frac{d\beta }{2\pi }},`$
$`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\beta ^2e^{\beta ^2b^2/2}}{D^2\beta ^4+\mathrm{\Omega }^2}}{\displaystyle \frac{d\beta }{2\pi }}{\displaystyle \frac{\sqrt{2}}{4D^{3/2}\sqrt{\mathrm{\Omega }}}}{\displaystyle \frac{1}{(1+(\mathrm{\Omega }\tau _b)^{3/4})^2}}.`$ (28)
This integral can be expressed through Lommel special functions, however, an approximation is used here, which works well for $`\mathrm{\Omega }\tau _b<1`$, where $`\tau _b=(\pi /4)^{1/3}b^2/D`$. Indeed, neglecting first item in denominator we obtain high–frequency approximation with dependence $`\mathrm{\Omega }^2`$ and neglecting exponential approximation for lower frequencies $`\mathrm{\Omega }^{1/2}`$ – dependence is obtained.

Finally:

$`S_{\delta \omega /\omega }(\mathrm{\Omega })`$ $`=`$ $`\left({\displaystyle \frac{dn}{dT}}\right)^2{\displaystyle \frac{\kappa T^2\sqrt{\mathrm{}}}{\pi ^{3/2}n^2R^2\sqrt{\lambda ^{}\rho C\mathrm{\Omega }}}}`$ (29)
$`\times `$ $`{\displaystyle \frac{1}{\sqrt{1b^2/d^2}}}{\displaystyle \frac{1}{(1+(\mathrm{\Omega }\tau _b)^{3/4})^2}}.`$
Numerical analysis which is omitted here, shows, that for the modes with $`lm`$, the power spectral density of the fluctuations is proportional to $`\sqrt{2(lm)+1}`$. This dependence is also confirmed by the fact, that azimuthal “width” of the modes, and hence their effective volume, have the same dependence.

To account for finite size of microspheres, calculations using (17) should be performed. As at room temperature the power radiated from the surface (Stephan–Boltzmann law) is much lower than the heat exchange due to thermal conductivity (Fourier law) the simplified boundary condition is used below:

$$\frac{u(r,\theta ,\varphi ,t)}{r}|_{r=R}=0,$$
(30)
$`\mathrm{\Phi }_{L,M,N}=C_{L,M,N}\mathrm{j}_L(\beta _{LN}r)\mathrm{P}_L^M(\mathrm{cos}\theta )\{{}_{\mathrm{sin}(M\varphi )}{}^{\mathrm{cos}(M\varphi )},`$ (31)
$`C_{L,M,N}^2`$ $`=`$ $`{\displaystyle \frac{2L+1}{\pi (1+\delta _{0M})}}{\displaystyle \frac{(LM)!}{(L+M)!}}`$ (32)
$`\times `$ $`{\displaystyle \frac{\xi _{L,N}^2}{R^3(\xi _{LN}^2L(L+1))\mathrm{j}_L^2(\xi _{LN})}},`$
$`G_{L,N}={\displaystyle \frac{C_{L,0,N}}{\pi bdR_0}}`$
$`\times `$ $`{\displaystyle e^{\frac{(rR_0)^2}{b^2}\frac{r^2\mathrm{cos}^2\theta }{d^2}}\mathrm{j}_L(\beta _{LN}r)\mathrm{P}_L(\mathrm{cos}\theta )r^2𝑑r\mathrm{sin}\theta d\theta }.`$
Although only functions with $`M=0`$ lead to nonzero integrals, it is appropriate to note, that functions with $`M=2\mathrm{}`$ could also be taken into account — they lead to the coupling between counter-propagating modes in the sphere and mode splitting.

To estimate this integral for small values of $`\mathrm{cos}\theta \psi =\theta \pi /2`$ near equator of microsphere, where the e.m. field is concentrated, the following approximation for the Legendre polynomial may be used for
$`L>0`$ ($`P_0(\mathrm{cos}\theta )=1`$):

$`\mathrm{P}_L(\mathrm{cos}\theta )`$ $``$ $`\sqrt{{\displaystyle \frac{2}{\pi L}}}\left(1{\displaystyle \frac{1}{4L}}\right)`$ (34)
$`\times `$ $`\mathrm{cos}\left((L+{\displaystyle \frac{1}{2}})\psi +{\displaystyle \frac{L\pi }{2}}\right),`$
and for the spherical Bessel functions:

$`j_L(z)`$ $`=`$ $`{\displaystyle \frac{1}{z}}\mathrm{sin}(zL\pi /2),`$
$`\xi _{LN}=\beta _{LN}R`$ $``$ $`{\displaystyle \frac{\pi (2N+L1)}{2}},`$
$`j_L(\xi _{LN})`$ $``$ $`{\displaystyle \frac{(1)^{N1}}{\xi _{LN}}}L=2K.`$ (35)
It is the crudeness of approximation for the roots of derivative of spherical Bessel functions that limits the applicability of the final sum we obtain below. In this way, the calculations below may be considered only as illustration.

$`|G_{L,N}|^2`$ $``$ $`{\displaystyle \frac{2(L,N)}{\pi ^2R^3}}e^{\frac{(L+1/2)^2d^2}{2R_0^2}\frac{\xi _{LN}^2b^2}{2R_0^2}};`$
$`(L,N)`$ $`=`$ $`{\displaystyle \frac{(11/4L)^2(1+1/2L)}{1L(L+1)/\xi _{LN}^2}}\text{if}L>0`$
$`(0,N)`$ $`=`$ $`\pi /4`$ (36)
In order to calculate these coefficients the approximation $`\beta _{LN}(RR_0)1`$ was used. And finally:

$`S_{\overline{u}}(\mathrm{\Omega })`$ $``$ $`{\displaystyle \frac{8kT^2}{\pi ^2\rho DCR}}`$
$`\times `$ $`{\displaystyle \underset{K=0}{\overset{\mathrm{}}{}}}{\displaystyle \underset{N=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\xi _{LN}^2e^{\frac{(L+1/2)^2d^2}{2R_0^2}\frac{\xi _{LN}^2b^2}{2R_0^2}}}{\xi _{LN}^4+\frac{\mathrm{\Omega }^2R^4}{D^2}}}(L,N).`$
To calculate the sums in high–frequency approximation, we can consider $`(L,N)=1`$. By making the following substitutions: $`x=(K+N)/\sqrt{2}`$, $`y=(KN)/\sqrt{2}`$ and $`\tau _R=R^2/D`$, we use integrals instead of sums:

$`S_u`$ $``$ $`{\displaystyle \frac{8\kappa T^2}{\pi ^2\rho CDR}}`$ (38)
$`\times `$ $`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \underset{x}{\overset{x}{}}}{\displaystyle \frac{2\pi ^2x^2}{4\pi ^4x^4+\mathrm{\Omega }^2\tau _R^2}}e^{\pi ^2x^2b^2/R^2}e^{(x+y)^2/\mathrm{}}𝑑y𝑑x`$
$``$ $`{\displaystyle \frac{8\kappa T^2}{\pi ^2\rho CDR}}\sqrt{{\displaystyle \frac{\mathrm{}}{8\pi }}}`$
$`\times `$ $`{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{t^2}{t^4+\mathrm{\Omega }^2\tau _R^2}}e^{t^2b^2/(2R^2)}\mathrm{erf}\left({\displaystyle \frac{\sqrt{2}t}{\pi \sqrt{\mathrm{}}}}\right)𝑑t.`$
This integral coincides with (27) obtained for infinite media, if erf–function, which is practically equal to unity for $`t>1/\sqrt{\mathrm{}}`$, is ignored. However, spectral densities described by integral (27) and sum (A) are different at frequencies for $`\mathrm{\Omega }<1/\tau _R`$ where (A) is finite for $`\mathrm{\Omega }0`$ (the term with two sums is close to unity).

## Appendix B Relation between regular calibration signal and noise spectrum

To find the link between the amplitude of harmonic frequency change, corresponding peak in calculated spectral density and noise spectral density, we consider the Fourier transform for regular and chaotic processes. Let $`x(t)`$ denotes chaotic process with correlation function of

$$B(\tau )=<x(t)x(t+\tau )>=\underset{\mathrm{}}{\overset{\mathrm{}}{}}S^\pm (\mathrm{\Omega })e^{i\mathrm{\Omega }\tau }\frac{d\mathrm{\Omega }}{2\pi },$$
(39)
where $`S^\pm (\mathrm{\Omega })`$ stands for “double–sided” spectral density which is symmetrical and two times less for positive frequencies than “one–sided” $`S(\mathrm{\Omega })`$. Let $`T`$ denotes the length of the “window” $`W(t)`$ in time domain, that is the duration of a set $`\{x(t)\}`$ (or $`\{x_i\}`$ in discrete case) for which the Fourier transform is carried out. Now we shall consider the transform with the square window $`w(t)=1`$ for $`t[T/2,T/2]`$. Fourier transform for $`x(t)`$ is:

$$X_j=\underset{T/2}{\overset{T/2}{}}W(t)x(t)e^{i\mathrm{\Omega }_jt}𝑑t,\text{where}\mathrm{\Omega }_j=\frac{2\pi j}{T}.$$
(40)
In the procedure of spectrum estimation the spectra $`X_j(\mathrm{\Omega })`$ are calculated for small intervals of data and then averaged, which in our notations corresponds to the calculation of the quantity $`|X_j|^2`$.

$`|X_j|^2`$ $`=`$ $`{\displaystyle \underset{T/2}{\overset{T/2}{}}}x(t)x(t^{})e^{i\omega _j(tt^{})}𝑑t𝑑t^{}`$ (41)
$`=`$ $`{\displaystyle \underset{\mathrm{}}{\overset{\mathrm{}}{}}}S(\omega ){\displaystyle \frac{4\mathrm{sin}^2[(\omega _j\omega )\frac{T}{2}]}{(\omega _j\omega )^2}}{\displaystyle \frac{d\omega }{2\pi }}.`$
Calibrating peak in our case corresponds to frequency $`\omega _j>>T^1`$. Then, the core in the latter integral has sharp maximum at $`\omega =\omega _j`$ and $`|X_j|^2=TS(\omega _j)`$.

Let now $`y(t)=Y_0\mathrm{cos}(\mathrm{\Omega }t)`$ represents a regular harmonic signal. The Fourier transform for such signal for square window will look as follows

$`Y_s`$ $`=`$ $`{\displaystyle \underset{T/2}{\overset{T/2}{}}}Y_0\mathrm{cos}(\mathrm{\Omega }t)e^{i\mathrm{\Omega }_jt}𝑑t,`$ (42)
leading to $`Y_{sj}=Y_0T/2`$. Now we compare the results of spectrum estimation procedure (41) and the results of Fourier transform for harmonic signal:

$$Y_{sj}^2=\frac{Y_0^2T^2}{4}TS^\pm (\omega _j)=|X_j|_s^2.$$

Here one can easily see the link between the amplitude of harmonic signal and corresponding peak in spectral density. If the harmonic signal has the dimension of frequency and the desired spectral density represents relative frequency changes, then the final relation, corresponding to one–sided spectral density of calibration spike with square window in Fourier transform, will be

$$\sqrt{S_{\delta f/f}}\frac{Y_0}{f_0}\sqrt{\frac{T}{2}}\left[\frac{\text{1}}{\sqrt{Hz}}\right],$$
(43)
where $`f_0`$ — is optical frequency. The same calculations for the Hann window give:

$$\sqrt{S_{\delta f/f}}\frac{Y_0}{f_0}\sqrt{\frac{T}{3}}\left[\frac{\text{1}}{\sqrt{Hz}}\right].$$
(44)
## Appendix C Formulas for WGMs identification

To identify the parameters of whispering–gallery modes from photographs the following approximation for half–width of WGM belt was used:

$$\mathrm{\Delta }\theta =\sqrt{2(\mathrm{}+1/2m)/\mathrm{}}.$$
(45)
To calculate $`\mathrm{}`$–index from the radius of microsphere we found by numerical approximation the following approximated formulas valid for $`R\lambda `$:

$`q=1:`$ $`l=0.5+t_{lq}2.287t_{lq}^{1/3}+0.1718t_{lq}^{1/3},`$ (46)
$`q=2:`$ $`l=0.5+t_{lq}4.617t_{lq}^{1/3}+0.6944t_{lq}^{1/3},`$
$`q=3:`$ $`l=0.5+t_{lq}6.895t_{lq}^{1/3}+1.518t_{lq}^{1/3},`$
$`q=4:`$ $`l=0.5+t_{lq}9.190t_{lq}^{1/3}+2.632t_{lq}^{1/3},`$
where

$`t_{lq}={\displaystyle \frac{2\pi n_iR}{\lambda _{lmq}}}+{\displaystyle \frac{pn_i}{\sqrt{n_i^21}}}.`$ (47)
These formulas allow to get an estimate of $`\mathrm{}`$ and $`\mathrm{}m`$ values from the radius of the microsphere and width of the belt of the mode visualized due to residual surface scattering.

###### Acknowledgements.

We are grateful to professor V.B.Braginsky, whose initial idea to measure thermorefractive noise in microspheres has embodied in this paper. Special thanks to professors S.P.Vyatchanin and F.Ya. Khalili for stimulating and fruitful discussions. This work was supported partially by NSF grant #PHY-0098715 and Russian Foundation for Basic Research grant #00-15-96706.


format: text

text:
# Narrow Resonances in Effective Field Theory

## Acknowledgments

UvK is grateful to the Nuclear Theory Group at the University of Washington for its hospitality, and to RIKEN, Brookhaven National Laboratory and the U.S. Department of Energy \[DE-AC02-98CH10886\] for providing the facilities essential for the completion of this work. This research was supported in part by the Director, Office of Energy Research, Office of High-Energy and Nuclear Physics, by the Office of Basic Energy Sciences, Division of Nuclear Sciences, of the U.S. Department of Energy under contract DE-AC03-76SF00098 (PFB), by a DOE Outstanding Junior Investigator Award (UvK), and by an Alfred P. Sloan Fellowship (UvK).


text:
# Breakdown of the Onsager Principle as a Sign of Aging

## I introduction

The Onsager principle onsager  is one of the basic tenets of statistical mechanics insofar as it establishes a connection between a property of equilibrium, the correlation function of a given variable $`A`$, and the regression to equilibrium of a macroscopic signal. For this reason, we judge the Onsager principle to be a fundamental step for the connection between dynamics and thermodynamics. It is important to stress, as clearly stated by Onsager himself onsager , that this principle holds true for aged systems, namely systems in contact with heat reservoirs that are supposed to be in a condition of thermal equilibrium. If the regression to equilibrium is very fast, it is not so important to ensure the equilibrium condition of the bath, at the moment when we begin measuring the regression to equilibrium of the system of interest. In fact, in the specific case where the bath is responsible for fluctuations that can be assumed to be white, the regression to equilibrium of the bath is instantaneous. The Onsager principle refers to a variable of interest whose dynamics are made stochastic by the interaction with a bath. Thus, when we discuss the process of regression to equilibrium, we have to specify if we are referring to the system of interest or to its bath. As earlier said, if we adopt the white noise approximation to describe the fluctuations that are responsible for the erratic motion of the variable of interest, then the regression to equilibrium of the reservoir is virtually instantaneous and we can easily fulfill the condition for the validity of the Onsager principle. The variable of interest is characterized by a stationary correlation function $`\mathrm{\Phi }_A(t_1,t_2)=A(t_1)A(t_2)`$, which only depends on $`|t_1t_2|`$. The stochastic behavior of the variable of interest is caused by the interaction between system of interest and bath, and this kind of process is often studied by means of the master equation method. A popular method to derive the master equation is the projection method by Zwanzig zwanzig . However, this method is easy and convenient to use, when special initial conditions can be adopted, with the total distribution expressed as the product of the relevant, or reduced distribution, and of the bath distribution grigo . These special initial conditions have the beneficial effect of annihilating the inhomogeneous term that makes the reduced equation of motion explicitly dependent on the initial condition grigo . If no initial condition of this kind is adopted, we have to consider also the inhomogeneous term and we have to wait a time comparable to the bath relaxation time for this term to disappear. This suggests that a non-Markov master equation with the same structure as the Zwanzig master equation, with no inhomogeneous term, does not necessarily refer to a stationary condition. It might refer to an initial condition with system of interest and bath uncorrelated one from the other, a non-stationary state. As we shall see, to create a master equation compatible with the Onsager principle, in the case when the bath is not infinitely fast, we must establish a delicate condition of entanglement between system and bath. This is not a trivial problem, since the departure from the Poisson condition generates an infinitely extended memory, and the system-bathentanglement is the result of a rearrangement process with an infinite time scale.

In the last few years there has been a growing interest for the formalism of Continuous Time Random Walk (CTRW) montrollweiss , because of its close connection with the adoption of fractional operators that, in turn, are revealing a powerful tool to describe cooperative processes in condensed matter bruce . On the other hand, the connection between the master equation formalism and the CTRW walk has been discussed over the years, beginning with the pioneering work of Ref. katja . According to the CTRW the walker travelling through an one dimensional path makes jumps from one site to the others after spending a given amount of time in the departure site. The distribution of sojourn times is assumed to be the same for all the sites, and is denoted by the function $`\psi (t)`$. The pioneering work of Ref. katja  proves that the Markovian master equation is compatible with the CTRW, if $`\psi (t)`$ is a Poissonian distribution.

The results of the pioneering work of Ref. katja  raise the important issue of the connection between CTRW, with non-Poissonian $`\psi (t)`$, and a non-Markovian master equation. Apparently, this important problem is solved by adopting the non-Markovian master equation proposed by the authors of Ref. kenkre . We consider this theoretical tool to be very important, and we denote with the term Generalized Master Equation (GME) all transport equations with the same structure. Recently the GME has been discussed by Metzler metzler  who argued that this equation unifies fractional calculus and CTRW. We want to remark, however, that here we focus our attention on how to make the GME stationary, and thus compatible with the Onsager principle. In fact, as earlier remarked, the absence of an inhomogeneous term does not guarantee the stationary nature of the transport equation.

We want to stress that another reason of interest of the GME of Kenkre, Montroll and Shlesinger is that it is formally identical to the master equation that the authors of Ref. mazza  derived from a quantum mechanical tight-binding model, with an erratic distribution of energy sites. The calculation was done along the lines established by Zwanzig zwanzig  using a projection operator similar to that used by Kenkre in an earlier work kenkre2 . This means that in principle understanding the connection between a master equation and CTRW might contribute to shedding light on the intriguing issue of unravelling quantum mechanical master equations. This problem becomes challenging in the case of non-Markovian processes. For some very recent references on this issue, we refer the reader to Refs. gambetta ; schreiber . The purpose of this paper is not so ambitious as to afford direct contributions to the settlement of this problem. However, we focus our attention on the case of a two-site system that seems to be closer to the problem of decoherence of a q-bit than to the random walk over an infinite path for which the CTRW theory was originally designed. We address a problem that is closely connected to that recently discussed by Sokolov and Metzler sokolovandmetzler , the discussion of the connection between the CTRW and a two-state non-markovian master equation. We do not discuss, however, the connection with the fractional calculus, an interesting issue addressed by the authors of Ref. sokolovandmetzler . Rather we focus on the problem of how to make both the CTRW and the equivalent GME compatible with the stationary condition and the Onsager principle. For instance, we define as GME the transport equation found by the authors of Ref. mazza , because it has the same structure as that proposed by Kenkre, Montroll and Shlesinger kenkre . However, the quantum transport equation of the type of Ref. mazza  might be non stationary, due to the fact it refers to an initial condition with the system totally disentangled from the bath. In this paper we shall show how to build up a genuinely stationary GME.
The Onsager principle is a property that makes it possible for us to derive the correlation function from the GME. If we set a condition of total equivalence between the GME and the CTRW, we find the apparently disconcerting result that only the Poisson statistics is compatible with the equivalence condition. The reason for this restriction is due to the fact that the departure from the Poisson statistics generates memory properties that make the GME incompatible with the Markov approximation. This means that the structure of the GME is dictated by the initial condition. If this is not stationary, the resulting GME is not a *bona fide* transport equation fox . On the other hand, if the waiting time distribution $`\psi (t)`$ is not exponential, there are aging effects. This means that we have to leave the system age till it reaches the condition where the Onsager principle holds true. In this condition it is possible to establish a *bona fide* GME. We do it, and in so doing we establish a complete equivalence between CTRW and GME. At same time, we prove what it is probably the most important result of this paper. This is that the transition from the Poisson to the non-Poisson statistics has the effect of creating a memory condition incompatible with the Markov approximation, at least in principle. In fact, the Markov approximation, as we shall see, means that the bath is infinitely fast, and this condition, in turn, is proven to force the waiting time distribution $`\psi (t)`$ to be exponential. This means that the Markov condition is incompatible with the non-Poissonian nature of the waiting time distribution.

## II The Poisson case and the Onsager Principle

In this section we prove that apparently the equivalence between the CTRW of Montroll and Weiss and the GME is complete only in the Poisson case. Let us express the GME under the following concise form:

$$\frac{d}{dt}𝐩(t)=_0^t\mathrm{\Phi }(tt^{})K𝐩(t^{})𝑑t^{},$$
(1)
where $`𝐩`$ is the $`m`$-dimensional population vector of $`m`$ sites, $`K`$ a transition matrix between the sites and $`\mathrm{\Phi }(t)`$ the memory kernel. The prescription of the CTRW montrollweiss  yields:

$$𝐩(t)=\underset{n=0}{\overset{\mathrm{}}{}}_0^t𝑑t^{}\psi _n(t^{})\mathrm{\Psi }(tt^{})M^n𝐩(0).$$
(2)
Note that $`\psi _n(t)`$ is the probability that $`n`$ jumps occurred and that the last took place at time $`t=t^{}`$. This means that

$$\psi _n(t)=_0^t\psi _{n1}(tt^{})\psi _1(t^{})𝑑t^{},$$
(3)
with

$$\psi _0(t)=\delta (t).$$
(4)
While $`M`$ is the transition matrix connecting the sites after one jump has occurred. The function $`\mathrm{\Psi }(t)`$ is the probability that no jump occurs in the time interval $`t`$, namely,

$$\mathrm{\Psi }(t)=_t^{\mathrm{}}\psi _1(t^{})𝑑t^{}.$$
(5)
It is evident that the time convolution structure of both Eq. (1) and Eq. (2) makes it straigthforward to derive in both cases the Laplace transform of $`𝐩(t)`$. Furthermore, Eq. (3) and Eq. (5) makes it possible to express the Laplace transform of $`𝐩(t)`$, in the second case, in terms of the function $`\psi _1(t)`$ only. This function is the earlier mentioned waiting time distribution, and from now on it will be denoted with the symbol $`\psi (t)`$. By comparing the Laplace transform of the GME to the Laplace transform of the CTRW we get

$$\widehat{\mathrm{\Phi }}(u)=\frac{u\widehat{\psi }(u)}{1\widehat{\psi }(u)}\frac{MI}{K},$$
(6)
where $`I`$ is the unity operator. We limit our discussion in this paper to the two-state case, where

$$M=\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right)$$
(7)
and

$$K=\left(\begin{array}{cc}1& 1\\ 1& 1\end{array}\right),$$
(8)
thereby turning Eq. (6) into

$$\widehat{\mathrm{\Phi }}(u)=\frac{u\widehat{\psi }(u)}{1\widehat{\psi }(u)}.$$
(9)
We note that the two-state master equation studied in this paper has to be related to the two-site CTRW. If we assign the value $`W`$ to the right site, and the value $`W`$ to the left site, and we adopt a discrete time representation, the motion of the random walker corresponds to a symbolic sequence $`\{\xi \}`$, with the form $`\{WWWWWWWWWWWWWWWWWWWWWWWWW\mathrm{}\}`$, which shows a significant time persistence of both states. The waiting time distribution $`\psi (t)`$ is the distribution of the patches filled with either $`W`$’s or $`W`$’s. We assume a symmetric condition. In the case when the first moment of $`\psi (t)`$ is finite, the process is stationary renewal ; renewaltheory  thereby making it possible to define the correlation function of the fluctuation $`\xi (t)`$, denoted with the symbol $`\mathrm{\Phi }_\xi (t)`$.

It is now the right time to use the Onsager principle. We set the Onsager principle in the form

$$\mathrm{\Phi }_\xi (t)=\frac{p_1(t)p_2(t)}{p_1(0)p_2(0)},$$
(10)
where $`p_1(t)`$ and $`p_2(t)`$ are the probabilities for the random walker to be, at time $`t`$, in the first and second state, respectively. In his original work Onsager referred himself to the case of a macroscopic fluctuation that is supposed to regress to the vanishing equilibrium value through a phenomemological equation of motion. We realize this macroscopic fluctuation selecting a large number $`N`$ of walkers, divided into two groups, with $`N_1`$ and $`N_2`$ walkers belonging to the first and second state, respectively. Then we relate this choice to the probabilities $`p_1(t)`$ and $`p_2(t)`$, by setting $`p_1=N_1/N`$ and $`p_2=N_2/N`$. The regression to equilibrium of this macroscopic fluctuation does not fit any phenomenological law. For this reason, we can refer ourselves to Eq. (10) as a generalized version of the Onsager principle. We plan to derive the mathematical expression of the law of regression to equilibrium from the use of this generalized version of the Onsager’s regression principle, which, as shown in Section 3, is exact in the dichotomous case. Due to its nature, the GME can afford information only on $`p_1(t)`$ and $`p_2(t)`$. There is no direct information on the correlation function $`\mathrm{\Phi }_\xi (t)`$, and the only possible way is through the Onsager’s principle, which, as shown by Eq. (10), implies that an initial off-equilibrium condition, $`p_1(0)p_2(0)0`$, is established. Then, the regression to equilibrium is expected to take place in such a way as to establish a connection with the correlation function in the form of Eq. (10). This assumption, applied to Eq. (1), yields:

$$\frac{d}{dt}\mathrm{\Phi }_\xi (t)=2_0^t𝑑t^{}\mathrm{\Phi }_\xi (tt^{})\mathrm{\Phi }(t^{}).$$
(11)
This equation is the non-Markovian counterpart of the phenomenological regression to equilibrium of the original work of Onsager onsager . Using the time convoluted expression of Eq. (11) and the central result of Eq. (9), we establish a connection between the correlation function $`\mathrm{\Phi }_\xi (t)`$ and the waiting time distribution $`\psi (t)`$, through their Laplace transforms, as follows:

$$\widehat{\mathrm{\Phi }}_\xi (u)=\frac{1}{u+\frac{2u\widehat{\psi }(u)}{(1\widehat{\psi }(u))}}.$$
(12)
On the other hand, the correlation function is connected to $`\psi (t)`$ through properties established by the renewal theory renewal . To properly establish this connection, we have to notice, with Zumofen and Klafter klafter , that it is convenient to introduce another type of waiting time distribution that we call $`\psi ^{}(t)`$. What’s the connection between $`\psi (t)`$ and $`\psi ^{}(t)`$? The waiting time distribution $`\psi (t)`$ is the experimental waiting time distributon. It could be evaluated experimentally observing the sequence $`\{\xi \}`$, and recording the time length of the laminar regions occupied only by $`W`$ or by $`W`$. However, we can imagine that a theoretical waiting time distribution exists, denoted by $`\psi ^{}(t)`$, and that the sequence $`\xi (t)`$ is obtained as follows. We select randomly a number $`t_1`$ from the distribution $`\psi ^{}(t)`$. We toss a coin, and assign to the first laminar region, of length $`t_1`$, the symbol $`W`$ or $`W`$, according to the coin tossing prescription. At the end of this laminar region, first we select, again from the distribution $`\psi ^{}(t)`$, a number $`t_2`$. This is the length of the second laminar region. We toss the coin again to decide the sign of it. It is evident that there is $`50\%`$ probability of getting the same sign as the earlier laminar region. We proceed in the same way with the length and the sign of third laminar region, and so on. We adopt this rule to create the sequence $`\{\xi (t)\}`$. Thus, from the renewal theory renewal  we obtain the following important result:

$$\mathrm{\Phi }_\xi (t)=\frac{1}{t}_t^{\mathrm{}}(t^{}t)\psi ^{}(t^{})𝑑t^{},$$
(13)
where $`t`$ is the mean waiting time of the $`\psi ^{}(t)`$-distribution. It is interesting to notice that this equation means that the second derivative of the correlation function is proportional to $`\psi ^{}(t)`$,

$$\frac{d^2}{dt^2}\mathrm{\Phi }_\xi (t)=\frac{\psi ^{}(t)}{t}.$$
(14)
Zumofen and Klafter klafter  explained with clear physical arguments the connection between $`\psi (t)`$ and $`\psi ^{}(t)`$. They prove that their Laplace transforms are related one to the other by

$$\widehat{\psi }^{}(u)=\frac{2\widehat{\psi }(u)}{1+\widehat{\psi }(u)}.$$
(15)
Using Eq. (14) and Eq. (15) we derive a further connection between the Laplace transform of $`\mathrm{\Phi }_\xi (t)`$ and the Laplace transform of $`\psi (t)`$. In fact, by Laplace transforming Eq. (14) and using Eq. (15), we derive a new expression for the Laplace transform of $`\mathrm{\Phi }_\xi (t)`$. By equating this new expression to that of Eq. (12), we find, after some algebra, the following form for the Laplace transform of $`\psi (t)`$,

$$\widehat{\psi }(u)=\frac{\frac{1}{2}\mathrm{\Phi }_\xi ^{^{}}(0)}{u+\frac{1}{t}+\frac{1}{2}\mathrm{\Phi }_\xi ^{^{}}(0)}.$$
(16)
Using Eq. (13) we prove that

$$\mathrm{\Phi }_\xi ^{^{}}(0)=\frac{1}{t}.$$
(17)
This equation allows us to write Eq. (16) under the form

$$\widehat{\psi }(u)=\frac{1}{2t}\frac{1}{(u+\frac{1}{2t})}.$$
(18)
This means that $`\psi (t)`$ is exponential and the explicit forms of $`\psi (t)`$ and $`\mathrm{\Phi }_\xi (t)`$ are given by

$$\psi (t)=\frac{\gamma }{2}\mathrm{exp}(\frac{\gamma }{2}t)$$
(19)
and

$$\mathrm{\Phi }_\xi (t)=\mathrm{exp}(\gamma t),$$
(20)
respectively.

This is the first important result of this paper. It is to some extent disconcerting since it seems to restrict the complete equivalence between GME and CTRW to the exponential case, with the effect of making the GME useless, since in the exponential case, it turns out to be an ordinary, memory-less, master equation. We want to notice that Sokolov and Metzler sokolovandmetzler  have recently discussed the derivation of a fractional transport equation totally equivalent to the CTRW. However, they did not discuss the intimately related problem of the formal equivalence between CTRW and GME, which is the main goal of the present paper.

The second important result of this paper is the discovery of a GME that is equivalent to the CTRW, with no restriction to the Poisson case. This second important result will be illustrasted in Section 5. To derive this important result we have to mention a fact recently observed by Barkai barkai : the Montroll-Weiss CTRW refers to a non-stationary condition. This means that the Onsager’s regression principle is invalidated. Therefore, we shall conclude that the disconcerting result of this section is due to the fact that the case of Poisson statistics is the only one where aging does not exist. At this stage, this is only a conjecture. We can notice, however, that this is a plausible conjecture. In fact, the exponential condition of Eq. (19) when applied to Eq. (9) makes $`\widehat{\mathrm{\Phi }}(u)`$ independent of $`u`$, thereby implying that instantaneous regression to the vanishing value of the memory kernel $`\mathrm{\Phi }(t)`$. This means that the fluctuations of the variable $`\xi `$ take place with its bath at equilibrium.

We give a further support to the conjecture that aging is the reason of the apparently disconcerting result of this section by referring ourselves to the work of Ref. mazza . The authors of this paper found a transport equation with the same form as that of Eq. (1) , with no inhomogeneous term. The adoption of the Zwanzig projection method makes this interesting result possible, provided that special initial conditions are chosen. These conditions might depart from the entanglement between the system of interest and its bath, necessary to ensure the stationary condition. On the other hand, the Hamiltonian treatment implies that the system exerts an action on its bath. Thus, the choice of an initial condition annihilating the inhomogeneous term might have the effect of perturbing the bath, with an ensuing long lasting re-equilibration process, and so a persistent non-stationary condition.

## III Onsager’s regression for dichotomous signals

We note that the constraint posed in the earlier section on the waiting time distribution $`\psi (t)`$, to ensure the equivalence between the Master equation and the CTRW, has been based on the assumption that the Onsager principle is valid. The fact that the equivalence between CTRW and GME is restricted only to the Poissonian case, and so to the case of ordinary master equation, might generate the false impression that the Onsager principle is valid only in the Poissonian case. It is not so. It is already known zambon  that in the case of Gaussian statistics the Onsager’s regression hypothesis is exact, and it holds true for initial excitation of any intensity. In this section, having in mind the two-state master equation, we focus our attention on the case of dichotomous statistics. In the dichotomous case, we show that the Onsager’s regression hypothesis turns out to be exact again, for initial departures from equilibrium of whatsoever intensity, and with any type of relaxation process, with an exponential form, and with an inverse power law form as well, provided that the system is ergodic. In short, in this section we prove that in the dichotomous case the Onsager regression principle holds true provided that the system is ergodic.

The ergodic property implies that we can create a Gibbs ensemble in two different, but equivalent, ways. The first refers to the original idea of Gibbs, that the same system is repeated infinitely many times. This means that we have to generate infinitely many sequences using the same physical prescription. The second way is based on the adoption of only one sequence, denoted by $`\xi (t)`$, which is supposed to be infinitely long. Then we define a generic trajectory $`\xi ^{(s)}(t)`$ with the following prescription

$$\xi ^{(s)}(t)=\xi (t+s).$$
(21)
In the practical case of a numerical treatment, time is discrete, time unity being, for instance $`t=1`$. In this case the superscript $`s`$ is an integer number.

We denote by $`P_+`$ and $`P_{}`$ the probability that the variable $`\xi `$ gets the value $`W`$ and $`W`$, respectively. These probabilities are frequencies that can be evaluated using either the first or the second ensemble of sequences. We divide our ensemble of sequences into two groups, characterized by the the initial condition $`\xi (0)=1`$, and $`\xi (0)=1`$, respectively. We define two distinct averages for the two groups, denoted by $`_+`$ and $`_{}`$, respectively. Of course, we have that:

$$=P_+_++P_{}_{},$$
(22)
with

$$P_++P_{}=1.$$
(23)
Although it is not essential for the main goal of this paper, we make our discussion as general as possible. Thus, we do not assume the two probabilities to be equal. We do not rule out that a bias might exists, given by

$$P_+P_{}=c.$$
(24)
In the ergodic, and stationary, condition that we are assuming, the stationary autocorrelation function exists, and it is defined as

$$\mathrm{\Phi }_\xi (t)\frac{\xi (0)\xi \xi (t)\xi }{\xi ^2\xi ^2}=\frac{\xi (0)\xi (t)c^2}{1c^2}.$$
(25)
Expressing the total average in terms of the averages over the two groups, according to the prescription of Eq. (22), we make Eq. (25) become

$`\mathrm{\Phi }_\xi (t)`$ $`=`$ $`{\displaystyle \frac{(+1)\xi (t)_++(1)\xi (t)_{}c^2}{1c^2}}`$ (26)
$`=`$ $`{\displaystyle \frac{\xi (t)_+\xi (t)_{}c^2}{1c^2}}.`$
Note that on the same token the property $`\xi =c`$ can be written as follows

$$P_+\xi _++P_{}\xi _{}=c.$$
(27)
We also note that

$$P_+=\frac{c+1}{2}.$$
(28)
Using Eqs.(27) and (28) we express the correlation function of Eq. (26) as follows:

$$\mathrm{\Phi }_\xi (t)=\frac{2P_+\xi _+cc^2}{1c^2}=\frac{\xi _+(c+1)c(c+1)}{1c^2}.$$
(29)
Finally, by dividing both numerator and denominator of the ratio corresponding to the last term of Eq. (29) by $`1+c`$, we obtain

$$\mathrm{\Phi }_\xi (t)=\frac{\xi _+c}{1c}.$$
(30)
Eq. (30) can be interpreted as follows. The group of trajectories corresponding to the condition $`\xi (0)=1`$ can be thought of as a way of creating an out of equilibrium condition, $`\xi _+`$. This out of equilibrium condition undergoes a process of regression to equilibrium that is proportional to the equilibrium correlation function. Using Eq. (27), we can express the correlation function in terms of $`\xi _{}`$, and we reach the same conclusion. It is evident that the same conclusion would be reached using an arbitrary mixture of $`\xi _+`$ and $`\xi _{}`$, departing from the vanishing equilibrium value. In conclusion, the Onsager’s regression principle is fulfilled.

We notice that the waiting time distributions that we are considering in this paper have an inverse power law nature with the power index $`\mu `$. If this power index fits the condition $`\mu >2`$ the ergodic condition is ensured. We shall focus our attention on this condition, and we shall prove that we can make the CTRW compatible with the GME if we set the constraint that the system is aged enough as to make the Onsager principle valid. As we shall see, this means the adoption of a form of CTRW different from that of Montroll and Weiss montrollweiss , which corresponds to a condition very far from the stationary state. It is evident that the ergodic condition might become hard to fulfill in practice, with relaxation processes described by inverse power laws. This is in fact the case when aging becomes important.

## IV Aging in renewal processes

In Section 2 we found that the GME and the CTRW are equivalent only in the Poisson case. Since the use of the Onsager’s regression assumption is the key ingredient used to establish this equivalence, we might be tempted to conclude that the Onsager principle does not hold true in the non-Possonian case. In Section 3, we found that it is not so, and that the Onsager principle is valid, provided that $`\mu >2`$. Here we shall prove that the non-Poisson case results in aging effects, which has significant consequences also in the case $`\mu >2`$. If the CTRW used does not refer to the stationary (aged) condition, it cannot be compatible with the Onsager regression principle. This is the second result of this paper. The problem of aging within the context of intermittent processes has been discussed recently in a very attractive paper by Barkai barkai . The work of Barkai focuses on the non-stationary condition, $`\mu <2`$, where no invariant measure exists, and aging lasts forever. In this case the Onsager regression principle is violated for reasons explicitly pointed out by Onsager himself. In fact, Onsager repeatedly mentioned onsager  that his principle has to be applied to aged systems. The condition of aged systems implies a state of thermodynamic equilibrium. In the case $`\mu <2`$ equilibrium is never reached, not even asymptotically. Thus, in this case the Onsager’s principle is always violated. In Section 4A, using the renewal theory renewaltheory , we shall study the case $`\mu >2`$, and we shall discuss two limiting conditions, corresponding to the birth and the death of the system, respectively. It is interesting to remark that death, meaning an ordinary thermodynamical condition, takes place eventually at the end of an infinitely long aging process.

### IV.1  Young and aged systems

First of all, we should define, as Barkai does barkai , the waiting time distribution $`\psi _{t_a}(t)`$. This means that we establish non-stationary conditions at time $`t_a<0`$, and we begin our observation at time $`t=0`$. The waiting time distribution depends on $`t_a`$. The naive conviction that the waiting time distribution is given by $`\psi (t)`$, actually is based on the assumption that $`t_a=0`$. This means, in fact, that we are considering a set of random walkers and that at time $`t=0`$ all of them begin their sojourn in the laminar region. This leads immediately to $`\psi _{t_a=0}(t)=\psi (t)`$. It is straighforward to evaluate also the distribution of sojourn times corresponding to $`t_a=\mathrm{}`$. We denote this distribution with the symbol $`\psi _{\mathrm{}}(t)`$.

This is the stationary case, corresponding to the following procedure. In the stationary condition, the probability of selecting a laminar zone of length $`T`$, by a random choice, is $`\frac{T}{\tau }\psi (T)dT`$, where $`\tau `$ is the mean length of a laminar zone. The probability density of observing the first change of laminar phase after a time $`t`$, being in a laminar zone of length $`T`$, is $`\mathrm{\Theta }(Tt)\frac{1}{T}`$. Consequently for the probability density of having the first change of laminar phase at time $`t`$ (i.e. $`\psi _{\mathrm{}}(t)`$ ), we have, integrating over all possible $`T`$:

$`\psi _{\mathrm{}}(t)`$ $`=`$ $`{\displaystyle \frac{1}{\tau }}{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}dTT\psi (T){\displaystyle \frac{1}{T}}\mathrm{\Theta }(Tt)`$ (31)
$`=`$ $`{\displaystyle \frac{1}{\tau }}{\displaystyle \underset{t}{\overset{\mathrm{}}{}}}dT\psi (T).`$
Note that due to the renewal theory, the distribution $`\psi _{\mathrm{}}(t)`$ concerns only the time that we have to wait to detect the first event. After the first event we have a total rejuvination. In fact, measuring the time at the moment of the first jump is equivalent to beginning the measurement process at the precise moment when the walker enters the laminar region. Thus the Laplace transform of $`𝐩(t)`$ is given by

$$\widehat{𝐩}(u)=(\widehat{\mathrm{\Psi }}_{\mathrm{}}(u)I+\frac{\widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)M}{1\widehat{\psi }(u)M})𝐩(0).$$
(32)
This result is derived from Eq. (2) as follows. First of all, we evaluate the Laplace transform of Eq. (2). We obtain

$`\widehat{𝐩}(u)`$ $`=`$ $`(\widehat{\mathrm{\Psi }}(u)+\widehat{\mathrm{\Psi }}(u)\widehat{\psi }(u)M+\widehat{\mathrm{\Psi }}(u)\widehat{\psi }^2M^2+\mathrm{})𝐩(0)`$
$`=`$ $`\widehat{\mathrm{\Psi }}(u){\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}(\widehat{\psi }M)^n𝐩(0)={\displaystyle \frac{\widehat{\mathrm{\Psi }}(u)}{1\widehat{\psi }(u)M}}𝐩(0).`$
Then, we replace the probability of occurrence, or of non occurrence, of the first event, called $`\psi (t)`$ and $`\mathrm{\Psi }(t)`$, respectively, with the corresponding aged quantities. These are called $`\psi _{\mathrm{}}(t)`$ and $`\mathrm{\Psi }_{\mathrm{}}(t)`$, respectively. The function $`\psi _{\mathrm{}}(t)`$ is given by Eq. (31) and the function $`\mathrm{\Psi }_{\mathrm{}}(t)`$ is given by

$$\mathrm{\Psi }_{\mathrm{}}(t)=_t^{\mathrm{}}\psi _{\mathrm{}}(t^{})𝑑t^{}.$$
(33)
All this yields Eq. (32).

We expect that the inverse Laplace transform of Eq. (IV.1), is a function of $`t`$ that asymptotically will become equivalent to the asymptotic value of the inverse Laplace transform of Eq.(32). What about the Onsager prescription of Eq. (10)? If we adopt the GME corresponding to the young condition, the adoption of Eq. (10) would be equivalent, as far as the left-hand-side of this equation is concerned, to adopting the non-stationary condition. We have to remark, in fact, that the GME corresponding to the Montroll-Weiss CTRW is not a *bona fide* master equation. The importance of a *bona fide* transport equation has been pointed out in a remarkable paper by Fox fox . A master equation is a *bona fide* master equation when it can be used with any initial condition (like the *bona fide* Fokker-Planck equation of Ref. fox ) . The master equation of Section 3, on the contrary, implies the choice of only one initial condition. The next subsection is devoted to sheding more light into this important aspect.

### IV.2 A theory for systems of any age

The purpose of this subsection is to shed further light into the aging problem, and into the ensuing conflict with the Onsager postulate. At the same time we shall derive an expression for the distribution of the first exit times valid for any age, and not only for $`t_a=0`$ and $`t_a=\mathrm{}`$, which denote observation taking place at the moment of birth and death of the dynamic process under study, respectively. This discussion is based on the following dynamic model. A variable $`x`$ moves in the interval $`I[0,1]`$ according to the following prescription

$$\frac{dx}{dt}=\lambda x(t)^z.$$
(34)
When it reaches the point $`x=1`$ it is injected back with uniform distribution, thereby producing another laminar region. We assign alternated signs to the sequel of laminar regions. The sojourn in a laminar region with a given sign is equivalent to sojourn in one of the two states discussed in Section 2. The resulting waiting time distribution coincides with the one earlier called $`\psi (t)`$. The explicit expression for $`\psi (t)`$ is obtained by expressing the exit time $`t`$ as a function of the initial condition $`x_0`$. Then we have to assume

$$p(x_0)dx_0=\psi (t)dt.$$
(35)
The choice of a uniform of back injection process implies $`p(x_0)=1`$. This is the condition behind the CTRW of Montroll and Weiss. All this results into

$$\psi (t)=(\mu 1)\frac{T^{(\mu 1)}}{(T+t)^\mu },$$
(36)
with

$$\mu =\frac{z}{(z1)}$$
(37)
and

$$T=\frac{1}{\lambda (z1)}.$$
(38)
We propose for the aging process a calculation procedure different from that adopted by Barkai barkai . As we shall see, this procedure yields the same results as those of Barkai. As done in Ref. rosa , we discuss the equation of motion for the probability density $`\rho (x,t)`$. This is given by

$$\frac{}{t}\rho (x,t)=\frac{}{x}(\lambda x^z\rho (x,t))+C(t).$$
(39)
The time evolution of the distribution function is described by the following formula

$$\rho (x,t)_{t_a}^t\frac{C(\tau +t_a)}{[1+1/\alpha (t\tau )x^{1/\alpha }]^{1+\alpha }}𝑑\tau ,$$
(40)
where we have used the following new parameter:

$$\alpha =\frac{1}{(z1)}.$$
(41)
We have assumed that the flat initial distribution is assigned at time $`t=t_a`$. We have neglected the contribution $`[1+1/\alpha (t+t_a)^{(1/\alpha )}]^{(1+\alpha )}`$, necessary to recover at $`t=t_a`$, the flat initial distribution. On the basis of the approach that we shall detail hereby to derive the distribution of the first exit times, it is straightforward to check that this term yields negligible contributions. Using the property

$$\rho (x,t=0)dx=\psi _{t_a}(t)dt,$$
(42)
where $`\rho (x,t=0)`$ is given by Eq. (40), with $`t=0`$, we can evaluate the distribution function of the first exit times in general, and not only for $`t_a=0`$ and $`t_a=\mathrm{}`$. It is interesting to remark that the infinitely old distribution $`\psi _{\mathrm{}}(t)`$ of Section 4 A, is easily obtained by noticing that

$$\rho (x,\mathrm{})=(2z)/x^{z1}.$$
(43)
Using Eqs.(43) and (42), after some algebra, we rederive the distribution of Eq. (31).

We can also establish a connection with the case discussed by Barkai barkai , namely the case $`z>2`$. In the case of $`z>2\alpha <1`$, the long time behavior of the function $`C(t)`$ is

$$C(t)\frac{\mathrm{sin}(\alpha \pi )}{\alpha ^\alpha \pi }t^{(1\alpha )}.$$
(44)
We now use the prescription of Eq. (42). A close inspection of Eqs. (40) and (44) reveals that it is not possible to give a close formula for the function $`\psi _{t_a}(t)`$. Then, we must study the two cases $`t_at`$ and $`t_at`$, separately.

* $`t_at`$:

$`\rho (x,t=0)`$ $``$ $`{\displaystyle _{t_a}^0}{\displaystyle \frac{(\tau +t_a)^{(1\alpha )}}{[11/\alpha \tau x^{1/\alpha }]^{1+\alpha }}}𝑑\tau `$ (45)
$`=`$ $`t_a^1{\displaystyle _1^0}{\displaystyle \frac{(1+y)^{(1\alpha )}}{[1/t_ay/\alpha x^{1/\alpha }]^{1+\alpha }}}𝑑y`$
$`=`$ $`t_a^1{\displaystyle _0^1}{\displaystyle \frac{(1y)^{(1\alpha )}}{[1/t_a+y/\alpha x^{1/\alpha }]^{1+\alpha }}}𝑑y`$
$``$ $`t_a^1{\displaystyle _0^1}{\displaystyle \frac{1}{[1/t_a+y/\alpha x^{1/\alpha }]^{1+\alpha }}}𝑑y`$
$`=`$ $`t_a^\alpha {\displaystyle _0^1}{\displaystyle \frac{1}{[1+y/\alpha t_ax^{1/\alpha }]^{1+\alpha }}}𝑑y`$
$``$ $`{\displaystyle \frac{t_a^{(1\alpha )}}{x^{1/\alpha }}}.`$

In Eq. (45), we have carried out the operation $`lim_{t_a\mathrm{}}`$, while keeping $`x0`$ fixed. The next step is based on the use of Eq. (42), of $`x=(1+t/\alpha )^\alpha `$ and of $`\left|\frac{dx}{dt}\right|=(1+t/\alpha )^{(\alpha +1)}`$. This allows us to obtain

$$\psi _{t_a}(t)|_{t_at}\frac{t_a^{(1\alpha )}}{t^\alpha }.$$
(46)
* $`t_at`$:

$`\rho (x,t=0)`$ $``$ $`{\displaystyle _{t_a}^0}{\displaystyle \frac{(\tau +t_a)^{(1\alpha )}}{[11/\alpha \tau x^{1/\alpha }]^{1+\alpha }}}𝑑\tau `$ (47)
$``$ $`{\displaystyle _{t_a}^0}(\tau +t_a)^{(1\alpha )}𝑑\tau `$
$`=`$ $`t_a^\alpha .`$
Thus, in this condition we have

$$\psi _{t_a}(t)|_{t_at}\frac{t_a^\alpha }{t^{1+\alpha }}.$$
(48)

Let us note that Eqs. (46) and (48) correctly reproduce the behavior predicted by Ref. barkai ) in the corresponding limits by means of the formula:

$$\psi _{t_a}(t)\frac{\mathrm{sin}(\alpha \pi )}{\pi }\frac{t_a^\alpha }{t^\alpha (t+t_a)},$$
(49)
where $`\alpha \frac{1}{z1}`$. In conclusion, our way of proceeding fits the conclusions of Barkai barkai . As earlier pointed out, in the condition $`\mu <2`$ the Onsager principle is always violated. This is expected, since it is known bettin  that the non-stationary condition $`\mu <2`$ is incompatible with the existence of the stationary correlation function $`\mathrm{\Phi }_\xi (t)`$.

## V The Stationary Master Equation

In this section we prove that the GME can be made equivalent to the CTRW, with no restriction to the Poisson statistics, if we set the stationary condition in both of them. The deep difference between stationary and non-stationary condition has been discussed years ago by Zumofen and Klafter in two pioneering papers a ; d . This means that we have to compare the GME of Eq. (1) to the stationary CTRW, which yields

$$\widehat{𝐩}(u)=(\widehat{\mathrm{\Psi }}_{\mathrm{}}(u)I+\frac{\widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)M}{1\widehat{\psi }(u)M})𝐩(0),$$
(50)
where, according to the stationary prescription chemphys :

$$\psi _{\mathrm{}}(t)=\frac{1}{\tau }_t^{\mathrm{}}𝑑t^{}\psi (t^{}),$$
(51)
as also explained in Section 4A, and $`\tau `$ is used again to denote the mean waiting time of the distribution $`\psi (t)`$. For $`M`$ we make the same choice as in Eq. (7). For calculation convenience, we adopt for $`K`$ the more general symmetric form:

$$K=\left(\begin{array}{cc}x& y\\ y& z\end{array}\right).$$
(52)
From (1) we get:

$$\widehat{𝐩}(u)=\frac{1}{uI+\widehat{\mathrm{\Phi }}(u)K}𝐩(0).$$
(53)
By equating (53) and (50) we derive

$$\frac{1}{uI+\widehat{\mathrm{\Phi }}(u)K}=\widehat{\mathrm{\Psi }}_{\mathrm{}}(u)I+\widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)M\underset{n=0}{\overset{\mathrm{}}{}}\widehat{\psi }(u)^nM^n.$$
(54)
By multiplying both sides by $`uI+\widehat{\mathrm{\Phi }}(u)K`$ and using the properties $`M^{2n}=I`$ and $`M^{2n+1}=M`$, we obtain:

$$\begin{array}{cc}\hfill I=& \widehat{\mathrm{\Psi }}_{\mathrm{}}(u)(uI+\widehat{\mathrm{\Phi }}(u)K)+\widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)(\frac{u}{1\widehat{\psi }(u)^2}M\hfill \\ & +\frac{\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}\mathrm{𝑀𝐾}+\frac{\widehat{\psi }(u)u}{1\widehat{\psi }(u)^2}I+\frac{\widehat{\psi }(u)\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}K).\hfill \end{array}$$
(55)
With some algebra we extract from Eq. (55) the following set of equations for $`x,y`$ and $`z`$:

$$\begin{array}{cc}\hfill 1=& \widehat{\mathrm{\Psi }}_{\mathrm{}}(u)(u+\widehat{\mathrm{\Phi }}(u)x)+\hfill \\ & \widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)\left(\frac{\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}y+\frac{\widehat{\psi }(u)u}{1\widehat{\psi }(u)^2}+\frac{\widehat{\psi }(u)\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}x\right)\hfill \\ \hfill 0=& \widehat{\mathrm{\Psi }}_{\mathrm{}}(u)\widehat{\mathrm{\Phi }}(u)y+\hfill \\ & \widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)\left(\frac{\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}z+\frac{u}{1\widehat{\psi }(u)^2}+\frac{\widehat{\psi }(u)\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}y\right)\hfill \\ \hfill 0=& \widehat{\mathrm{\Psi }}_{\mathrm{}}(u)\widehat{\mathrm{\Phi }}(u)y+\hfill \\ & \widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)\left(\frac{\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}x+\frac{u}{1\widehat{\psi }(u)^2}+\frac{\widehat{\psi }(u)\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}y\right)\hfill \\ \hfill 1=& \widehat{\mathrm{\Psi }}_{\mathrm{}}(u)(u+\widehat{\mathrm{\Phi }}(u)z)+\hfill \\ & \widehat{\psi }_{\mathrm{}}(u)\widehat{\mathrm{\Psi }}(u)\left(\frac{\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}y+\frac{\widehat{\psi }(u)u}{1\widehat{\psi }(u)^2}+\frac{\widehat{\psi }(u)\widehat{\mathrm{\Phi }}(u)}{1\widehat{\psi }(u)^2}z\right)\hfill \end{array}$$

The solution of this set of equations is given by $`x=z,y=x`$ and

$$x=\frac{u(1\widehat{\psi }(u))}{\widehat{\mathrm{\Phi }}(u)(2+u\tau +(2+u\tau )\widehat{\psi }(u))}.$$
(56)
So, the GME and the stationary CTRW are compatible if $`K`$ has the same form as that of Eq. (8) and

$$\widehat{\mathrm{\Phi }}(u)=\frac{u(1\widehat{\psi }(u))}{(2+u\tau +(2+u\tau )\widehat{\psi }(u))}.$$
(57)
Using the Onsager’s principle in the form of Eq. (10) and consequently the time convoluted equation of Eq. (11), we obtain

$$\widehat{\mathrm{\Phi }}_\xi (u)=\frac{1}{u+2\widehat{\mathrm{\Phi }}(u)}.$$
(58)
Finally, using for $`\widehat{\mathrm{\Phi }}(u)`$ the expression of Eq. (57) and for $`\widehat{\mathrm{\Phi }}_\xi (u)`$ the expression afforded by the renewal theory, we get an identity in $`\widehat{\psi }(u)`$ (the Laplace transform of the waiting times distribution). This means that the aged CTRW is totally equivalent to the GME with the memory kernel given by Eq. (57).

## VI Concluding remarks

This paper affords a significant contribution to understanding the role of aging in systems departing from Poisson statistics. The Onsager regression principle rests on the second principle of thermodynamics, and it applies to aged systems. However, non-Poissonian statistics might imply a very slow regression to equilibrium, so as to create conditions, correctly depicted by the Montroll-Weiss CTRW, but incompatible with the Onsager regression postulate. This might be a condition different from the out of equilibrium thermodynamics, and represent instead a state of matter intermediate between the dynamic and the thermodynamic condition. This new condition has been recently proposed to explain the emergency of life on earth and has been denoted with the appealing term of *living state of matter* livingstateofmatter .

The master equation of Section 5 corresponds to the aged systems where the Onsager principle holds true. The Montroll-Weiss CTRW and the corresponding master equation refers to the dynamics of a young system. From a technical point of view one might wonder if it is possible to produce a master equation, or a CTRW, reflecting correctly the process of aging. The young and the old sequence reflect two different initial conditions of the system. Moving from the young condition the system will eventually relax to the equilibrium, but this relaxation will not follow the correlation function behavior.

Eventually, the system will settle at equilibrium. However, to establish, at least ideally, the validity of the Onsager regression principle, we have to set some suitable out of equilibrium initial conditions. If we do select initial conditions that do not fit the stationary prescription a ; d , the correlation function $`\mathrm{\Phi }_\xi (t)`$ cannot be defined, and consequently we are forced to determine the memory kernel of the GME of Eq. (1) using the non-stationary prescription of Eq. (9). Neverthless, this equation will produce an exact time evolution. Let $`𝐩(0)`$ evolve in time till it reaches the new condition $`𝐩(T)`$. If we adopt this as a new initial condition, namely we set the new origin of time at $`t=T`$, and we apply again Eq. (1) to predict the ensuing time evolution, we get a wrong prediction. In fact, the GME equivalent to the CTRW of Montroll and Weiss is not a *bona fide* master equation fox . We should build up a new GME equation with the distribution of the first jumps given by $`\psi _{t_a=T}(t)`$. If, on the contrary, we select an initial out of equilibrium condition, as done in Section 3, which is compatible with the bath being at equilibrium, then the Onsager principle is fulfilled, the GME with the kernel of Eq. (57) is a *bona fide* master equation and can be used with any kind of initial conditions concerning $`𝐩`$.

To derive the results here illustrated, we significantly benefitted from the important work of Refs. sokolovandmetzler ; barkai . However, we want to point out that the material derivative of Ref. sokolovandmetzler , although more advantageous than the GME to study the influence of an external force field, is equivalent to the GME of Ref. kenkre , with the condition of Eq. (6) rather than the stationary condition of Eq. (57), and, consequently, as proved in this paper,it violates the Onsager principle. As to the work of Barkai on aging barkai , we would like to note that it is confined to the case $`\mu <2`$, a condition incompatible with the stationary master equation found in this paper, namely the master equation with the memory kernel of Eq. (57). Barkai reach the conclusion that the aging effect for $`\mu >2`$ are “confined to a certain time window”. We think that this paper has the effect of suggesting the plausible conjecture that actually this time window is of infinite size. Let us see why. We notice that in a sense the memory kernel of Eq. (57) affords a way to measure the aging process. In fact, if the waiting time distribution $`\psi (t)`$ is exponential, the memory kernel is equivalent to that produced by a bath that would regress to equilibrium instantanesously, with no aging effect. The aging effect, from within our perspective, can also be interpreted as a regression to equilibrium of a system of interest, under the influence of a bath that is regressing, too, to equilibrium. In the Poisson case, the bath regression to equilibrium would be instantaneous, thereby yielding no aging. In the case of $`\mu <\mathrm{}`$ we have a memory lasting forever. This is evident when $`\mu <3`$, since in this case the adoption of Eq. (57) proves that the Laplace transform of the memory kernel, $`\widehat{\mathrm{\Phi }}(u)`$, tends to vanish for $`u0`$. This is an indication of infinitely persistent oscillations that are clearly incompatible with the Markov approximation.

Apparently, the condition $`\mu >3`$ seems to be compatible with the Markov approximation, given the fact that in that case, $`\widehat{\mathrm{\Phi }}(u)`$ tends to a finite value for $`u`$ tending to $`0`$. We note that, while the analytical expression of $`\mathrm{\Phi }(t)`$ is at the moment not yet known to us, the analytical form of $`\mathrm{\Phi }_\xi (t)`$ is under our total control, and it is a an inverse power law with index $`\beta =\mu 2`$. The free diffusion process resulting from the fluctuations of $`\xi (t)`$ is known. The pioneering work of Ref. klafter  proves that even in the case $`\mu >3`$ the Gaussian condition is not exactly realized. The central part of the diffusion pdf is the ordinary Gaussian diffusion, but at large distances slow tails with an inverse power law appear. These tails become weaker and weaker upon increase of $`\mu `$ and might be annihilated by even weak external fluctuations, with the ensuing annihilation of aging. However, there are deep physical motivations to fully understand the theoretical issues generated by the aging process. A relevant example of this theoretical request is given by the physics of glassy systems, which are characterized by aging and consequently by a violation of the ordinary fluctuation-dissipation condition, which is intimately related to the validity of the Onsager principle crisanti . In other words, we have to do further research work to understand how robust the non-Poissonian nature of $`\psi (t)`$ is, and to establish the intensity of environemental fluctuations that are expected to produce a truncation of the inverse power law nature of $`\psi (t)`$. This can be done along the lines of Ref. elena . This research work can be carried at a much deeper level, as a consequence of the fundamental result of Eq. (57). In fact, this equation defines the memory kernel $`\mathrm{\Phi }(t)`$ of the GME, in terms of the waiting function $`\psi (t)`$. We assume that the function $`\psi (t)`$ can be derived from the experimental observation as in the case of the blinking quantum dots of Ref. chemphys . Thus, with the help of Eq. (57) we can draw from the experimental function $`\psi (t)`$ information on the bath responsible for the intermittent character of the process under study. This information will help the foundation of a dynamic model, and consequently will make possible to discuss the robustness against noise of the intermittent process under study, from a perspective more realistic than in the earlier work, for instance, the work of Ref. elena .

Let us summarize the main results of this paper. The functions $`\mathrm{\Phi }(t)`$ and $`\psi (t)`$ are the main mathematical properties behind the GME and CTRW, respectively, the former affording indications on the complex dynamics of the irrelevant variables and the latter illustrating the nature of the intermittent process that can be observed as a result of a special experimental detection chemphys . The deviation of $`\psi (t)`$ from the exponential condition is equivalent to forcing $`\mathrm{\Phi }(t)`$ depart dramatically from the white noise condition. The non-Poisson nature of $`\psi (t)`$ generates memory infinitely extended in time, and aging with no finite time scale as well, no matter how big the scaling parameter $`\mu `$, provided that $`\mu <\mathrm{}`$. In fact, if a finite time scale for both aging and memory existed, at a much larger time scale, the function $`\mathrm{\Phi }(t)`$ would be equivalent to a delta of Dirac, the bath fluctuations would be white noise, the bath would be at equilibrium, and, of course, $`\psi (t)`$ would be exponential. This is equivalent to proving that the non-exponential nature of $`\psi (t)`$ yields aging effects that are not limited to the non-stationary case $`\mu <2`$, even if in this case these aging effects seem to be more natural. The condition $`\mu >3`$ is compatible with the stationary condition. However, a bath generating this condition, through the key result of Eq. (57), is characterized by dynamic properties with memory infinitely extended in time.

###### Acknowledgements.

We thank Professor Bill Schieve for his remark that we are using a generalized form of Onsager’s principle. We thankfully acknowledge the financial support of ARO, through Grant DAAD19-2-0037.