text:
# MS-TP-03-27hep-lat/0311032 Chiral perturbation theory for lattice QCD with a twisted mass term

(November 21, 2003
(revised: February 26, 2004))

## Abstract

Quantum Chromodynamics on a lattice with Wilson fermions and a chirally twisted mass term for two degenerate quark flavours is considered in the framework of chiral perturbation theory. The pion mass and decay constant are calculated in next-to-leading order including terms linear in the lattice spacing $`a`$.

Numerical simulations of Quantum Chromodynamics with dynamical quarks are facing a difficulty with light quark masses. Due to a severe slowing down of the common algorithms for small quark masses it is presently not possible to implement realistic values for the masses of the u- and d-quarks. In this context chiral perturbation theory  is a very useful tool. It amounts to an expansion around the chiral limit where the quarks are massless. Meson masses and other physical quantities are then expanded in powers of quark masses (modified by logarithms). In its range of applicability, chiral perturbation theory can be used to extrapolate the results of Monte Carlo simulations of QCD into the region of small quark masses. On the other hand, numerical simulations of lattice QCD can provide information about low energy parameters of chiral perturbation theory, the Gasser-Leutwyler coefficients, see  for a review.

When chiral perturbation theory is applied to lattice QCD, the lattice spacing $`a`$ enters as a second expansion parameter. Correspondingly, the low energy effective Lagrangean contains additional terms proportional to powers of $`a`$, representing the lattice artifacts . Physical quantities appear in double expansions in quark masses and in $`a`$. A comparison of numerical results with chiral perturbation theory in next-to-leading order has been made in .

An attractive framework for numerical simulations of QCD with Wilson fermions is the lattice formulation with a chirally twisted quark mass matrix . It allows to avoid problems with quark zero modes and simplifies the renormalization procedure. In view of Monte Carlo calculations in twisted mass lattice QCD it is desirable to extend chiral perturbation theory to this case.

In this letter we consider lattice QCD with $`N_f=2`$ quark flavours and a twisted mass term. For simplicity we restrict ourselves to the case of degenerate quark masses ($`m_u=m_d=m`$). The pion mass and the pion decay constant are calculated in chiral perturbation theory in next-to-leading (one-loop) order, including lattice terms linear in the lattice spacing $`a`$. A new feature of the twisted mass case is the fact that the minimum of the effective Lagrangean is shifted by an amount proportional to $`a`$ and the expansion has to be rearranged correspondingly.

In lattice QCD a chirally twisted mass term can be introduced in the form $`\overline{q}M(\omega )q`$ with

$$M(\omega )=m^{}\mathrm{e}^{\mathrm{i}\omega \gamma _5\tau _3}=m^{}\mathrm{cos}(\omega )\mathbf{\hspace{0.17em}1}+\mathrm{i}m^{}\mathrm{sin}(\omega )\gamma _5\tau _3m+\mathrm{i}\mu \gamma _5\tau _3,$$
(1)
where

$$m=m^{}\mathrm{cos}(\omega ),\mu =m^{}\mathrm{sin}(\omega ),$$
(2)
and $`\tau _b`$ are Pauli matrices. In the continuum the chiral rotation of the mass term can be removed by a chiral transformation of the quark fields according to

$$q=\mathrm{e}^{\mathrm{i}\omega \gamma _5\tau _3/2}q^{},$$
(3)
which leaves the kinetic term invariant. On the lattice, however, it is not possible to remove the chiral twist due to the presence of lattice artifacts.

Chiral perturbation theory is based on the low energy effective Lagrangean, which describes the dynamics of pions in terms of SU(2)-valued matrices

$$U(x)=\mathrm{exp}\left(\frac{\mathrm{i}}{F_0}\pi _b(x)\tau _b\right),$$
(4)
where $`\pi _b(x)`$ are the pion fields. In the continuum it is given to leading order by

$$_2=_{2,kin}+_m=\frac{F_0^2}{4}\text{Tr}\left(_\mu U^{}^\mu U\right)\frac{F_0^2}{4}\text{Tr}\left(\chi U^{}+U\chi ^{}\right).$$
(5)
The chiral symmetry breaking mass term contains the matrix

$$\chi =2B_0M=2B_0m\mathbf{\hspace{0.17em}1}\chi _0\mathbf{\hspace{0.17em}1},$$
(6)
with

$$B_0=\frac{1}{F_0^2}0|\overline{u}u|0=\frac{1}{F_0^2}0|\overline{d}d|0.$$
(7)
To lowest order in $`m`$ the pion mass is given by

$$m_\pi ^2=\chi _0=2B_0m.$$
(8)
Lattice QCD can be treated in chiral perturbation theory by adding terms which represent the lattice artifacts. The chiral Lagrangean to lowest order is extended to

$$_2=_{2,kin}+_m\frac{F_0^2}{4}\text{Tr}\left(\rho U^{}+U\rho ^{}\right),$$
(9)
where

$$\rho =2W_0a\mathbf{\hspace{0.17em}1}$$
(10)
and the coefficient $`W_0`$ parametrizes the lattice artifacts.

Now let us consider the effective Lagrangean for QCD with a twisted mass term. The twisting of the mass term

$$\chi \chi (\omega )=2B_0\mathrm{e}^{\mathrm{i}\omega \tau _3/2}M\mathrm{e}^{\mathrm{i}\omega \tau _3/2}=2B_0m\mathrm{e}^{\mathrm{i}\omega \tau _3}$$
(11)
can be undone in the kinetic and mass terms by a chiral transformation of the axial type

$$U=\mathrm{e}^{\mathrm{i}\omega \tau _3/2}U^{}\mathrm{e}^{\mathrm{i}\omega \tau _3/2}.$$
(12)
The lattice term, however, is transformed to

$$\text{Tr}\left(\rho (\omega )U^{}+U^{}\rho (\omega )^{}\right)$$
(13)
with the matrix

$$\rho (\omega )=\mathrm{e}^{\mathrm{i}\omega \tau _3}\rho .$$
(14)
In this way the twist can be transferred to the lattice term. Omitting the primes on the chiral field $`U`$, the leading order effective Lagrangean reads

$$_2=\frac{F_0^2}{4}\text{Tr}\left(_\mu U^{}^\mu U\right)\frac{F_0^2}{4}\text{Tr}\left((\chi +\rho (\omega ))U^{}+U(\chi +\rho (\omega ))^{}\right)$$
(15)
with

$$\chi =2B_0m^{}\mathbf{\hspace{0.17em}1}.$$
(16)
In next-to-leading order the following terms contribute :

$`_4`$ $`=`$ $`{\displaystyle \frac{F_0^2}{4}}\text{Tr}\left(_\mu U_\mu U^{}\right){\displaystyle \frac{F_0^2}{4}}\text{Tr}\left(\chi U^{}+U\chi ^{}\right){\displaystyle \frac{F_0^2}{4}}\text{Tr}\left(\rho (\omega )U^{}+U\rho (\omega )^{}\right)`$
$`L_1\left[\text{Tr}\left(_\mu U_\mu U^{}\right)\right]^2L_2\text{Tr}\left(_\mu U_\nu U^{}\right)\text{Tr}\left(_\mu U_\nu U^{}\right)`$
$`L_3\text{Tr}\left(\left[_\mu U_\mu U^{}\right]^2\right)+L_4\text{Tr}\left(_\mu U_\mu U^{}\right)\text{Tr}\left(\chi U^{}+U\chi ^{}\right)`$
$`+W_4\text{Tr}\left(_\mu U_\mu U^{}\right)\text{Tr}\left(\rho (\omega )U^{}+U\rho (\omega )^{}\right)+L_5\text{Tr}\left(_\mu U_\mu U^{}\left[\chi U^{}+U\chi ^{}\right]\right)`$
$`+W_5\text{Tr}\left(_\mu U_\mu U^{}\left[\rho (\omega )U^{}+U\rho (\omega )^{}\right]\right)L_6\left[\text{Tr}\left(\chi U^{}+U\chi ^{}\right)\right]^2`$
$`W_6\text{Tr}\left(\chi U^{}+U\chi ^{}\right)\text{Tr}\left(\rho (\omega )U^{}+U\rho (\omega )^{}\right)L_7\left[\text{Tr}\left(\chi U^{}U\chi ^{}\right)\right]^2`$
$`W_7\text{Tr}\left(\chi U^{}U\chi ^{}\right)\text{Tr}\left(\rho (\omega )U^{}U\rho (\omega )^{}\right)L_8\text{Tr}\left(\chi U^{}\chi U^{}+U\chi ^{}U\chi ^{}\right)`$
$`W_8\text{Tr}\left(\chi U^{}\rho (\omega )U^{}+U\rho (\omega )^{}U\chi ^{}\right)+𝒪(a^2).`$
Spurion analysis shows that the twisting of the mass term does not produce further terms in $`_4`$, but amounts to replacing $`\rho `$ by $`\rho (\omega )`$.

In contrast to the untwisted case, the minimum of the effective action is not located at vanishing fields $`\pi _b=0`$, corresponding to $`U=\mathrm{𝟏}`$, but at the point

$$\stackrel{~}{\pi }_3=F_0\frac{W_0a}{B_0m^{}}\mathrm{sin}\omega \left[1\frac{8\chi _0}{F_0^2}(4L_6+2L_82W_6W_8)\right]+𝒪(a^2),\stackrel{~}{\pi }_1=0,\stackrel{~}{\pi }_2=0.$$
(18)
Therefore a shift

$$\pi _3=\pi _3^{}+\stackrel{~}{\pi }_3$$
(19)
has to be performed before expanding the effective action.

Using the resulting expression we calculated the pion propagator in next-to-leading order. From it the pion mass and wave function renormalization are obtained. For the pion mass we get

$`m_\pi ^2`$ $`=`$ $`\chi _0+\rho _0+8{\displaystyle \frac{\chi _0^2}{F_0^2}}(4L_6^r+2L_8^r2L_4^rL_5^r)`$ (20)
$`+8{\displaystyle \frac{\chi _0\rho _0}{F_0^2}}(4W_6^r+2W_8^r2W_4^rW_5^r2L_4^rL_5^r)`$
$`+{\displaystyle \frac{(\chi _0+\rho _0)^2}{32\pi ^2F_0^2}}\mathrm{ln}\left({\displaystyle \frac{\chi _0+\rho _0}{\mathrm{\Lambda }^2}}\right),`$
where

$$\chi _0=2B_0m^{},\rho _0=2W_0a\mathrm{cos}\omega ,$$
(21)
$`L_k^r`$ are the renormalized chiral parameters and $`\mathrm{\Lambda }`$ is the renormalization scale. The wavefunction renormalization in dimensional regularization is given by

$`Z`$ $`=`$ $`1{\displaystyle \frac{8}{F_0^2}}\chi _0(2L_4+L_5){\displaystyle \frac{8}{F_0^2}}\rho _0(2W_4+W_5)`$ (22)
$`{\displaystyle \frac{1}{24\pi ^2F_0^2}}(\chi _0+\rho _0)\left[{\displaystyle \frac{2}{ϵ}}+\mathrm{ln}(4\pi )\gamma +1\mathrm{ln}\left({\displaystyle \frac{\chi _0+\rho _0}{\mathrm{\Lambda }^2}}\right)\right].`$
Setting the twist angle $`\omega `$ to zero, the expression for $`m_\pi ^2`$ is consistent with the result of . Apart from numerical factors which depend on $`N_f`$, to this order of chiral perturbation theory the chiral twist amounts to adding the factor $`\mathrm{cos}\omega `$ to $`\rho `$. It should be noted, however, that this rule will not apply to higher orders, where the above mentioned shift in the pion fields introduces new vertices.

In the case of maximal twist, $`\omega =\pi /2`$, the lattice artifacts vanish to the order of the expansion considered here, as has been observed for lattice QCD in general in .

The second physical quantity we calculated is the pion decay constant $`F_\pi `$, given by

$$0|J_A^{\mu ,a}|\pi _b(p)=\mathrm{i}F_\pi p^\mu \delta _{ab},$$
(23)
where $`J_A`$ is the axial current. The Noether procedure yields an expression for the axial current, which depends on the coefficients $`L_k`$ and $`W_k`$ and contains products of the pion fields. A one-loop calculation gives

$$F_\pi =F_0\left(1+\frac{4}{F_0^2}[\chi _0(2L_4^r+L_5^r)+\rho _0(2W_4^r+W_5^r)]\frac{1}{16\pi ^2F_0^2}\left[(\chi _0+\rho _0)\mathrm{ln}\frac{\chi _0+\rho _0}{\mathrm{\Lambda }^2}\right]\right).$$
(24)
In the limit $`a0`$ this coincides with the result of .

The formulae given here can be used in the analysis of numerical results from unquenched simulations of twisted mass lattice QCD with two light flavours. A useful technique in numerical simulations of QCD is partial quenching of quarks. In view of applications to partially quenched simulations it suggests itself to generalize the results above to partially quenched chiral perturbation theory .

We thank I. Montvay, E. Scholz and F. Farchioni for discussions.


format: text

text:
# Single pulses from PSR B1641–45

## 1 Introduction

The integrated pulse profile from pulsars is one of the key characteristics of the emission process. Although single pulses vary widely both in intensity and shape, integrating a few thousand pulses produces a highly stable integrated profile. The integrated profile therefore yields clues as to the permanent structures associated with the pulsar magnetosphere. Kramer et al. (1994) showed that the integrated profile can be decomposed into a small number of Gaussian components. These components represent a 2-D Gaussian region of emission from a fixed height in the patchy pulsar beam \[Lyne & Manchester 1988\]. In single pulses, some fraction of the component is emitting when the beam sweeps the line of sight and therefore the single pulses tend to have narrower features (sub-pulses) than those of the integrated profile. The sub-pulse features may also rotate about the magnetic axis, producing the drifting sub-pulse phenomenology \[Deshpande & Rankin 1999\]. Integrated profiles show high levels of linear polarization and moderate circular polarization. Single pulses tend to be more highly polarized than the integrated profile. The smooth position angle swing in integrated profiles as a function of pulse longitude can occasionally be broken by 90° jumps which shows that the emission is present in two orthogonal modes of radiation. There has been considerable debate in the literature over whether these two orthogonal modes are emitted simultaneously or disjointly. Recently, McKinnon & Stinebring (1998) and Karastergiou et al. (2002) have shown that simultaneously occuring modes nicely explain many aspects of single pulse data.

PSR B1641$``$45 is one of the brightest pulsars in the sky at 1.4 GHz. It has a pulse period of 455 ms and its pulse profile consists of a single component preceeded by a low amplitude ‘precursor’. Previous observations of this pulsar at a frequency of 1.6 GHz by Manchester, Hamilton & McCulloch (1980) showed that it had moderate linear and circular polarization and that the leading ‘precursor’ had a high degree of linear polarization orthogonal to the polarization of the main pulse emission. Rankin (1990) and van Ommen et al. (1997) classify the pulsar as a core single even though it has a relatively shallow position angle swing across the pulse.

## 2 Observations and Data Reduction

As part of a project investigating single pulses from a large sample of pulsars, PSR B1641$``$45 was observed on 2000 March 16 and 17 using the 64-m Parkes radio telescope. The centre frequency, $`\nu `$, of the observations was 1413 MHz; at this frequency the system equivalent flux density is 26 Jy. The receiver consists of two orthogonal feeds sensitive to linear polarization. The signals are down-converted and amplified before being passed into the backend. The backend, CPSR, is an enhanced version of the Caltech Baseband Recorder (Jenet et al. 1997). It consists of an analogue dual-channel down-converter and digitizer which yields 2-bit quadrature samples at 20 MHz. The data stream is written to DLT for subsequent off-line processing allowing all 4 Stokes parameters to be computed. On the two occasions, 3102 and 3718 single pulses were recorded. Before each observation, a 90-s observation of a pulsed signal, directly injected into the receiver at a 45° angle to the feed, is made. This enables the gain of each channel and the phase difference between the channels as a function of frequency to be calibrated. Observations of the flux calibrator Hydra A were made which allows absolute fluxes to be obtained.

The data were processed off-line using a workstation cluster at the Swinburne Supercomputer Centre. Data reduction involves coherent de-dispersion (Hankins & Rickett 1975) and includes quantization error corrections as described by Jenet & Anderson (1998). The data are folded at the apparent topocentric period of the pulsar and the full Stokes profiles for each pulse are written to disk. Flux calibration and instrumental calibration are then carried out using information contained in the observation of the pulsed (calibration) signal. The data in each frequency channel is corrected for the rotation measure of the pulsar and all the frequency channels are then summed to produce the final profile. The final profiles consist of 1024 time-bins per pulse period for an effective time resolution of 0.44 ms.

## 3 Integrated Profile

The integrated profile of PSR B1641$``$45 is shown in Figure 1. The continuum flux density is 430 mJy and the peak of the profile has a flux density of 21 Jy. The pulsar shows a moderate amount of linear polarization ($``$20%) and a small amount of circular polarization ($``$6%), which is predominantly negative throughout. The position angle swings by about 50° throughout the main pulse. There is a ‘pedestal’ or ‘pre-cursor’ at the leading edge of the profile and the polarization of this feature is orthogonal to the rest of the pulse. The peak of the linear polarization occurs earlier than the peak in total intensity. These results are consistent with the observations of Manchester et al. (1980). There is scatter broadening of the pulse profile at this observing frequency although it is difficult to measure without knowing the underlying pulsar profile. As part of a separate study of this pulsar (Rickett et al. In preparation), we measured a pulse broadening time of 56 ms at an observing frequency of 660 MHz. This is an improvement over the only value currently in the literature, 40$`\pm `$10 ms at 750 MHz by Komesaroff et al. (1973). Extrapolating to $`\nu `$=1413 MHz by assuming that the broadening time scales as $`\nu ^{4.4}`$ yields a scattering time of $``$2.0 ms at the observing frequency used here. This time is small compared to the duration of the pulse, but likely manifests itself as the tail of the emission at times beyond 13 ms in Figure 1.

After accounting for the orthogonal mode jump, the rotating vector model (RVM) can be fitted to the position angle swing. As is commonly the case, the angles $`\alpha `$ (between the rotation and magnetic axes) and $`\beta `$ (the offset of the line of sight from the magnetic axis) are not well constrained, mostly because emission only occurs from a small fraction of pulse longitude. I adopt the parameters $`\alpha =30\mathrm{°}`$, $`\beta =8\mathrm{°}`$ with the magnetic pole crossing located at phase 0.931, just prior to the pulse peak. These parameters provide an adequate fit to the data, and their exact values are not crucial to the discussion below.

Kramer et al. (1994) demonstrated the method of decomposing the total intensity profiles of pulsars into Gaussian components and showed that typically fewer than five were necessary in the majority of cases. Applying this method to PSR B1641$``$45, I find that four Gaussian components provide the best fit to the data. The fit is excellent apart from the the trailing edge of the pulsar where scatter broadening of the profile is clearly seen. I define a Gaussian by

$$y=Ae^{\frac{1}{2}(\frac{x\overline{x}}{\sigma })^2}$$
(1)
and list the parameters of the fitted Gaussians in Table 1. The fitted components are shown in Figure 2. The smallest amplitude component is necessary to fit the precursor. The location of the second Gaussian is close to the peak of the linear polarization. I therefore attempted to fit both the linear and the circular polarization by varying the fractional polarization in each of the four components. The results are listed in Table 1. Note that the linear polarization of the precursor is negative which implies that it is orthogonally polarized to the rest of the emission. The position angle at a given longitude is computed from the parameters of the RVM fit given above. The appropriate values of Stokes $`Q`$ and $`U`$ are then computed from the position angle for each Gaussian. These are then summed to produce the integrated Stokes $`Q`$ and $`U`$ which are then re-converted back to position angle. The profile of the ‘simulated’ pulsar is shown in Fig. 3. The simulated profile can be compared directly with Fig. 1; in particular the simulated profile nicely reproduces most of the main features of the real profile, including the location of the orthogonal mode jump which occurs when the amplitude of Gaussian 1 drops below that of Gaussian 2. Also, both the linear and circular polarizations match well with the true values.

## 4 Single Pulses

### 4.1 Total Intensity

The single pulses from this pulsar are rather featureless. They generally show a similar morphology to that of the integrated profile and there is no evidence for microstructure or drifting sub-pulses as is the case for many other pulsars. The precursor component seen in the integrated profile is rarely seen in the single pulses as the noise dominates its emission.

Given the similarity between the single pulses and the integrated profile, I attempted to perform Gaussian decomposition on the single pulses. Most recently, Gupta & Gangadhara (2003) have applied a ‘window-thresholding’ technique to single pulse data. This allowed them to detect extra components in the pulse profile which they then Gaussian fitted to obtain their parameters. Here, I adopt a slightly different technique. The width and centroid of the gaussians are fixed to the values obtained from the fit to the integrated profile as it seemed apparent from visual inspection that the single pulses had the same pulse width as the integrated profile. However, the low amplitude Gaussian necessary to fit the precursor region of the integrated profile was not included as it is too weak to be seen in single pulses and thus could not be constrained in the fitting process. I therefore performed fitting to the single pulses with three free parameters (the amplitude of the three Gaussians listed as 2-4 in Table 1). This provided good fits (as judged by the $`\chi ^2`$ statistic) to the single pulse data in more than 95% of the cases. The majority of the failures came when the flux density of the pulsar was low. Often, the amplitude of Gaussian 2 was very low and the single pulses were just as adequately fit with two Gaussians. However, some pulses required a high amplitude for Gaussian 2.

The scatter broadening of the pulse is not taken into account in the fitting procedure in that symmetrical Gaussians rather than scatter broadened Gaussians are used. The scattering only affects the tail of the profile, where the signal to noise ratio is low. However, the goodness of fit is weighted towards the high signal to noise points (i.e. the pulse peak). The fitting procedure is therefore largely unaffected by the scattering. I return to the effects of scattering on the implications of the results in a later section. Figure 4 shows two examples of single pulses and the fitted Gaussians. The fits are excellent and show that the scatter broadening is not playing a significant role in determining the fitted parameters.

The question then arises as to the distribution of the amplitudes of the Gaussians. Figure 5 shows the distribution for Gaussian 4. The distribution is clearly log-normal with a mean amplitude of 16.7 Jy and a $`\sigma `$ of 0.27 (in the log). Similarly, Gaussian 3 is log-normally distributed with a mean of 38.7 Jy and a $`\sigma `$ of 0.13 (in the log). The statistics of Gaussian 2 are more problematical to determine as often this component is weak and hence dominated by noise. A reasonable fit is obtained with a log-normal distribution with a mean of 7.45 Jy and a $`\sigma `$ of 0.3 in the log. Note that the means obtained in the fitting process correspond nicely to the peak amplitudes obtained during the Gaussian decomposition process on the integrated profile. This gives confidence that the fitting procedure is working correctly. Remarkably, therefore, the single pulses from this pulsar can be reproduced by drawing three random numbers from the above distributions and summing the resultant Gaussians together.

### 4.2 Polarization

In order to examine the polarization of the single pulses, an analysis similar to that of Stinebring et al. (1984a,b) was performed. For the linear polarization, for example, a two dimensional grid is constructed with the pulse phase along the x-axis and 20 bins along the y-axis from 0 to 100%. For each pulse the appropriate cell is then updated by one depending on its linear polarization. The same procedure is carried out for the circular polarization and the position angle. The values in the cells are then converted to a linear grey-scale between white and black depending on the counts in that cell. Figure 6 shows the grey-scale representation of the data. The figure shows that although there is a large spread in the linear and circular polarization, rather few pulses have a net polarization greater than 50% apart from the leading edge of the profile. There is also little evidence of orthogonally polarised modes in the position angle distribution.

Figure 7 shows the polarization properties of the phase bin where the linear polarization is maximum. In this bin, one mode dominates virtually all the time, the position angle distribution is narrow around its mean value. The fractional linear polarization is high and the circular polarization shows large flux densities of either sign although there is a slight preference for negative values. This contributes to give a low net circular polarization.

Figure 8 shows the same information for a phase bin on the trailing edge of the profile. Both the total intensity and linearly polarized intensity are lower as is the fractional linear polarization. The PA distribution is much broader than seen in Fig. 7. This is partly due to increased uncertainity in the measurement of the PA (as the linear polarization is lower) and partly due to the presence of the orthogonal mode appearing near PA $`80\mathrm{°}`$. The circular polarization is more evenly distributed around zero.

## 5 Simulating Single Pulses

I have already demonstrated that replicating the single pulses in total intensity is a relatively simple task for this pulsar. It involves drawing three random numbers from a log-normal distribution. What further parameters need to be added in order to replicate the polarization of the single pulses?

As a starting point I use the model proposed by McKinnon & Stinebring (1998). In their model both polarization modes are emitted simultaneously and are 100% polarized. Then, the total power is the sum of the intensities of the two modes and the linear polarization is their difference. The amplitudes of the modes are highly correlated. They found that, under these assumptions, the single pulse data of Stinebring et al. (1994) could be replicated. McKinnon & Stinebring (1998) treated each phase bin in an independent fashion. In a further paper, they also showed that the integrated profile could be ‘split’ into the two orthogonal modes \[McKinnon & Stinebring 2000\]. This involved a simple mathematical manipulation based on the fractional linear polarization in each bin (see also Karastergiou et al. 2003).

In this instance, I assume that the Gaussians described in Section 3 are a coherent emitting unit, or patch, and that they emit radiation simultaneously in two, 100% polarized orthogonal modes. Ignoring for now the circular polarization, one can apply the technique of McKinnon & Stinebring (2000) to compute the mean amplitude of each of the modes for each of the Gaussians under these assumptions using the data from Table 1 and knowing the integrated fractional polarization. Furthermore, the distribution of the amplitudes of these Gaussians are also known from the fitting described above; they are all log-normal.

In order to simulate the polarization of the single pulses from PSR B1641$``$45, therefore, we now need to draw 8 random numbers (one pair for each of the Gaussians in the table) from a log-normal distribution with the mean and sigma as given in Table 2. Note that the random values for each mode are correlated according to the correlation coefficient given in the table \[McKinnon & Stinebring 1998\]. Also, the parameters for Gaussian 1 are not well constrained, but do not play a major role in the single pulse data as the flux density is low. The total intensity of each Gaussian is then the sum of the amplitudes of the two orthogonal modes and the linear polarization is the difference of their amplitudes. This accounts for the linear polarization but does not include any circularly polarized component. The mean circular polarization of each individual Gaussian is already known from the fits to the integrated profile. It is tempting to associate a handedness of circular polarization with a given orthogonal mode as has been proposed by e.g. Cordes, Rankin & Backer (1978). However, this seems be to true only at observing frequencies lower than $``$600 MHz. Recent observations of the Vela pulsar \[Kramer, Johnston & van Straten 2002\] and PSR B1133+16 \[Karastergiou, Johnston & Kramer 2003\] show that the situation is much more complex at frequencies above about 1 GHz. There is a lack of correlation between the handedness of circular polarization and the dominant orthogonal mode, but a rather large spread in the values of both hands of circular polarization. For each component (not each orthogonal mode), I therefore draw the circular polarization from a Gaussian distribution with a mean and standard deviation as listed in Table 2.

Single pulses simulated in such a way bear a strong resemblance to the true single pulse data. Figs 9 and 10 show the simulated data for the same phase bins as Figs 7 and 8. Comparison between the total intensity distributions show that the low flux density pulses are under-represented, however the high intensity data are well matched. The linear and circular polarizations are well replicated. The biggest discrepancy lies in the spread of the position angle distributions. In the simulated data, the spread is relatively narrow and the orthogonal modes are well separated. The spread in values is almost entirely due to receiver noise as I have assumed that the position angle is fixed for a given bin. In the real data the spread is much larger than a delta function convolved with the noise. This is similar to the results found by McKinnon & Stinebring (1998) and Karastergiou et al. (2002).

## 6 Discussion

I have shown that Gaussian decomposition of the integrated profile not only provides a good fit to the total intensity but can also fit the polarization profile, including the location of the orthogonal jump. In this pulsar, at least, these same Gaussians can be made to fit the single pulse data with only their amplitudes as a free parameter. Interestingly the distribution of these amplitudes is log-normal. Log-normal distributions are seen both in the time-averaged flux density of pulsars and in phase resolved statistics \[Cairns, Johnston & Das 2001, Cairns, Johnston & Das 2003\]. The implications are then that pulsar emission may be interpreted in terms of stochastic growth theory and that the emission mechanism is linear, involving either a direct linear instability or linear mode conversion of non-escaping waves driven by a linear instability as discussed in detail in Cairns et al. (2003).

Under the assumptions of 100% polarized orthogonal modes, the polarization distributions of the single pulse data can also be reproduced from a relatively small number of random variables. This extends the work of McKinnon & Stinebring (2000) by implying that the sub-pulses can be characterised in such a way (not just a particular phase bin) and that log-normal statistics are clearly the dominant statistic for the observed flux density distributions. Any physical model must therefore have as its basis an emission beam which is 100% polarized and which is subsequently split through e.g. bi-refringence into two orthogonally polarized modes. However, it is clear that propagation effects, not taken into account in the simulation process, strongly affect the position angle and the degree of circular polarization seen in single pulses. This situation has also been seen previously by a number of observers and discussed in the context of theoretical models by Petrova (2001). In her model, the polarization modes recombine in the ‘polarization limiting region’. This mode coupling then allows for significant deviations in the position angle and in the circular polarization, with greater position angle deviations leading to increased circular polarization. The data here do not entirely support this picture. Comparison between Figs 7 and 8 show that the former has a much narrower spread of position angles than the latter but that the percentage V distributions look very similar. The effects of refraction in the magnetosphere are also likely to be important as the two modes are affected differently. This is likely to result in a mode mixture at a given observer’s longitude and such an effect would also serve to broaden the position angle distribution.

The applicability of this method to other pulsars is not clear. For example, in the Vela pulsar, Krishnamohan & Downs (1983) showed that four Gaussian components fitted the total intensity data and could account for some of the single pulse data. More recent observations with much higher time resolution \[Kramer, Johnston & van Straten 2002\] showed that microstructure was prevalent in the single pulse data. Clearly, no amount of Gaussian fitting will reproduce microstructure. Other pulsars, such as PSR B0950+08 show a bewildering variation in the amplitude and location of single pulse components, and there is almost no resemblance between the single pulses and the integrated profile. Many pulsars, such as PSR B1133+16, show orthogonal mode jumps near the peak of a pulse component and again such a feature is hard to reproduce using simple Gaussian decomposition.

However, it may be that there is a class of pulsars similar to PSR B1641$``$45 for which this method will work very well. One possibility is that the scatter broadening in this pulsar blurs the features with time-scales less than the sub-pulse time. The same was also true of the low resolution observations of the Vela pulsar where the dispersion measure smearing also blurred out the microstructure. Processing techniques which involve smoothing over the microstructure are also possible. For example, van Leeuwen et al. (2002) fitted Gaussian components to PSR B0809+74 in order to study the drift bands and noted how similar the Gaussians were in amplitude and width.

It therefore appears to be possible to get a handle on processes on the scale of tens of milliseconds (the subpulses) without delving into the complexities which occur on timescales significantly shorter. The data seem to favour the idea of simultaneously occuring modes of emission, with propagation through the magnetosphere influencing the degree of circular polarization and the position angle variations.

## Acknowledgments

The Australia Telescope is funded by the Commonwealth of Australia for operation as a National Facility managed by the CSIRO. I thank W. van Straten for help with the data reduction, A. Karastergiou for useful discussions and the referee for improving the content of the paper.


text:
# Generic features of the fluctuation dissipation relation in coarsening systems

## Abstract

The integrated response function in phase-ordering systems with scalar, vector, conserved and non conserved order parameter is studied at various space dimensionalities. Assuming scaling of the aging contribution $`\chi _{ag}(t,t_w)=t_w^{a_\chi }\widehat{\chi }(t/t_w)`$ we obtain, by numerical simulations and analytical arguments, the phenomenological formula describing the dimensionality dependence of $`a_\chi `$ in all cases considered. The primary result is that $`a_\chi `$ vanishes continuously as $`d`$ approaches the lower critical dimensionality $`d_L`$. This implies that i) the existence of a non trivial fluctuation dissipation relation and ii) the failure of the connection between statics and dynamics are generic features of phase ordering at $`d_L`$.

After the groundbreaking work of Cugliandolo and Kurchan Cugliandolo93  on mean field spin glasses, the study of the out of equilibrium linear response function has been gaining an increasingly important role in the understanding of slow relaxation phenomena. The key concept is that of the fluctuation dissipation relation (FDR) Cugliandolo2002 . In terms of the integrated response function $`\chi (t,t_w)`$, i.e. the response to an external field acting in the time interval $`(t_w,t)`$, an FDR arises if $`\chi (t,t_w)`$ depends on time only through the unperturbed autocorrelation function $`C(t,t_w)`$. If this happens, there remains defined a function $`\chi =S(C)`$ which generalizes the fluctuation dissipation theorem into the out of equilibrium regime.

The existence of an FDR is important for several reasons Cugliandolo2002 . Here we focus on a specific aspect: to what extent the FDR shape is revealing of the mechanism of relaxation and of the structure of the equilibrium state. In particular, we aim at dispelling the common belief that relaxation by coarsening and a simple equilibrium state do necessarily imply a flat or trivial FDR, i.e. $`S(C)=1q_{EA}`$ when $`C`$ falls below the Edwards-Anderson order parameter $`q_{EA}`$.

To appreciate the relevance of the problem, consider that, by reversing the argument, the observation of a non flat FDR,would rule out relaxation by coarsening. This is a statement of far reaching consequences. For instance, an argument of this type plays a role in the discrimination between the mean field and the droplet picture of the low temperature phase of finite dimensional spin glasses Ruiz . In that case, the final conclusion may well be right, but for the argument to be sound, the behavior of the response function, when relaxation proceeds by coarsening, needs to be thoroughly understood.

As a contribution in this direction, we have undertaken a large program of systematic investigation of the FDR in the context of phase ordering systems Bray94 , which provide the workbench for the study of all aspects of relaxation driven by coarsening. We have considered pure ferromagnetic systems quenched from above to below the critical point. We have covered the whole spectrum of systems with non conserved (NCOP), conserved (COP), scalar ($`N=1`$) and vector ($`N>1`$) order parameter at different space dimensionalities $`d`$, where $`N`$ is the number of components of the order parameter. The manifold of the systems considered is displayed in Table I. Some of these (marked by a dot) have been studied before. The important novelty is that, with the new entries, the picture becomes rich enough to promote to generic the behavior previously observed in the particular case of the Ising model Lippiello2000 ; Corberi2001 ; Castellano ; Corberi2003  and in the large $`N`$ model Corberi2002 . Namely, that a flat FDR is obtained for $`d>d_L`$, while a non flat FDR is found for $`d=d_L`$, where $`d_L`$ is the lower critical dimensionality. The implication is that a flat FDR is not a necessary condition for coarsening.

To explain in more detail, let us recall Cugliandolo2002  that, quite generally, one can write $`\chi (t,t_w)=\chi _{st}(tt_w)+\chi _{ag}(t,t_w)`$. The first is the stationary contribution due to the fast degrees of freedom which rapidly equilibrate with the bath, while the second is the aging contribution coming from the slow out of equilibrium degrees of freedom. What one can also show, in general, is that a flat FDR is obtained if $`\chi _{ag}(t,t_w)`$ vanishes asymptotically Cugliandolo2002 ; Corberi2001 . Now, in phase ordering for large $`t_w`$ one expects the scaling behavior

$$\chi _{ag}(t,t_w)=t_w^{a_\chi }\widehat{\chi }(t/t_w)$$
(1)
from which it follows that the FDR is or is not flat according to $`a_\chi >0`$ or $`a_\chi 0`$ (here we will restrict to the case with $`a_\chi =0`$). Therefore, investigating the FDR shape requires the investigation of $`a_\chi `$.

Let us see what is the situation with this exponent. In the Ising model Lippiello2000 ; Corberi2001 ; Castellano ; Corberi2003  and in the large $`N`$ model Corberi2002 , we have found that $`a_\chi `$ depends on dimensionality according to

$$a_\chi =\{\begin{array}{cc}\delta \left(\frac{dd_L}{d_Ud_L}\right)\hfill & \text{for}d<d_U\hfill \\ \delta \text{with log corrections}\hfill & \text{for}d=d_U\hfill \\ \delta \hfill & \text{for}d>d_U\hfill \end{array}$$
(2)
where $`\delta `$ enters the time dependence of the density of defects and $`d_L`$, $`d_U`$ are the two special dimensionalities (with $`d_L<d_U`$) where $`a_\chi =0`$ and above which $`a_\chi =\delta `$, respectively. The density of defects goes like $`\rho (t)L(t)^nt^\delta `$, where $`L(t)t^{1/z}`$ is the typical defect distance, $`z`$ is the dynamic exponent and $`n=1`$ or $`n=2`$ for scalar or vector order parameter Bray94 . Hence, $`\delta =n/z`$.

Opposite to Eq. (2) stands a qualitative argument according to which $`\chi _{ag}(t,t_w)`$ ought to be simply proportional to the defect density Barrat98 ; Franz98 . Namely, $`a_\chi =\delta `$ at any dimensionality. There are no measurements or derivations of $`a_\chi `$ supporting such a statement. Furthermore, the argument is incompatible with the exact result $`a_\chi =0`$ Lippiello2000 ; Godreche2000  in the $`d=1`$ Ising model. The implication for the FDR is that it should be always flat, except for possible pathological cases. This point of view is diametrally opposed to the one underlying Eq. (2), where the vanishing of $`a_\chi `$ at $`d=d_L`$, far from being a pathology, comes out as a limiting behavior in the smooth dimensionality dependence of $`a_\chi `$ for $`d<d_U`$.

The understanding of the FDR in coarsening systems requires to clarify whether $`a_\chi `$ does or does not to depend on $`d`$. Here we present strong evidence supporting Eq. (2) as the generic pattern of behavior.

We have computed $`\chi (t,t_w)`$ for systems quenched from infinite to zero final temperature. This is computationally efficient and can be done without loss of generality, since all quenches below $`T_C`$ are controlled by the $`T=0`$ fixed point Bray94 . In all cases we have used the time dependent Ginzburg-Landau equation Bray94 , except for NCOP with $`N>1`$ and $`d>2`$ where the Bray-Humayun Bray90  algorithm has been used nota2 . The stationary response has been computed from equilibrium simulations. The aging part, then, has been obtained from $`\chi _{ag}(t,t_w)=\chi (t,t_w)\chi _{st}(tt_w)`$. To get $`a_\chi `$, one ought to extract the $`t_w`$ dependence of $`\chi _{ag}(t,t_w)`$ for fixed $`x=t/t_w`$ Corberi2003 . However, this is computationally very demanding and would make it impossible to get the vast overview we are aiming at. So, we have measured $`a_\chi `$ from the large $`t`$ behavior for a fixed $`t_w`$, assuming $`\chi _{ag}(t,t_w)t^{a_\chi }`$. This holds if $`\widehat{\chi }(x)x^{a_\chi }`$ for $`x1`$, which has been verified in the NCOP scalar case Corberi2001 ; Corberi2003  and it is an exact result in the soluble models Lippiello2000 ; Corberi2002 . The assumption is that it holds in general. The choice of $`t_w`$ is inessential provided it is larger than some microscopic time necessary for scaling to set in Corberi2003 .

The time dependence of $`\chi _{ag}(t,t_w)`$ is depicted in Figs. 12 and 3. We have extracted $`a_\chi `$ from the asymptotic power law decay and we have collected all results, old and new, in Table II. At $`d_L`$ we have used the parametric plot $`\chi _{ag}(C)`$ (insets of Figs. 12 and 3), showing more effectively the absence of asymptotic decay, due to $`a_\chi =0`$. In Table II, we have also reported the values of $`a_\chi `$ predicted by Eq. (2). The comparison with the computed values is quite good. For convenience, we have collected in Table III the values of all the parameters entering Eq. (2). Finally, Fig. 4 provides the pictorial representation of Table II, and it is the main result in the paper.

Let us now comment the results. From Fig. 4 it is evident that the pattern of behavior predicted by Eq. (2) is obeyed with good accuracy in the scalar cases, with $`d_L=1`$ and $`d_U=3`$. In the vector cases, given the great numerical effort needed, values of $`N`$ were chosen according to the criterion of the best numerical efficiency, together with the requirement to simulate both systems with ($`N<d`$) and without ($`N>d`$) stable topological defects. The overall behavior of the data in Fig. 4 shows that Eq. (2) well represents the dimensionality dependence of $`a_\chi `$ also in the vector case with $`d_L=2`$ and $`d_U=4`$. Finally, the insets in Figs 12 and 3 (together with the analogous figures for the $`d=1`$ Ising model in Ref.s Lippiello2000 ; Castellano  and in the large $`N`$ model Corberi2002 ) show quite clearly that $`a_\chi =0`$ and a non flat FDR are common features in phase ordering kinetics at $`d_L`$.

At this stage Eq. (2) is a phenomenological formula. Apart from the exact solution of the large $`N`$ model Corberi2002 , there is no derivation of Eq. (2). Here we propose an argument for the dependence of $`a_\chi `$ on $`d`$ in the scalar case. It is based on two simple physical ingredients: a) the aging response is given by the density of defects $`\rho (t)`$ times the response of a single defect Corberi2001  $`\chi _{ag}(t,t_w)=\rho (t)\chi _{ag}^s(t,t_w)`$ and b) each defect responds to the perturbation by optimizing its position with respect to the external field in a quasi-equilibrium way. In $`d=1`$ this occurs via a displacement of the defect. In higher dimensions, since defects are spatially extended, the response is produced by a deformation of the defect shape.

We develop the argument for a 2-d system, the extension to arbitrary $`d`$ being straightforward. A defect is a sharp interface separating two domains of opposite magnetization. In order to analyse $`\chi _{ag}^s(t,t_w)`$ we consider configurations with a single defect as depicted in Fig. 5. The corresponding integrated response function reads Corberi2001  $`\chi _{ag}^s(t,t_w)=1/(h^2^{d1})𝑑x𝑑y\overline{S(x,y)h(x,y)}`$, where $`S(x,y)`$ is the order parameter field which saturates to $`\pm 1`$ in the bulk of domains. $`h(x,y)`$ is the external random field with expectations $`\overline{h(x,y)}=0`$, $`\overline{h(x,y)h(x^{},y^{})}=h^2\delta (xx^{})\delta (yy^{})`$ and $``$ is the linear system size. The overbar and angular brackets denote averages over the random field and thermal histories, respectively. With an interface of shape $`z_w(y)`$ at time $`t_w`$ (Fig. 5), we can write $`\chi _{ag}^s(t,t_w)=1/(h^2^{d1})\overline{_{\{z\}}E_hP_h(\{z(y)\},t)}`$, where $`P_h(\{z(y)\},t)`$ is the probability that an interface profile $`\{z(y)\}`$ occurs at time $`t`$ and $`E_h=_0^{}𝑑y_{z_w(y)}^{z(y)}𝑑xh(x,y)\text{sign}[z(y)z_w(y)]`$ is the magnetic energy. We now introduce assumption b) making the ansatz for the correction to the unperturbed probability $`P_0(\{z(y)\},t)`$ in the form of a Boltzmann factor $`P_h(\{z(y)\},t)=P_0(\{z(y)\},t)\mathrm{exp}(\beta E_h)P_0(\{z(y)\},t)[1\beta E_h]`$. Then $`\chi _{ag}^s(t,t_w)=1/(h^2^{d1})\overline{_{\{z\}}E_h(1\beta E_h)P_0(\{z(y)\},t)}`$. Taking into account that the term linear in $`E_h`$ vanishes by symmetry and neglecting $`z_w(y)`$ with respect to $`z(y)`$ for $`tt_w`$, we eventually find $`\beta ^1\chi _{ag}^s(t,t_w)=^{1d}_{\{z\}}_0^{}𝑑y|z(y)|P_0(\{z(y)\},t)`$. This defines a length which scales as the roughness of the interface Henkel2003  given by $`W(t)=[^{1d}_{\{z\}}𝑑yz(y)^2P_0(\{z(y)\},t)]^{1/2}`$. The behavior of $`W(t)`$ in the coarsening process can be inferred from an argument due to Villain Abraham89 . In the case $`d3`$, when interfaces are rough Rough , for NCOP one has $`W(t)t^{(3d)/4}`$, while for COP $`W(t)t^{(3d)/6}`$, with logarithmic corrections in both cases for $`d=3`$. For $`d>3`$ interfaces are flat and $`W(t)const.`$ Finally, multiplying $`\chi _{ag}^s`$ by $`\rho (t)L(t)^1`$ Eq. (2) is recovered note2  and $`d_U`$ is identified with the roughening dimensionality $`d_R=3`$. The relevance of roughening in the large time behavior of $`\chi (t,t_w)`$ has been independently pointed out by Henkel, Paessens and Pleimling in Ref. Henkel2003 . The crucial difference with these authors is that they believe roughening to be unrelated to aging behavior, while we claim the opposite.

In summary, we have investigated the scaling properties of the response function over a large variety of systems designed to bring forward the generic features when relaxation is driven by coarsening. The primary result is that the exponent $`a_\chi `$ depends on dimensionality and that it vanishes smoothly as $`dd_L`$. This implies that a non trivial FDR is not exceptional, rather is the rule for coarsening systems at $`d_L`$. Another important consequence is that the failure of the connection between statics and dynamics at $`d_L`$ Corberi2001  is also a generic feature of coarsening. The connection between the FDR and the overlap probability function is derived Franz98  under the assumptions of stochastic stability and that $`\chi (t,t_w)`$ goes to the equilibrium value as $`t\mathrm{}`$. The latter assumption does not hold at $`d_L`$ due to the existence of a non flat FDR (insets of Figs 1,2,3), which makes the limiting value of $`\chi (t,t_w)`$ to rise above the equilibrium value. Obviously, the important and, as of yet, unanswered question is why all this happens at $`d_L`$. The scaling behavior of the response function reported in this Letter adds to the many already existing challenges posed by a theory of phase ordering kinetics.

This work has been partially supported from MURST through PRIN-2002.


text:
# Opinion formation in networked societies with strong leaders

## I Modeling opinion formation in agent based computer models

Computer simulations are quickly becoming recognized and accepted tool in social sciences, viewed as a way of testing hypotheses and predictions and revealing some simple mechanisms underlying complex behaviour. The studies are especially useful when the studied social phenomena can be quantified, and when enough empirical data exists to compare the numerical models with reality.

One of such areas is opinion formation in societies, measurable both on the large scale (for example through election results or polls) or on small scale direct experiments. Computer simulation models of opinion forming date back to 1960s, however, the particular model we base our research on is the one of Nowak et al. , developed further by Nowak and Lewenstein  and Kacperski and Hołyst , Hołyst et al. . The general approach analyses the formation of public opinion through interactions between individual members of the society, taking into account differences in receptiveness, strength of influence and preexisting biases. The original work Nowak et al.  has shown that interesting macroscopic behaviour can result from simple microscopical properties of computer agents.

The basis for the models in Nowak et al. , Nowak and Lewenstein , Kacperski and Hołyst , Hołyst et al. , which will be generally followed here, is the following:

* A set of $`N`$ interacting agents form a society. Interactions take form of one to one activities. The dynamics of the total system is assumed to take place through discrete time steps, with updates done globally.
* Each agent $`i`$ has, at a given time, his ‘opinion’ $`\sigma _i`$. The changes of aggregates of ‘opinion’ within the society or certain subgroups is the topic of the research.
* Each agent is characterized by the strength of his possible influence on other agents. This allows to model situations of uneven distribution of influence, such as presence of strong leader(s).
* As a special case, one of the agents (the leader) is assumed to have the strength of influence much greater than the rest of the agents.
* Interaction between particular agents is governed by the strength of agents as well as the remoteness of the interacting agents. We use the term ‘remoteness’ in abstract sense, to denote the measure of the social separation between members of society and its effect on transmission of opinions. The overall social impact of agent $`i`$ on agent $`j`$ is given by combination of acting agent strength and the separation between $`i`$ and $`j`$.
* The agents interact and influence each other in turns, changing their opinion after each full turn of interactions take place. Agents are allowed to interact with themselves, this mimics the phenomenon of self-support, or inclination to hold agent’s present opinion.
* The model allows for extra-social influence or bias, applying it uniformly to all agents.
* The models may allow for the noise in communication and changing individual opinion by adding an equivalent of temperature to the simulations.

It is extremely interesting that the results of the works cited above have significant dependence on the way the ‘remoteness’ between the agents is introduced. The works of Nowak et al. Nowak et al. , Nowak and Lewenstein  have used simple spatial two-dimensional (2D) model of agent distribution. It has been argued that

> ‘People are more likely to interact with neighbors, that is those who live close to them in physical space’ …‘Our choice of a 2D lattice represents quite well the physical distribution of people on flat surfaces. The results of studies conducted in Boca Raton, Warsaw and Shanghai, have shown that the probability of social interactions is decreasing as a square of physical distance’ (Nowak and Lewenstein ).

Results obtained for the 2D lattice based society have shown that there are localized ‘bubbles’ of uniform opinion around strong leaders, growing and merging or diminishing with changes of external influence. The geometrical model on a 2D disc has been used by Kacperski and Hołyst , with the measure of decrease of the influence with distance studied mostly using linear relation.

In their later works of Kacperski, Hołyst and Schweitzer Kacperski and Hołyst , Hołyst et al.  have argued that the social distance, which does not have to fulfill the same conditions as geometrical one (such as, e.g. triangle inequality) should be modelled by more general and flexible model. The authors introduce the notion of social immediacy $`m_{ij}`$ (between agents $`i`$ and $`j`$ ). The immediacy $`m_{ij}`$ does not need to be equal $`m_{ji}`$ and the distribution of values of $`m_{ij}`$ may be arbitrary. In Kacperski and Hołyst  three immediacy distributions were analysed: uniform random distribution $`0m_{ij}2\overline{m}`$, exponential distribution and discrete multimodal distribution. The results show rapid jumps in majority-minority proportion of opinion and presence of local clusters of opinion, which in abstract social distance space correspond to Nowak’s bubbles of opinion.

## II Opinion formation in networked societies

### II.1 Network systems – an overview

Recent studies of various networked systems (see for example Strogatz , Albert and Barabási , Dorogovtsev and Mendes \[2002a\], Dorogovtsev et al. , Dorogovtsev and Mendes \[2002b\], Newman , Newman and Park ), including biological and social systems, the Internet and World Wide Web have shown several universal and interesting effects. The most important among them are ‘small world effect’, degree distributions exhibiting power-law or truncated power-law behaviour and significant clustering effects. Large number of results indicate that the formation and properties of networked systems found in nature — including social networks — shows certain regularities and rules. In this paper we study the same basic model of opinion formation as Nowak, Hołyst, Kacperski and co-workers, but we place the agents on nodes of several specific networks, differing in topology of connections. The different network topologies result in significant differences in network characteristics, which might influence the process of opinion formation. The main networks used in this work are described below.

In random (RAND) networks the links are distributed randomly, thus we have a meshed network of links, with agents differing in number of the neighbours, and no general structure. On the average, the maximum distance between agents in random network grows as $`\mathrm{log}N`$ and the number of nieghbours per agent is given by the Poisson distribution.

The nearest neighbor (NN) nets are formed by linking together a fixed number of closest neighbors. The easiest way to visualize such network is to place the agents on an imaginary circle and connecting each agent to $`n`$ neighbors. Interesting property of NN networks is that for small filling factors, agents on the opposite points at the circle to communicate must go through many intermediaries. For 4000 agents and number of neighbors set at 10, the longest ‘distance’ is 200 ‘hops’. Due to such large separation any change of behavior of the agent $`i`$ is seen immediately only by his closest neighbors but the majority of the society is affected only after filtering by intermediate agents located far from $`i`$.

The small world Strogatz-Watts (SW) networks, introduced and popularized in recent years, reproduce a curious fact observed in many natural and human-produced networks, namely that the distance between any two nodes of the network, measured in number of links needed to connect them, is usually much smaller than that in nearest neighbor or even random networks of the same filling ratio Milgram . The name of the network category comes exactly from such observation. One of the ways to realize a SW network is through a simple reworking of the NN model: one cuts a small (even very small!) fraction of the links from between nearest neighbors and applies them instead between random agents. Keeping in mind the visualization of NN networks as nodes along a circle, this corresponds to adding connections that criss-cross the circle at random. Due to such shortcuts, even if their number is very small, the average distance between any two nodes drops dramatically. Thus we have a network that for each agent, locally is very similar to NN model (as most of the neighbors are, in fact, the same), but globally the communication through the network is much faster.

The scale-free Albert-Barabási network (AB) reproduce another effect found in natural and artificial systems. In AB networks the dispersion of the number of immediate neighbors (degree distribution) scales according to power law rather than exponential or Gaussian. This type of networks results from interplay of two processes: growth of the network and preferential attachment (the ‘rich get richer’ principle). As a result, the AB networks exhibit much greater presence of highly connected agents than other types of networks considered here.

We study also two types of hierarchical networks, simple one (HS), with each agent having one link up (to its ‘parent’ on higher hierarchy level) and fixed number of links down (to descendants), and clustered hierarchical network (HCL), in which in addition to links up and down, the agents having the same parent are all linked together (within the same level). While the simple hierarchy shows no clustering at all, the clustered version preserves the division of the society into separate groups and levels, while providing high clustering ratio within groups

For completeness, in addition to the abstract space networks, we study also two traditional spatial geometries, 2D square array (2DSQ) and 3D cubic array (3DCU). In both cases agents may be linked with all immediate neighbors, including diagonal ones (thus the 2D agent can have 8 neighbours, and in 3D case up to 26 neighbours). We have used periodic boundary conditions.

In our simulations we have tried to compare results for different networks but the same average number of neighbors $`N_N`$. For the regular 2D and 3D networks, where the geometrical number of neighbours is fixed, to achieve the value of $`N_N`$ correspondig to other networks, we have randomly removed appropriate number of links. The only exceptions were the hierarchical structures, where to preserve the essential top-down asymmetry of the network, we have had to accept the lower number of connections at the lowest level (edge of the network due to finite size), where only links up (HS) or links up and in horizontal cluster (HCL) are present. The hierarchical networks have therefore smaller average number of links, but preserve the $`N_N`$ within the core of the network.

Figures 1 and 2 compare the two main properties of the discussed networks relevant to the present research. The results shown were calculated for the average number of neighbours $`N_N=6`$. The first figure compares the distribution of immediacies calculated according to Equation 1 ($`m_{ij}=1/d_{ij}`$, i.e. with $`\alpha =1`$) for various networks. In all cases one sees very pronounced maxima in distribution of $`m_{ij}`$ (notice the lin-log scale), but the most probable values differ: they are the highest for the RAND and AB networks, smaller for SW, 3DCU and HS networks, and very small for NN, 2DSQ and HCL. Increasing $`\alpha `$ would shift the distributions of $`m_{ij}`$ to even smaller values. It is worth noting that in the NN system case due to very large distances most $`m_{ij}`$ values are very small indeed.

Examining the distributions of the number of neighbours one sees that most networks show some limited dispersion around the average $`N_N=6`$, the only exception being the scale-free network of Albert-Barabási. Out of the ‘normal’ networks, the dispersion is largest for RAND system, while for NN network the number of neighbors is fixed $`N_N6`$. The AB network contains quite a number of agents with high connectivity, even up to 250 in this particular case. This turns up as the origin of one of the interesting results of our simulations.

### II.2 Use of networks in modelling opinion formation

Even if we concentrate only on studies of strictly social phenomena, there is growing body of data that the communication processes in human societies are best modelled by complicated network arrangements (for examples see Albert and Barabási , Newman and Park ). The ‘small world’ theorem, well established experimentally, stating that the separation between randomly chosen members of society is much smaller than expected from geometrical and NN models, and even from randomly linked model, can in obvious way influence the process of communication between members of society, and thus the formation of opinions, or the ‘distance’ from leaders to ‘normal’ people. On the other hand, the presence of extremely highly connected nodes (seen in Albert-Barabási models) gives particular influence to some individual members of the society.

The aim of the current work is to present some results of the same type of simulations used before on geometrical or abstract connections topologies applied here on various types of networks. The differences between the general network properties might shed some light on the way the opinion shift occurs, such as the Nowak’s ‘bubbles of opinion’ in 2D model.

## III Details of the model

As mentioned above, our model followed almost exactly that of Hołyst and Kacperski. There are $`N`$ agents in the simulated system (we have used $`N=4096`$). Each agent $`i`$ is described by strength $`s_i`$ (unchanged during the simulation run) and the value of opinion, described by $`\sigma _i=\pm 1`$. The process of establishing the opinion in the society is modelled in discrete time steps. At each step the opinion of any given agent is updated in accordance with the combined influence of the other agents, agent’s own self-influence and possible external conditions. The interaction between agents, for example agent $`j`$ influence on agent $`i`$ is given by combination of agent $`j`$ strength $`s_j`$ ‘remoteness’ of the agents, described by immediacy $`m_{ij}`$. In our work we set a very natural model for $`m_{ij}`$ in networked system

$$m_{ij}=\frac{1}{d_{ij}^{}{}_{}{}^{\alpha }}$$
(1)
where $`d_{ij}`$ is the distance between agents $`i`$ and $`j`$ measured in number of network links and the exponent $`\alpha `$ determines the ratio of decrease of immediacies with increasing distance (we have simulated societies setting $`\alpha =1,2,3`$. As $`d_{ij}=1,2,3\mathrm{}`$ the immediacies are always $`1`$. The self-influence terms $`m_{jj}`$ were set at 2 for all agents, which corresponds to relatively high tendency to preserve one’s own opinion.

The strengths $`s_i`$ were chosen randomly from interval $`[0,2\overline{s}]`$. The only exception was a single agent, called the leader, for whom the strength $`s_L`$ was much greater than the average $`\overline{s}`$. Moreover, we have fixed $`\sigma _L=1`$ and set $`m_{LL}1`$ so that the leader’s opinion $`\sigma _L`$ would remain unchanged during the simulation runs.

The states of the agents are updated synchronously in discrete time steps according to the following formula:

$$\sigma _i(t+\mathrm{\Delta }t)=\{\begin{array}{cc}1\hfill & \text{with probability }\frac{\mathrm{exp}(I_i/T)}{\mathrm{exp}(I_i/T)+\mathrm{exp}(I_i/T)}\hfill \\ & \\ 1\hfill & \text{with probability }\frac{\mathrm{exp}(I_i/T)}{\mathrm{exp}(I_i/T)+\mathrm{exp}(I_i/T)}\hfill \end{array}$$
(2)
where $`T`$ is the measure of randomness (‘social temperature’) and

$$I_i=\underset{j=1}{\overset{N}{}}s_jm_{ij}\sigma _j+h.$$
(3)
The value of $`h`$ measures uniform external influence on the system.

It is useful to introduce here the rescaled values. The scaling is given by

$`s_L^R`$ $`=`$ $`{\displaystyle \frac{s_L}{N\overline{s}}},`$ (4)
$`h^R`$ $`=`$ $`{\displaystyle \frac{h}{N\overline{m}\overline{s}}},`$ (5)
$`T^R`$ $`=`$ $`{\displaystyle \frac{T}{N\overline{m}\overline{s}}},`$ (6)
with averages $`\overline{m},\overline{s}`$ excluding the leader. Additionally, we define the reference background field $`B=N\overline{s}\overline{m}`$, which is the maximum value of the background influence of all non-leader agents if they all have $`\sigma _j=1`$. Now the social impact on agent $`i`$, in terms of the rescaled values is given by

$$I_i=s_L^R\frac{m_{iL}}{\overline{m}}B+\underset{jL}{}s_j^R\frac{m_{ij}}{\overline{m}}\sigma _jB+h^RB.$$
(7)
The fact, that through the Equation 1 the network itself establishes values of the immediacies $`m_{ij}`$ allows to connect topological properties of the networks (which in turn may be closely related to conceptual characteristics of the society) to the opinion formation process. Thus, our work is much closer to possible experimental and practical uses than the previous research.

## IV Results

In our simulations we were particularly interested in searching for conditions under which the leader can convince a significant part or the whole society to his views. We have studied the fraction $`f`$ of leader’s supporters as function of unfavorable external influence ($`h^R<0`$) and leader’s strength (rescaled) $`s_L^R`$. We have compared results of simulations performed with three types of starting conditions: assuming that initially all agents have the same opinion as the leader ($`\sigma _i1`$), random initial distribution of $`\sigma _i`$ and assuming that all agents disagree initially with the leader ($`\sigma _i1`$). Results of the three types of initial conditions are presented as rows in figures.

In some of the simulations the leader was chosen randomly, with the number of immediate neighbors close to average, in other cases the highest connected (HC) agent was chosen to be the leader. Results are presented in the left and right columns of figures 37 (the only exception being the NN network where all agents have exactly the same connectivity, so only one set of results is presented). The number of connections of random and HC leader are rather similar in most cases — with the exception of the scale free AB network, where there may be a difference of 1–2 orders of magnitude.

Figures 37 present contour plot of fraction of leader supporters $`f`$ for some of the networks studied. The contours correspond to the following values of the fraction of leader’s supporters: 0.01, 0.10, 0.25, 0.50, 0.75, 0.90, 0.99 (from red to green). The social temperature $`T^R`$ was set at 0.3. The results, are in most cases remarkably similar. The exceptions, which will be discussed below are the NN network and the AB system.

Figure 8 presents $`f`$ as function of $`s_L^R`$ for a fixed value of $`h^R=1.5`$, for a set of network types and initial conditions. Left column corresponds to extremely disfavourable initial conditions $`\sigma _j1`$, right column to a random distribution with average $`\overline{\sigma }=0`$. Two sets of lines correspond to leaders placed at random (red) and in a HC position (green).

Let’s discuss first a typical situation, such as the one found in RAND and SW networks, as well as the AB network when the leader is chosen randomly. For low $`s_L^R`$ and significant negative $`h^R`$ there are no supporters of leader’s opinion. Increasing $`s_L^R`$ leads at first to establishment of small cluster of supporters (manifested in the area where $`0.01<f<0.25`$. In figure 8 this is clearly visible as plateaus of $`f`$. Further increase leads to rapid transition of the whole system to polarized, supporting state of $`f1`$. The transition region is approximately given by equation

$$s_L^{}=h^R+C,$$
(8)
where the value of the constant $`C`$ depends on the temperature and initial conditions. This equation simply reflects the strength of the leader needed to overcome the average background field of initial configuration. In low $`T^R`$ limit $`C1`$ for the $`\sigma _i=1`$ starting condition and $`C=1`$ for the opposite case; for random starting condition $`C=0`$. These results are quite straightforward: when the leader’s strength overcomes the combined effect of the background influence and external conditions, the agents start to change opinion. With the increase of their number the conditions become even more favourable and transition to supportive state occurs.

That the range of values for which the local cluster of supporters is present is greater for the favourable starting conditions than in the random or disfavourable case is also quite evident: due to relatively large value of self-influence terms $`m_{ii}`$ in Equation 7, it is more difficult for the leader to convert an agent than to maintain the supporting agent’s attitude.

The results are presented for moderate values of $`h^R`$ and $`s_L^R`$, comparable to 1. For much larger values the external influence and the leader’s influence dominate any local statistical fluctuations due to quenched disorder. The dynamics of the system is then trivial, given, for any agent $`i`$ by competition between the leader term $`s_L^R\frac{m_{iL}}{\overline{m}}`$ and $`h^R`$. Due to discrete nature of distances $`d_{iL}`$ and therefore $`m_{iL}`$, agents with the same $`m_{iL}`$ are convinced (for given $`h^R`$) at discrete values of $`s_L^R`$, and $`f`$ shows simple multi-step behaviour.

We return now to the special cases. We’ll start with the NN network. Due to the presence of very large separations between agents in NN system, leader’s influence on distant agents is extremely weak, and the only practical way to supportive state is the gradual ‘conviction’ of the growing range of neighbours. As a result we see that the region of $`f`$ grows linearly with increasing $`s_L^R`$. The growth of $`f`$ is slow, and very large values of $`s_L^R`$ are needed to ensure the wholly supportive state. The NN society may be called strongly localized, and the presence of very high separation between agents clearly influences the opinion formation. It should be noted that our NN network is, in fact, a 1-dimensional (1D) network. The 2D and 3D networks, for which the average separation grows not linearly with $`N`$ but rather as $`\sqrt{N}`$ and $`\sqrt[3]{N}`$ respectively should weaken the effect, but still preserve it. Our simulation size of 4096 agents is, however, too small to make it visible, as the furthest distance on 2D network is 32 and for 3D case only 8. On the other hand, for random networks, where the largest distance grows as $`\mathrm{log}N`$ and for SW and AB networks where it is even smaller, the conditions for rapid transition are present and the such transition is observed.

The second interesting case is the AB network with the leader in high connectivity position. Such leader has a lot of neighbours with $`d_{ij}=1`$ and thus with $`m_{ij}=1`$. The influence on those agents is relatively strong, and as these agents form a significant fraction of the society’s size (in our case $`250/40966`$%) it is much easier to achieve the supportive state

$$s_L^{}K(h^R+C),$$
(9)
with $`K<1`$.

Figure 7 presents results for the AB network with $`\alpha =2`$. The decrease of $`m_{ij}`$ with distance being much faster, there are significant differences to the $`\alpha =1`$ case. The transition region is here more gradual, with visible spread of the contours especially for large negative $`h^R`$ region, also for the leader in random position. While the contour $`f=0.90`$ falls into similar position as in the $`\alpha =1`$ case, the $`f=0.99`$ contour is shifted to much higher values of $`s_L^R`$ in the case of disfavourable and random initial conditions. This is due to the difficulty of converting even moderately remote agents due to rapid decrease of $`m_{ij}`$.

The results presented in this work indicate that simple models of networked societies, based on network types actually found in nature and created by humans, together with a simple formula for turning social distance to influence strength one is able to derive the picture of the conditions necessary to bring the consensus of opinion due to the influence of a single strong individual. For most of the networks the simulated behaviour is quite simple, but for some of the network topologies (especially for the scale-free Albert-Barabási networks) the results are significantly different. The conditions necessary to convert the population to leader’s opinion depend crucially on leader’s position and connectivity in the network. Also, for strongly localized nearest neighbour networks the achievement of fully supportive state requires far greater leader strengths than for small-world or even random networks.

As the Strogatz-Watts and Albert-Barabási networks are found to be relevant in many human societies and activities (such as friendship and acquaintance networks, sexual contacts, scientific collaborations, scientific citations, internet WEB pages links and even telephone calls structure) the presented model may be found of importance in the analysis of the processes of opinion formation in some of these cases.


text:
# Power and Heat Fluctuation Theorems for Electric Circuits

(March 31, 2004)

## Abstract

Using recent fluctuation theorems from nonequilibrium statistical mechanics, we extend the theory for voltage fluctuations in electric circuits to power and heat fluctuations. They could be of particular relevance for the functioning of small circuits. This is done for a parallel resistor and capacitor with a constant current source for which we use the analogy with a Brownian particle dragged through a fluid by a moving harmonic potential, where circuit-specific analogues are needed on top of the Brownian-Nyquist analogy. The results may also hold for other circuits as another example shows.

Nanotechnology is quickly getting within reach, but the physics at these scales could be different from that at the macroscopic scale. In particular, large fluctuations will occur, with mostly unknown consequences. In this letter, we will investigate properties of electric circuits concerning the fluctuations of power and heat within the context of the so-called *Fluctuation Theorems* (FTs). These theorems were originally found in the context of non-equilibrium dynamical systems theory. Surprisingly, these can be applied also to electric circuits, as we will show, and thus give further insight into their behavior.

Let us first give a brief introduction to the FTs. First found in dynamical systemsEvansetal93 ; GallavottiCohen95a  and later extended to stochastic systemsKurchan98 , these conventional FTs give a relation between the probabilities to observe a positive value of the (time averaged) “entropy production rate” and a negative one. This relation is of the form $`P(\sigma )/P(\sigma )=\mathrm{exp}[\sigma \tau ]`$, where $`\sigma `$ and $`\sigma `$ are equal but opposite values for the entropy production rate, $`P(\sigma )`$ and $`P(\sigma )`$ give their probabilities and $`\tau `$ is the length of the interval over which $`\sigma `$ is measured. In these systems, the above mentioned FT is derived for a mathematical quantity $`\sigma `$, which has a form similar to that of the entropy production rate in Irreversible Thermodynamics.

Apart from an early experiment in a turbulent flowotherexperiment , for quite some time, the investigations of the FTs were restricted to theoretical approaches and simulations. In 2002, Wang et al. performed an experiment on a micron-sized Brownian particle dragged through water by a moving optical tweezer. In this experiment, a *Transient Fluctuation Theorem* (TFT) was demonstrated for fluctuations of the total external work done on the system in the transient state of the system, i.e., considering a time interval of duration $`\tau `$ which starts immediately after the tweezer has been set in motionWangetal02 . In contrast, a *Stationary State Fluctuation Theorem* (SSFT), which was not measured, would concern fluctuations in the stationary state, i.e., in intervals of duration $`\tau `$ starting at a time long after the tweezer has been set in motion. While the work fluctuations satisfy the conventional TFT and SSFTMazonkaJarzynski99 ; VanZonCohen02b , the heat fluctuations satisfy different, extended FTs due to the interplay of the stochastic motion of the fluid with the deterministic harmonic potential induced by the optical tweezerVanZonCohen03a ; VanZonCohen03b . Given the possible problems with identifying the entropy productionWangetal02 ; MazonkaJarzynski99 , in this letter we prefer to consider the work and the heat.

We remark that the *conventional* FTs hold in time-reversible, chaotic dynamical systemsEvansetal93 ; GallavottiCohen95a  and in finite stochastic systems if transitions can occur forward and backwardKurchan98 . However, the general condition for an *extended* FT is unknown.

In view of the well-known analogy of Brownian motion (as in the experiment of Wang et al.) and Nyquist noise in electric circuitsVanKampenMazo , one could ask whether the TFT and SSFT based on the Langevin equation also apply to electric circuits. Electric circuits are interesting as they are directly relevant to nanotechnology and because they lend themselves easily to experiments. Indeed, it turns out the conventional FTs hold for *work* and the extended FTs for *heat*, as we will show here. We emphasize that to connect with previous papers on which this one is based, we will also use the term *work* for the electric circuits, which is nothing but the time integral of the *power*.

To exploit the analogy of the Langevin descriptions for electric circuit and that of the Brownian particle, we first recall the form of the Langevin equation for the Brownian particle in the experiment of Wang et al.VanZonCohen02b :

$`m{\displaystyle \frac{d^2x_t}{dt^2}}=\alpha {\displaystyle \frac{dx_t}{dt}}+\xi _tk(x_tv^{}t).`$ (1)
Here, $`m`$ is the mass of the Brownian particle, $`x_t`$ is its position at time $`t`$, $`\alpha `$ is the (Stokes’) friction coefficient, $`k`$ is the strength of the harmonic potential induced by the optical tweezer, and $`v^{}`$ is the constant speed at which the tweezer is moved. $`\xi _t`$ represents a Gaussian white noise which satisfies

$$\xi _t=0;\xi _t\xi _t^{}=2k_BT\alpha \delta (tt^{}),$$
(2)
where $``$ denotes an average over an ensemble of similar systemsVanKampenMazo , $`T`$ is the temperature of the water and $`k_B`$ is Boltzmann’s constant. Usually, the velocity relaxes quickly, so one can set $`m=0`$ in Eq. (1).

Next we consider in Fig. 1 an electric circuit in which a resistor with resistance $`R`$ and a capacitor with capacitance $`C`$ are arranged *in parallel* and are subject to a constant, nonfluctuating current source $`I`$. Energy is being dissipated in the resistor and according to the fluctuation-dissipation theorem this means there are fluctuations too. To a good approximation, one can use a Gaussian random noise term $`\delta V_t`$ to describe these fluctuationsVanKampenMazo , which is depicted in the figure by a voltage generator. In addition, we define $`q_t`$ as the charge that has gone through the resistor, $`i_t`$ as the current that is going through it (so $`i_t=dq_t/dt`$) and $`q_t^{}`$ as the charge on the capacitor, all at time $`t`$. Using that $`q_t^{}=_0^t𝑑t^{}(Ii_t^{})=Itq_t`$ and that $`V_{ab}=q_t^{}/C`$, standard calculations for electric circuits give

$`0=R{\displaystyle \frac{dq_t}{dt}}\delta V_t{\displaystyle \frac{q_tIt}{C}},`$ (3)
in which $`\delta V_t`$ satisfiesVanKampenMazo

$$\delta V_t=0;\delta V_t\delta V_t^{}=2k_BTR\delta (tt^{}).$$
(4)
We see that Eqs. (3) and (4) are of a very similar form as those for a one-dimensional massless Brownian particle dragged through a fluid by means of a harmonic potential in Eqs. (1) and (2)footnote11 . Table 1 gives all the analogues, i.e., both the well-known Brownian motion-Nyquist noise ones ($`\xi _t`$, $`\alpha `$, $`T`$ vs. $`\delta V_t`$, $`R`$, $`T`$) as well as additional circuit-specific ones.

Given this analogy, we turn to the *heat* fluctuations in this circuit. This heat is developed in the resistor. Thus, the dissipated heat over a time $`\tau `$ is given by the time integral of the voltage over the resistor, $`V_{ab}`$, times the current through it, $`i_t`$, i.e.,

$`Q_\tau ={\displaystyle _0^\tau }𝑑ti_t[i_tR+\delta V_t]={\displaystyle _0^\tau }𝑑ti_t{\displaystyle \frac{q_tIt}{C}},`$ (5)
where we used Eq. (3). This is precisely the quantity found in Refs. VanZonCohen03a ; VanZonCohen03b  for the Brownian particle, when we use the analogies in Table 1. Hence $`Q_\tau `$ in this parallel RC circuit behaves completely analogous to the heat for the Brownian particle and thus we know that it satisfies the extended FT. That is, defining a *fluctuation function* by

$$f_\tau ^Q\frac{k_BT}{Q_\tau }\mathrm{ln}\left[\frac{P(+Q_\tau )}{P(Q_\tau )}\right]$$
(6)
and a scaled heat fluctuation by $`p_QQ_\tau /Q_\tau `$, one has for large $`\tau `$

$$f_\tau ^Q(p_Q)=\{\begin{array}{cc}p_Q+O(\tau ^1)\hfill & \text{if }0<p_Q<1\hfill \\ p_Q\frac{1}{4}(p_Q1)^2+O(\tau ^1)\hfill & \text{if }1<p_Q<3\hfill \\ 2+O\left(\sqrt{(p_Q3)/\tau }\right)\hfill & \text{if }p_Q>3.\hfill \end{array}$$
(7)
Here, we only gave the orders of magnitude of the finite-$`\tau `$ correction terms. Their detailed forms — which differ in the transient and the stationary state — require an involved calculation using the saddle-point method which can be found in Ref. VanZonCohen03b . Note that these calculations need not be redone for the current case but that we can make the substitutions in Table 1 and consider the one-dimensional case as can be obtained from footnote 24 of Ref. VanZonCohen03b .

Next, we discuss the *work* fluctuations in this circuit. The total work done in the circuit is the time integral of the power. The power is the current through the circuit, $`I`$, times the voltage over the whole circuit, $`V_{cd}`$, so that

$`W_\tau ={\displaystyle _0^\tau }𝑑tI\left[{\displaystyle \frac{q_tIt}{C}}\right].`$ (8)
This happens to be precisely the form of the work as we found in the Brownian caseVanZonCohen02b , if we use Table 1. This is somewhat surprising because those analogues were based on Eq. (3) which in principle only involves the current through the resistor, while $`W_\tau `$ is the work done on the whole circuit. Since the work fluctuations for the Brownian particle satisfy the conventional FT for $`\tau \mathrm{}`$MazonkaJarzynski99 ; VanZonCohen02b , by analogy we know the same will hold here, i.e.,

$`f_\tau ^W{\displaystyle \frac{k_BT}{W_\tau }}\mathrm{ln}\left[{\displaystyle \frac{P(+W_\tau )}{P(W_\tau )}}\right]={\displaystyle \frac{p_W}{1\epsilon (\tau )}},`$ (9)
where $`p_W=W_\tau /W_\tau `$. Here, for the TFT, $`\epsilon =0`$ while for the SSFT (for details see Ref. VanZonCohen02b )

$`\epsilon (\tau )={\displaystyle \frac{\tau _r(1e^{\tau /\tau _r})}{\tau }},`$ (10)
where $`\tau _r=RC`$. Note that $`\epsilon (\tau )0`$ for $`\tau 0`$.

We note that the behavior of the heat fluctuations in Eq. (7) differs from that of the work fluctuations in Eq. (9) due to exponential tails of the distribution of the heat fluctuationsVanZonCohen03a ; VanZonCohen03b , while those of $`P(W_\tau )`$ are Gaussian.

A second example of a circuit that satisfies the extended heat fluctuation theorem is depicted in Fig. 2. Here, the resistor with resistance $`R`$ and the capacitor with capacitance $`C`$ are arranged *in series* and are subject to a linearly increasing voltage source $`V(t)=At`$. Again, there is a thermal noise generator next to the resistance. The definitions of $`i_t`$ and $`q_t`$ are still that they are the current and charge through $`R`$ at time $`t`$, respectively, but $`q_t^{}=q_t`$ here, since all charge that runs through $`R`$ ends up in $`C`$. We find the appropriate Langevin equation as follows. The imposed voltage, $`V(t)`$, is equal to the potential difference $`V_{ad}=i_tR+\delta V_t+q_t/C`$ (cf. Fig. 2). Using $`i_t=dq_t/dt`$ and $`V(t)=At`$, we find

$`0=R{\displaystyle \frac{dq_t}{dt}}\delta V_t{\displaystyle \frac{q_tCAt}{C}}.`$ (11)
This equation is the same as that of Eq. (3) except that $`I`$ is replaced by $`CA`$. This means that as far as the current through the resistance is concerned, a constant current source with a capacitor in parallel to it is equivalent to a linearly increasing voltage with a capacitor in series with it. It also means that this case is analogous to the Brownian particle dragged through a fluid by a harmonic potential as well. Quantities for this circuit and the Brownian model can therefore be translated into each other using table 1, except for the work, as we will explain below.

Consider now the *heat* developed in the resistor during a time $`\tau `$, which is given by the time integral of the voltage over the resistor, $`V_{ab}`$, times the current through it, $`i_t`$, so

$`Q_\tau ={\displaystyle _0^\tau }𝑑ti_t[i_tR+\delta V]={\displaystyle _0^\tau }𝑑ti_t{\displaystyle \frac{q_tCAt}{C}},`$ (12)
where we used Eq. (11). Note that using $`ICA`$, this form for $`Q_\tau `$ is the same as in Eq. (5) for the parallel case. Hence, the heat in the serial case behaves precisely the same as in the parallel circuit, because we have seen that $`q_t`$ and $`i_t`$ behave the same in both circuits. Thus, we known that the heat in the serial RC circuit, in the parallel RC circuit and in the Brownian system all behave analogously, and all satisfy the extended FTVanZonCohen03a ; VanZonCohen03b  in Eq. (7).

We will now consider the work fluctuations in the serial RC circuit. The work is the time integral of the total current, $`i_t`$, times the total voltage, $`V(t)=At`$, hence

$`W_\tau ^{}={\displaystyle _0^\tau }𝑑ti_tAt.`$ (13)
Even when we use $`ICA`$, this is not the same form as $`W_\tau `$ in Eq. (8) (which is why we added a superscript ), and likewise, using table 1, it is not of the same form as in the Brownian case. Clearly, we cannot use the same results for $`W_\tau ^{}`$ as we obtained for $`W_\tau `$ from the Brownian case or the parallel circuit: we in fact need a additional calculation. For this, we use the method in Ref. VanZonCohen02b . For a Gaussian $`P(W_\tau ^{})`$ one has

$`f_\tau ^W{\displaystyle \frac{k_BT}{W_\tau ^{}}}\mathrm{ln}\left[{\displaystyle \frac{P(+W_\tau ^{})}{P(W_\tau ^{})}}\right]={\displaystyle \frac{p_W^{}}{1\epsilon (\tau )}},`$ (14)
where $`p_W^{}=W_\tau ^{}/W_\tau ^{}`$ and $`\epsilon (\tau )=1V/(2M)`$ with $`M=W_\tau `$ and $`V=(W_\tau W_\tau )^2`$. The quantities $`M`$ and $`V`$ in turn are calculated using the definition of $`W_\tau ^{}`$ in Eq. (13) and the relations $`q_t=CA[t\tau _r(1e^{t/\tau _r})]`$ and $`(q_tq_t)(q_t^{}q_t^{})=k_BTCe^{|tt^{}|/\tau _r}`$ (from Ref. VanZonCohen02b , using Table 1). This yields

$`\epsilon (\tau )=2{\displaystyle \frac{\tau _r^2(1e^{\tau /\tau _r})\tau \tau _re^{\tau /\tau _r}}{\tau ^2}}.`$ (15)
Note that this goes to zero asymptotically as $`1/\tau ^2`$, which is faster than the $`1/\tau `$ decay of $`\epsilon (\tau )`$ in Eq. (10).

These result so far made no explicit reference to nanoscale circuits. Indeed, these results are valid for a circuit of any size. We will now turn to their relevance for nanocircuits. The extended FT shows that large heat fluctuations are more likely to occur than according to the conventional FT, due to the exponential tails of the distribution of heat fluctuations. This could be important for the design of nanostructures because of the rise of temperature due to heat development. Either the *average* current through $`R`$ could be too large or a *large heat fluctuation* might occur. Fluctuations of the energy exist already in an equilibrium system in contact with a heat reservoir, i.e. a circuit for which $`I=0`$. In that case, the energy fluctuations of an $`N`$ atomic resistor will be of the order of $`\sqrt{N}k_BT`$. To explore the properties of nanostructures, composed of just a few atoms, we start with $`N=1`$. For that case, energy fluctuations are of the order of $`k_BT`$ which, if they were not removed, could amount to a significant increase or decrease of the resistor’s temperature. An estimate of this temperature change can be made using the law of Dulong and Petit that (at room temperatures) the specific heat per atom of a solid is $`3k_B`$. Thus, the change in temperature could be as large as $`\mathrm{\Delta }T=k_BT/(3k_B)=T/3=100\mathrm{K}`$, for $`T=300\mathrm{K}`$.

In our Langevin theory for the circuits, we have assumed a constant temperature. To still be able to use the theory, we need that the temperature does not vary too much. A similar condition is needed if the resistor is not to fail. To satisfy this condition, the heat developed will have be transported away at a fast enough time scale $`\tau _T`$. If the heat $`Q_{\tau _T}`$ developed in this time is large enough to significantly increase the material’s temperature, failure may occur. The induced temperature difference, given by $`\mathrm{\Delta }T=Q_{\tau _T}/C_V`$, is insignificant if $`\mathrm{\Delta }T/T1`$, i.e., if

$$\frac{Q_{\tau _T}}{C_VT}1.$$
(16)
This condition should hold both for the average and the fluctuations of $`Q_{\tau _T}`$. Let us first consider the average (i.e., what one would do for macroscopic circuits). In the stationary state of the parallel circuit, the average heat is $`Q_{\tau _T}=IV\tau _T`$, where $`V=IR`$, so

$$\frac{Q_{\tau _T}}{C_VT}=\frac{I^2R\tau _T}{C_VT}=\frac{\tau _T}{\tau _H},$$
(17)
where $`\tau _H=C_VT/I^2R`$. (This in fact also holds for the serial circuit if we replace $`I`$ by $`CA`$.) Equation (16) now becomes

$$\tau _T\tau _H.$$
(18)
That is, the time scale of temperature relaxation given by $`\tau _T`$ has to be fast compared to the heating time $`\tau _H`$. However, even if the average of the heat is well-behaved, the heat fluctuations might still damage the circuit. We need the typical size of these fluctuations. From our previous work on the FTs, we found that for long times $`\tau `$, $`\sqrt{(Q_\tau Q_\tau )^2}\sqrt{2k_BTQ_\tau }`$VanZonCohen03a ; VanZonCohen03b . Thus, for $`\tau =\tau _T`$, $`\sqrt{2k_BTI^2R\tau _T}`$ is a typical value for $`Q_{\tau _T}`$ which we can insert into Eq. (16), giving

$$\frac{2k_BI^2R\tau _T}{TC_V^2}\frac{\tau _T}{\tau _F}1$$
(19)
where the time scale for the heat fluctuations is $`\tau _F=TC_V^2/2k_BI^2R=(C_V/2k_B)\tau _H`$. For an $`N`$ atomic solid at room temperature, for which Dulong-Petit is valid, we know that $`\tau _F/\tau _H=C_V/2k_B`$ = $`\frac{3}{2}N`$. Thus, the factor between $`\tau _F`$ and $`\tau _H`$ is always bigger than one, so that Eq. (19) in fact follows from Eq. (18). As a result, if on average the circuit will not fail, then the fluctuations will not make it fail either. However, for systems below room temperature, for which quantum effects become relevant, the ratio $`C_V/k_B`$ can be much smaller than one (i.e. Debye’s law), opening up the possibility that in that case the requirement on $`\tau _F`$ in Eq. (19) could be stricter than that on $`\tau _H`$ in Eq. (18).

Although reassuring, Eq. (19) concerns only typical fluctuations. However, one also needs to be concerned with large fluctuations. Compared to the Gaussian distributed work (or power) fluctuations, large fluctuations for heat are much more likely due to the exponential tails of its distribution function. Furthermore, if the condition in Eq. (16) is not met, temperature variations will occur and the theory then needs to include a coupling to the heat diffusion equation. This could be relevant for future experiments.

In conclusion, we studied work and heat fluctuations in electric circuits using analogies to Brownian systems with non-universal additions to the Nyquist noise-Brownian motion analogy. Our analogy links the work and heat fluctuations in a parallel RC circuit to those in a Brownian system for which the work and heat fluctuations are known to satisfy the conventional FT and extended FT respectively. For the serial circuit, the analogy also works for the heat fluctuations, but not for the work fluctuations. However, a short calculation \[below Eq. (13)\] shows they still satisfy the conventional FT.

RVZ and EGDC acknowledge the support of the Office of Basic Engineering of the US Department of Energy, under grant No. DE-FG-02-88-ER13847.


text:
# 𝑆𝑈(3)_𝐹 Meson Mass Formula from Random Phase Approximation

## Acknowledgements

This project was supported by the National Natural Science Foundation of China under Grant 10375003, Ministry of Education of China, FANEDD and SRF for ROCS, SEM.