text:
# On The Spectral Energy Dependence of Gamma-Ray Burst Variability

## 1. Introduction

Until a few years ago, gamma-ray bursts (GRBs) were known only as brief, intense flashes of high-energy radiation, with no observable traces at other wavelengths. The study of fading afterglows has enabled the measurement of redshift distances and the identification of host galaxies, establishing that GRBs are extremely luminous events detectable to much larger distances than quasars or galaxies. Consequently, GRBs can provide novel information about early epochs in the history of the universe.

Present distance estimates - which rely on optical line features in the afterglow spectrum or emission lines in the spectrum of the host galaxy - are relatively rare. There are now $`20`$ GRBs with optical spectroscopic redshifts (all in the range 0.43 $`z`$ 4.5).111See Jochen Greiner’s page at http://www.aip.de/ jcg/grb.html for an excellent compilation of information on GRBs with redshifts. Hence, until recently it appeared that using GRBs to map the high-$`z`$ universe would have to wait for dedicated localization and follow-up programs like Swift222http://swift.gsfc.nasa.gov and NGST333http://ngst.gsfc.nasa.gov. However, the discovery of two recent correlations between the degree of variability of the $`\gamma `$-ray light curve and the GRB luminosity (Ramirez-Ruiz & Fenimore 1999; Fenimore & Ramirez-Ruiz 2001), and between the differential time lags for the arrival of burst pulses at different energies and the luminosity (Norris, Marani & Bonnell 2000) offer the possibility of deriving independent estimates of the redshift of a GRB. Interestingly, in a large sample of BATSE bursts, the time lags for the arrival of burst pulses at different energies and the degree of variability appear to be strongly related (Schaefer, Ming & Band 2001), lending credence to each correlation. While these correlations are still tentative, they seem to be a natural consequence of the variation in energy per unit solid angle (or bulk Lorentz factor) of the emitting region (Salmonson 2000; Ioka & Nakamura 2001; Kobayashi, Ryde & MacFadyen 2001; Plaga 2001; Ramirez-Ruiz & Lloyd-Ronning, 2002). That is, an increase in energy per unit solid angle through, for example, an increase in the relativistic source expansion velocities can lead to more luminous bursts as well as shorter observed timesales in accordance with the observed correlations (see Ramirez-Ruiz & Lloyd-Ronning 2002).

Here we show that - when transformed to the cosmological rest frame using redshifts derived either from spectroscopic observations or from the luminosity-variability ($`L-V`$) relation - bursts with highly variable light curves have greater typical peak energies. The paper is organized as follows: In §2, we present the intrinsic peak energy-variability correlation. In §3, we discuss possible observational selection effects that may affect this correlation and show that - even when assuming the most conservative, severe data truncation - the correlation still holds to high statistical significance. In §4, we briefly discuss how a wide variety of burst phenomenology may be attributable to the existence of this correlation (in conjunction with the $`L-V`$ relation). We also suggest that this result not only provides useful insight into the physics of the GRB mechanism, but also may support the validity of the $`L-V`$ relation as a reliable luminosity indicator. Finally, we present our conclusions in §5. Throughout our analysis, we assume $`H_0=65\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$, a matter density $`\mathrm{\Omega }_{\mathrm{matter}}=0.3`$, and a vacuum energy density $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$.

## 2. GRB Spectra and Variability

GRB temporal profiles are so enormously varied and complicated that, at first sight, their behavior obeys no simple rule. Many bursts have a highly variable temporal profile with a timescale of variability that is significantly shorter than the overall duration. Several studies have suggested the possibility of relating properties of the time structure with the burst luminosity (Stern, Poutanen & Svennson 1997; Beloborodov et al. 2000; Norris et al. 2000; Fenimore & Ramirez-Ruiz 2001; Reichart et al. 2001). In particular, Fenimore & Ramirez-Ruiz (2001) explored the possibility of using the “spikiness” of the time structure, combined with the observed flux, to obtain distances, much as the Cepheid relationship gives distances from the pulsation period. Several hundred long and bright bursts were amenable for their analysis, producing a large sample of events with derived redshifts and luminosities. Besides using this large sample to understand both the intrinsic GRB luminosity function and GRB formation rate (e.g. Fenimore & Ramirez-Ruiz 2001; Lloyd-Ronning, Fryer & Ramirez-Ruiz 2002), these estimates offer the possibility of studying the physical nature of the luminosity indicator itself. To that effect, we investigate the dependence of the burst spectra on variability. From the set of 220 bright, long, BATSE bursts that Fenimore & Ramirez-Ruiz (2001) analyzed, we use all 159 that have time resolved fits from 16 channel spectral data (R.S. Mallozzi444Deceased., private communication). The observed spectra are phenomenologically well characterized by the “Band” function (Band et al. 1993), defined by a low-energy spectral index, $`\alpha `$, a high-energy spectral index, $`\beta `$, and the peak of the $`\nu F_\nu `$ distribution, $`E_p`$555The parameter $`E_p`$ corresponds to the peak of the spectrum in $`\nu F_\nu `$ only if $`\beta `$ is less than -2. Otherwise, the spectral peak is given by the lower boundary in energy of the high-energy power-law component characterized by $`\beta `$ (Preece et al. 2000). (at which the source is observed to emit the bulk of its luminosity). We use the spectral parameters from the time when each burst’s photon flux is maximum (i.e.the burst’s “peak”), but find qualitatively similar results if time averaged spectra are used. Table 1 lists the data we have used in our analysis.

Figure 1 shows the GRB peak energy in the cosmological rest frame, $`E_p^{}=E_p(1+z)`$, versus the observed variability, $`V`$. The filled circles are BATSE bursts with secure redshifts, high-resolution light curves and resolved spectral fits (Fenimore & Ramirez-Ruiz 2001; Jimenez, Band & Piran 2001), while the open circles are bursts with redshift estimates derived from the $`L-V`$ indicator itself. We find a significant ($`5\sigma `$) positive correlation between $`E_p^{}`$ and $`V`$ that extends for about two orders of magnitude; this correlation - taken at face value - can be parameterized as $`E_p^{}V^{0.8\pm 0.2}`$. However, we caution that the quantitative estimate of this correlation can be affected by selection effects (indicated by the shaded regions and solid line in Figure 1). We discuss this in detail in §3 below. Of the bursts in our sample, GRB 980425 is unique because of its possible association with SN 1998bw (Galama et al. 1998). The fact that it is consistent with the observed trend (see Figure 1) suggests that this event and the cosmological bursts may share a common physical origin. This speculation is made more intriguing by a recent discovery that, at least in some bursts, a supernova may be involved (Bloom et al. 1999; Reichart 1999; Lazzati et al. 2001) which may have contributed to an otherwise unexplained bump and reddening in the optical light curve (but see Esin & Blandford 2000 and Ramirez-Ruiz et al. 2001 for alternative explanations).

It is also important to note that each burst’s redshift used to transform the peak energy into the cosmological rest frame is derived from the luminosity indicator itself (apart from the bursts with secure $`z`$ that are used to calibrate the correlation). As more GRBs with independent spectroscopic redshifts are obtained, the existence of this (and the $`L-V`$) correlation will be more definitively tested. So far, however, those bursts with secure redshift estimates seem to fall well along this trend as seen in Figure 1. If we fit a power law to the correlation for just the 7 (or 8) bursts with measured redshifts, we find $`E_p^{}V^{0.75\pm 0.3}`$ and $`E_p^{}V^{0.45\pm 0.15}`$, including and excluding GRB 980425 respectively. The existence of the $`E_p^{}-V`$ correlation in our sample of 159 GRBs, therefore, may provide some confidence in the validity of the redshifts derived from the $`L-V`$ luminosity indicator.

As an illustration, in Figure 2 we show a histogram of the rest frame peak energy $`E_p^{}`$, for those 159 bursts which have spectral fits and redshifts from the $`L-V`$ relation. The superposed dotted histogram shows the observed peak energy for reference. Note that the distribution of $`E_p^{}`$ peaks at about $`1\mathrm{M}\mathrm{e}\mathrm{V}2m_ec^2`$. If one believes the redshifts from this luminosity indicator, then we can use this intrinsic distribution of GRB spectral peak energies to gain insight into the relevant particle acceleration processes and emission mechanisms present in GRBs. These possibilities and their consequences for the predicted prompt and afterglow emissions are investigated in Ramirez-Ruiz & Lloyd-Ronning (2002). One should keep in mind that the observed $`E_p`$ distribution in Figure 2 is for a limited sample of bright bursts observed by BATSE and in particular does not include the increasing number of bursts with values of $`E_p`$ that fall below the BATSE threshold, which possibly account for up to 1/3 of all GRBs (e.g. Kippen et al. 2001, Heise et al. 2001). As of yet, these so-called “X-ray Flashes” (XRFs) have no quantitative variability measurements (due to their very low fluxes, at least in the BATSE data) although it has been qualitatively claimed that they exhibit rapid variation in their time profiles, representative of the “typical” BATSE population (Kippen et al. 2001). It will be interesting to see if these bursts follow the trend exhibited in Figure 1. We note that all analysis to date of these bursts indicates that they tend to at least marginally exhibit the same trends as the bulk of the BATSE bursts (see, e.g., Kippen et al. 2001). However, once variability measurements of these bursts are made, this correlation should be re-examined, including this sample.

## 3. The Role of Selection Effects

We realize some of the limitations that are intrinsic to our procedure. Before drawing any conclusions from the $`E_p^{}-V`$ correlation, it is essential to understand the role selection effects play in determining this trend, as well as the uncertainties imposed by the scatter in the $`L-V`$ relation. We discuss each in turn below.

### 3.1. Truncation due to Flux Limit

The $`L-V`$ sample is a selected sample of bursts above a flux threshold of $`1.5\mathrm{ph}/\mathrm{cm}^2/\mathrm{s}`$, and a duration threshold $`>20`$ s. This latter selection criterion is unlikely to play an important role in our analysis since $`L`$ and $`E_p^{}`$ are relatively independent of duration. On the other hand, the flux threshold has the effect of causing a very strong truncation in the $`L-z`$ plane, in the sense that low luminosity bursts at high redshift are not “observed” (see Figure 2 in Lloyd-Ronning et al. 2001). Because $`L`$ and $`V`$ are correlated, this translates to a truncation in the $`V`$-$`z`$ plane. In other words, the flux selection criterion excludes bursts with low variability at high redshift. We would like to understand where these observationally “missing” bursts might fall in the $`E_p^{}-V`$ plane. Now, because our observed $`E_p`$ distribution is relatively narrow and uncorrelated with redshift, when we transform into the cosmological rest frame $`E_p^{}=E_p(1+z)`$, we find (on average) higher $`E_p^{}`$ values at higher redshifts. Therefore there is a selection against low variability, high peak energy bursts (the upper left of the $`E_p^{}-V`$ plane). We investigate the role - if any - this selection effect plays in producing the observed correlation.

A quantitative formulation of the truncation in the $`E_p^{}-V`$ plane from the luminosity-redshift selection effect is not straightforward (primarily because the scatter in the observed $`L-V`$ relation does not allow for an unambiguous change of variables between $`L`$ and $`V`$). However, we can make an estimate of this truncation if we assume a perfect correlation between $`L`$ and $`V`$, $`LV^{3.3}`$ (this is the best-fit relation derived by Fenimore & Ramirez-Ruiz 2001). From the flux threshold selection criterion, a truncation is produced in the luminosity-redshift plane: $`L_{lim}f_{p,lim}d_z^2`$, where $`d_z=(1+z)_0^z𝑑z/\sqrt{\mathrm{\Omega }_\mathrm{\Lambda }+\mathrm{\Omega }_m(1+z)^3}`$, and $`f_{p,lim}`$ is the peak flux limit that Fenimore and Ramirez-Ruiz chose - $`1.5`$ ph/cm2/s. We can then estimate the limiting variability as a function of $`E_p^{}`$ by changing variables from $`L`$ to $`V`$ (using $`LV^{3.3}`$) and the fact that $`E_p^{}/E_p=(1+z)`$. We find $`V_{\mathrm{lim}}(E_p^{})[E_p^{}_1^{E_p^{}/E_p}𝑑x(\mathrm{\Omega }_\mathrm{\Lambda }+\mathrm{\Omega }_\mathrm{m}x^3)^{1/2}]^{2/3.3}`$. This estimate of the truncation is in fact quite shallow as shown by the solid line in Figure 2. I has no effect on the quantitative results reported in section 2 (which did not account for selection effects). This point is further illustrated in Figure 3, which shows the intrinsic peak energy $`E_p^{}`$ vs. V for increasing flux bins (the squares, triangles, pentagons, hexagons, and horizontal lines are for $`f_p=1.5-2.2`$ ph/cm2/s, $`f_p=2.2-2.7`$ ph/cm2/s, $`f_p=2.7-4.0`$ ph/cm2/s, $`f_p=4.0-7.0`$ ph/cm2/s, and $`f_p>7.0`$ ph/cm2/s, respectively). If the flux limit caused a severe bias against observing bursts in the upper left hand corner of the $`E_p^{}-V`$ plane, we expect to see these points migrate toward the right lower corner of the $`E_p^{}-V`$ plot, with increasing flux bin. This is clearly not the case, and in fact the points are scattered throughout one another in the plot.

In the analysis that follows, we in fact take the most conservative approach by imposing a steeper, more severe truncation than what the analytical estimate in the preceeding paragraph predicts. We feel this more conservative approach not only places our results on much firmer statistical ground, but also implicitly accounts for the effects of scatter in the $`E_p^{}-V`$ plane.

### 3.2. Accounting for the Selection Effects

As mentioned above, it is important to quantify the role of selection effects in the data. Fortunately, there are firmly established non-parametric statistical methods that have been developed to deal with precisely such selection effects (e.g. Lynden-Bell 1971; Efron & Petrosian 1992). These techniques use a well-defined truncation criterion (and the assumption that the observed sample is the most likely to be observed) to estimate the correlation between (and underlying parent distributions of) the relevant variables. For each data point indexed by $`i`$, an “eligible set” is defined based on those points that fall within the observational limits of the $`i`$’th data point at hand. This amounts to making a truncation parallel to the axes for each data point; a weight is then assigned to each point given the number of points in its eligible set. For example, an eligible set for the $`i`$’th data point in our sample consists of all bursts indexed by $`j`$, where $`V_{lim,j}<V_i`$ and $`V_j>V_i`$. Correlations are then computed via non-parametric rank statistics (such as a Kendell’s $`\tau `$ test), where rank comparisons are made only among those within the eligible sets. Further details of these techniques can be found in Efron and Petrosian (1992) and in the appendix of Lloyd, Petrosian, & Mallozzi (2000).

As mentioned in §3.1, to estimate how the peak flux selection criterion affects our results, we take the most conservative view that all of the lack of bursts in the upper left corner is the result of an observational selection effect. We in fact employed two estimates of this truncation for our data - shown by the two shaded regions in Figure 1. The first (dark shaded region) fits as closely to the data as possible without eliminating any points; the latter, more severe truncation (light shaded region) artificially eliminates “high scatter” points, to allow the truncation to fit tightly to the majority of the data. For the first truncation, we find a significant $`>5\sigma `$ correlation between $`E_p^{}`$ and $`V`$. The functional form of this correlation (accounting for the truncation of course) can be expressed as $`E_p^{}V^{0.8\pm 0.15}`$. For the second truncation, we again find a $`>5\sigma `$ correlation between $`E_p^{}`$ and $`V`$, which can be parameterized as $`E_p^{}V^{0.7\pm 0.15}`$. The existence of the correlation, therefore, is likely to primarily be due to the lack of bursts with high variability and low $`E_p^{}`$ (i.e. the lower right hand corner of Figure 1) which we believe is real and is not a result of any selection effect. We again emphasize that the flux limit most likely imposes much less severe truncation than we have assumed here (see solid line in Figure 1).

### 3.3. Consequences of Uncertainties in the $`L-V`$ Relation

We have also computed the correlation between $`E_p^{}`$ and $`V`$ given the luminosities and redshifts derived from the upper and lower limits to the fitted $`L-V`$ relation. This correlation is parameterized as $`LV^\beta `$ with $`\beta =2.2,3.3`$, and $`5.8`$ for the lower limit, best-fit, and upper limit respectively (see Fenimore & Ramirez-Ruiz 2001). As for the best fit relation (see §3.2 and Figure 1), we chose two truncations for each data set – the first fitting as closely as possible to the data without eliminating any data points (truncation 1), and the second eliminating “high scatter” points (truncation 2). In all cases, we find a $`>4\sigma `$ correlation between $`E_p^{}`$ and $`V`$. For $`\beta =2.2`$, the correlation can be parameterized by $`E_p^{}V^{0.5\pm 0.15}`$ (truncation 1) and $`E_p^{}V^{0.4\pm 0.15}`$ (truncation 2). For $`\beta =5.8`$, the parameterization is $`E_p^{}V^{1.15\pm 0.15}`$ (truncation 1) and $`E_p^{}V^{1.05\pm 0.15}`$ (truncation 2). All of these results are summarized in Table 2.

It is also important to notice that we are relying on redshifts that are obtained from the $`L-V`$ relation, which we have extrapolated to higher variabilities than those from which the correlation was derived. We have therefore computed the correlation between $`E_p^{}`$ and $`V`$ for only those bursts in our sample which fall in the variability range corresponding to the observed bursts that were used to calibrate the $`L-V`$ correlation. Again, accounting for selection effects, we still find a highly signficant correlation ($`>5\sigma `$) between $`E_p^{}`$ and $`V`$ for bursts in this limited range (both including and excluding the variability of GRB 980425 in defining the lower limit of our range), with $`E_p^{}V^{0.9\pm 0.2}`$. In addition, we have computed the correlation between $`E_p^{}`$ and $`V`$ for only those bursts in our sample which fall in the redshift range corresponding to the observed bursts that were used to calibrate the $`L-V`$ correlation, $`0.8<z<3.4`$ (note that we exclude GRB 980425; including this burst in our redshift range will increase the significance of our results). Once again, we find a $`4\sigma `$ correlation between $`E_p^{}`$ and $`V`$ for the $`80`$ bursts that made this redshift cut, with $`E_p^{}V^{0.75\pm 0.2}`$. We conclude from this analysis that if the $`L-V`$ relation holds true (even over a limited range), then highly variable bursts have greater intrinsic peak energies.

## 4. On the relationship between temporal and spectral structure in GRBs

The existence of the $`E_p^{}-V`$ correlation also implies the (not surprising) fact that the GRB luminosity and intrinsic peak energy are correlated, as directly follows from the $`L-V`$ and $`E_p^{}-V`$ relations. We indeed find such a correlation in our data at high ($`>5\sigma `$) statistical significance. This result is qualitatively consistent with the findings of Lloyd, Petrosian & Mallozzi (2000), who suggested that the observed hardness-intensity correlation in GRBs is an intrinsic (and not cosmological) effect - namely, a correlation between GRB rest frame peak energy and luminosity.

Moreover, several mutually reinforcing trends have been found in the past that support the validity of our results. One key trend involves the tendency for pulses or “peaks” in GRB time histories to be narrower at higher energies. This was first noted by Fishman et al. (1992), Link, Epstein & Priedhorsky (1993) and quantitatively explored by Fenimore et al. (1995). The latter showed that the average pulse width has a power-law dependence on energy with an index of about -0.4. A visual inspection of the pulses fitted to GRBs by Norris et al. (1996) shows that the low-amplitude pulses (within a single burst) tend to be wider. Finally, Ramirez-Ruiz & Fenimore (2000) found a quantitative relationship between pulse amplitude and pulse width: the smaller amplitude peaks tend to be wider with the pulse width following a power law with an index of about -2.8. Therefore, it is not surprising that we find that more variable profiles contain a larger number of high energy pulses, which are intrinsically narrower and brighter.

## 5. Conclusions

We have shown there exists a correlation between the characteristic photon energy in the cosmological rest frame and the gamma-ray burst variability, as well as the GRB luminosity. While these correlations are still tentative, it is reassuring that they are consistent with several mutually reinforcing trends found between spectral and temporal properties in a diversity of GRB profiles. These relationships can help shed light on the relevant physical mechanisms responsible for the observed properties of a gamma-ray burst - namely the structure of the ultra relativistic outflow, the microphysics of shock acceleration, and the magnetic field generation (see Ramirez-Ruiz & Lloyd-Ronning 2002).

We would like to thank D. Band, E. Fenimore, D. Lamb and V. Petrosian, and J. Salmonson for very useful discussions and comments. We also thank the referee for comments that led to improvements of this paper. ERR thanks CONACYT, SEP and the ORS for support.


format: text

text:
# Astrophysical magnetic field reconstruction and spectroscopy with ultra high energy cosmic rays

## 1 Introduction

The origin and nature of ultra high energy cosmic rays (UHECRs), after several decades of study, is still puzzling . If they are charged particles, such as protons or heavier nuclei, their trajectories are affected by the magnetic fields along their path from the sources to the observer. For energies below $`Z\times 10^{18}`$eV ($`Ze`$ is the CR electric charge) they are confined by Galactic magnetic fields. For larger energies the galactic field is not strong enough to confine them, which combined with the lack of any significant excess from the galactic disk suggests that UHECRs are most probably of extragalactic origin. However, due to the effects of the magnetic fields present along their path, they do not point to their birth places, complicating the task of identifying their sources. The small angular scale clustering observed in the AGASA data  already hints to the existence of UHECR point sources , but the small present statistics does not allow to solve the source identification question. With the next generation of detectors, like AUGER  and EUSO , larger clusters are expected, having each several dozens of events. If this were the case, it would not only become possible to reconstruct the source positions, but also the data could be used to obtain information about the magnetic field along the CR path .

In this paper we discuss some strategies which can be adopted to reconstruct the main parameters of Galactic magnetic fields from a set of events originating from an extragalactic point source, and test them using simulated data. The knowledge of the magnetic field responsible for the deflections allows to better reconstruct the locations of the sources and, moreover, can also be helpful to do spectroscopy of CRs, making possible the measurement of their charges. The method we develop here profits from several specific features imprinted on the events by the magnetic fields, which not only deflect the trajectories but also lead to strong lensing phenomena, including the formation of multiple images and energy dependent magnifications or demagnifications of the CR fluxes .

## 2 Magnetic lensing

The magnetic field in the Galaxy is known to have a regular and a turbulent component . The regular field in the disk follows the spiral arms structure, with reversals taking place from arm to arm. The local value is $`B_{reg}2\mu `$G. The possible existence of an extended halo field (of few kpc scale height) is not settled down yet. The typical deflection in the direction of propagation of a CR particle of charge $`Ze`$ and energy $`E`$ produced by a regular magnetic field is

$$\delta 8.1^{}\frac{40\mathrm{EeV}}{E/Z}\left|_0^L\frac{\mathrm{d}𝐬}{3\mathrm{kpc}}\times \frac{𝐁}{2\mu \mathrm{G}}\right|,$$
(1)
where $`1\mathrm{EeV}=10^{18}`$ eV.

The turbulent galactic magnetic field component has a root mean square amplitude larger than the regular field, $`B_{rms}1`$$`3B_{reg}`$, and the largest turbulence scale is $`L_{max}100`$ pc, while the minimum one is believed to be much smaller ($`L_{min}L_{max}`$). It is often modeled as a Gaussian random field with a Kolmogorov spectrum  (i.e. such that the energy density satisfies d$`E/\mathrm{d}kk^{5/3}`$). The mean deflection of particles moving in the turbulent field vanishes, while the root mean square deflection is given by

$$\delta _{rms}1.4^{}\frac{40\mathrm{EeV}}{E/Z}\frac{B_{rms}}{4\mu \mathrm{G}}\sqrt{\frac{L}{3\mathrm{kpc}}}\sqrt{\frac{L_c}{50\mathrm{pc}}},$$
(2)
where $`L_c`$ is the coherence length. For a narrow band ($`L_{min}L_{max}`$) spectrum, one has $`L_cL_{max}/2`$, while for a broad band ($`L_{max}L_{min}`$) Kolmogorov spectrum it is $`L_cL_{max}/5`$.

We have written the random and regular field deflections in terms of the characteristic values of the Galactic fields, thus the numerical value of the deflections are typical expected values. We see that the regular field generally produces the dominant effect on the deflection, which leads to a coherent displacement of the apparent position of the source that is inversely proportional to the energy. Compared to this, the random field produces rather small changes in the arrival direction of the CRs. However, as the charged particles propagate through the magnetic field, they are not only deflected, but their flux is also focused or defocused due to the differential deflections of neighboring paths , with the effects being larger at smaller energies. As lower energies are considered, CR particles from one source can have more than one path leading to the detector, each traversing an uncorrelated magnetic field patch, and this leads to the appearance of multiple images of the source. This typically occurs when $`\delta _{rms}L_c/L`$ , what allows to define a critical energy $`E_c`$, around and below which the formation of multiple images is very likely, through111 The regular field can also lead to multiple images, with its characteristic scale of homogeneity (few kpc) playing the role of the coherence length.

$$\delta _{rms}\frac{E_c}{E}\frac{L_c}{L}.$$
(3)
Its numerical value is given by222 If the amplitude of the turbulent field were not constant, but had instead some modulation along the Galaxy, the factor $`B_{rms}\sqrt{L}`$ in eq. (2) should be replaced by $`|_0^LdsB_{rms}^2(s)|^{1/2}`$, and $`B_{rms}L^{3/2}`$ in eq. (4) should be replaced by$`\left[3_0^Ldxx_0^xdsB_{rms}^2(s)\right]^{1/2}`$.

$$E_cZ60\mathrm{EeV}\frac{B_{rms}}{4\mu \mathrm{G}}\left(\frac{L}{3\mathrm{kpc}}\right)^{3/2}\sqrt{\frac{50\mathrm{pc}}{L_c}}.$$
(4)
The new images appear in pairs, at an angular distance $`L_c/L`$ from the position of the original source image and with a large magnification of their fluxes. The peak in the amplification when new images appear leads to small angular scale clustering and to an excess of events in an energy bin close to $`E_c`$ .

The overall picture can be summarized as follows: the magnetic fields present along the CR trajectories change the mapping between the observed arrival directions and the source ones. At the highest energies, the effects are small and hence the mapping is close to the identity. Going down in energy, deflections arise due to the turbulent and the regular magnetic field components. The turbulent component leads to small rms deflections, but it gives rise to a network of lines in the sky around which the fluxes of potential sources would be significantly magnified. This network corresponds to the directions where folds in the mapping (caustics) develop. The caustics appear at energies near (and below) $`2E_c`$, giving rise to multiple images of the sources lying close to those directions . The regular field, which is responsible for larger deflections, leads to a global displacement of this network of caustics relative to the source directions. As a result of this drift caused by the regular field, a given CR source will experience successive magnifications and demagnifications of its flux at different energies caused by the motion of the caustic network across the sky. For decreasing energies, the continuous appearance of new images also leads to additional peaks in the spectrum, and strong lensing effects associated to the regular field component will also show-up.

## 3 Magnetic field reconstruction

The properties of the magnetic deflections and lensing effects can be used to determine the magnetic field parameters from the study of clustered UHECR events having their origin in isolated point-like sources. For this purpose one should start by defining the regions containing the excess events (which have to be larger than at least a few times the angular resolution of the detector but yet not too large to avoid excessive background events coming from other independent nearby sources). Different strategies must be considered, depending on the relative strength of the effects due to regular and random components, and in what follows we develop several ideas to perform this reconstruction.

Consider first the likely case in which the deflection caused by a regular component (either galactic or extragalactic) is larger than the rms deflection imprinted by the turbulent galactic field. In this case, the overall angular motion of the source images as a function of energy can be used to reconstruct the integral along the line of sight to the source of the perpendicular component of the regular magnetic field. In order to do this, one first finds the overall direction of motion of the images caused by the regular field by fitting a straight line ($`\alpha _y=a+b\alpha _x`$) to the event coordinates $`\alpha _{x,y}`$. The integral along the line of sight of the perpendicular component of the regular magnetic field will just be at a right angle with respect to this line, and its absolute value can be obtained choosing now the coordinates $`\alpha _{}`$ (in the direction along the deflections produced by the regular field) and $`\alpha _{}`$ (orthogonal to it), and plotting $`\alpha _{}`$ vs. $`1/E`$ for the set of selected events. One should have, ignoring at this step the random field component, that

$$\alpha _{}^{}{}_{i}{}^{}\alpha _0+\frac{Z_i}{E_i}K_{}$$
(5)
with

$$K_{}\pm |e_0^Ld𝐬\times 𝐁^{reg}|\pm 1.3^{}(40\mathrm{EeV})\left|_0^L\frac{\mathrm{d}𝐬}{\mathrm{kpc}}\times \frac{𝐁^{reg}}{\mu \mathrm{G}}\right|.$$
(6)
If the CRs have all identical electric charge $`Z`$, the events would fit nicely to just one straight line, $`\alpha _{}=\alpha _0+c/E`$. The slope $`c`$ would give the value of $`ZK_{}`$, which essentially provides the magnitude and direction of the regular magnetic field orthogonal to the l.o.s. multiplied by its relevant scale length, modulo the CR electric charge. The value of $`\alpha _0`$ just gives the original location of the source ($`\alpha _{}=\alpha _0`$, $`\alpha _{}=0`$). After the source location is estimated, it proves convenient to change variables to $`\alpha _{}\mathrm{sign}(K_{})(\alpha _{}\alpha _0)`$, defined such that it tends to zero at large energies and increases as $`E`$ decreases. In what follows we will use the name $`\alpha _{}`$ just to refer to this shifted variable.

We illustrate this procedure using simulated data, obtained through a ray shooting technique, in which a large number of antiparticles were thrown isotropically from Earth, their trajectories numerically integrated across a distance $`L`$ within a homogeneous magnetic field plus a turbulent component, and were finally recorded only if they pointed within a small cone (of aperture $`\delta _{rms}/10`$) from the direction to a fixed ‘source’. The source was assumed to inject CRs with a differential energy spectrum $`\mathrm{d}N/\mathrm{d}EE^{2.7}`$. About 100 simulated events with energy above $`40`$ EeV were selected, assuming the detector efficiency to be energy-independent. We display the ‘theoretical’ results from the simulation (with black dots) as well as more realistic ones (open circles) obtained by adding noise to the energy (corresponding to a detector energy resolution of 10%) and to the angles (corresponding to an angular resolution of $`0.2^{}`$ in each direction).

In the case illustrated in Figure 1 the homogeneous magnetic field had strength $`2\mu `$G along the $`\alpha _x`$ direction, the turbulent component had $`B_{rms}=4\mu `$G with just one turbulence scale $`L_{max}=L_{min}=100`$ pc (and thus $`L_c=50`$pc), and the distance traversed by the CRs (assumed to be protons) within these fields was $`L=3`$kpc. The source actual position was in the direction ($`\alpha _x=0,\alpha _y=0`$). The first panel displays the simulated data, with the straight line being the reconstructed overall direction of motion of the images and the star being the reconstructed source location. The second panel plots instead $`\alpha _{}`$ vs. $`1/E`$, while the third one will be discussed farther below.

In Figure 2 we display the results for another simulation with the same parameters for the regular field and same root mean square amplitude of the turbulent field but with a Kolmogorov spectrum with $`L_{max}=10L_{min}=100`$ pc (corresponding to $`L_c=25`$ pc). The result for the best fit of the numerical reconstruction was for both cases $`Z|_0^Ld𝐬\times 𝐁|6.1\mu `$G kpc, in good agreement with the values used in the simulations.

In the case of a mixed composition, the fit of $`\alpha _{}`$ vs. $`1/E`$ to just one straight line should fail (excessive $`\chi ^2`$), since one expects to see the superposition of several straight lines, all with the same value of $`\alpha _0`$ but different slopes, with ratios fixed by the ratio of the different electric charges $`Z_i`$ involved333The value of $`\alpha _0`$ in the case of mixed composition was just determined by eye, looking for the value of $`\alpha _{}`$ towards which the different lines converged.. In this case one can separate the events in groups, above and beyond limiting slopes, and then fit each group independently to an expression of the form given in Eq. (5). This allows to obtain information on both $`K_{}`$ and the ratio of charges $`Z_i`$ involved, besides the relative amounts of different nuclei from the relative number of events on the different lines. Notice that the reconstructed fraction of different nuclei reflects the composition upon arrival to Earth, which may differ from the fraction injected by the source within the same energy range because the flux magnification due to magnetic lensing depends on the combination $`E/Z`$. Thus, the flux of different nuclei is magnified by a different factor at a given energy .

In Figure 3 we exemplify this procedure for the same values of the fields and Kolmogorov spectrum with 60% of protons and 40% of Lithium ($`Z=3`$). (The actual fraction of Lithium injected by the source was 50% in this case.) The two slopes in $`\alpha _{}`$ vs. $`1/E`$ are clearly separated and the reconstructed values are $`Z_1|_0^Ld𝐬\times 𝐁|5.9\mu `$G kpc, and $`Z_2|_0^Ld𝐬\times 𝐁|18.6\mu `$G kpc, implying that $`Z_2/Z_13.1`$.

In Figure 4 we show the analogous results for 80% of protons and 20% of Helium, for which we obtain $`Z_1|_0^Ld𝐬\times 𝐁|6.2\mu `$G kpc, and $`Z_2|_0^Ld𝐬\times 𝐁|12.3\mu `$G kpc, implying that $`Z_2/Z_12`$. The reconstructed direction of the regular magnetic field was off by only a few degrees in all cases. Let us mention that we performed several other simulations with different realizations of the random field component, obtaining reconstructed values of similar quality.

Having determined the overall motion of the images due to the regular magnetic field component, we turn to analyze the effects of the turbulent field. At high energies the main effect is the formation of multiple images, with their associated peaks manifesting in the observed flux. The critical energy of this process, $`E_c`$ (given by eq. (4)), can be determined from the location of the first lensing peak, which is expected to lie at energies somewhere between $`2E_c`$ and $`E_c/2`$. We then take bins in $`\alpha _{}`$ covering the range of angles where clustered events are observed above a certain energy threshold. The bins should be larger than the angular resolution of the detector and such that a sizeable number of events per bin result (taking e.g. ten angular bins if the number of clustered events is of the order of one hundred). Now one can plot the number of events in each bin vs. $`\alpha _{}`$. If there were only deflections but no flux (de)magnification due to lensing effects, the regular field would just lead to $`\mathrm{\Delta }N/\mathrm{\Delta }\alpha _{}\alpha _{}^\beta `$, where $`\beta `$ characterizes the slope of the CR energy spectrum, with $`\mathrm{d}N/\mathrm{d}EE^{2\beta }`$.

The advantage of plotting $`\mathrm{\Delta }N`$ in terms of $`\alpha _{}`$ (and not e.g. vs. $`E^1`$) is that regardless of the CR composition the first lensing peak always appears at the same value of $`\alpha _{}`$. Referring to this value as $`\alpha _{}^{peak}`$, one has then (up to an overall factor of two)

$$\frac{E_c}{Z}\frac{|K_{}|}{\alpha _{}^{peak}}$$
(7)
The third panels of Figures 1-4 illustrate the results of this method applied to the simulated data. In Figure 1 the first lensing peak associated to multiple image formation appears at $`\alpha _{}^{peak}6.4^{}`$ and corresponds to the sharp peak. It leads to a critical energy value $`E_c/Z49`$EeV, in good agreement with the parameters used in the simulation (which correspond to $`E_c60`$ EeV). A shallower and more symmetric excess of events is also evident at smaller deflection angles, $`\alpha _{}^{peak}2.5^{}`$, and it is due to a high magnification region crossing the source location at high energies (larger than $`2E_c`$, when folds were not yet formed).

In Figure 2 the first peak appears at $`\alpha _{}^{peak}4.3^{}`$, corresponding to $`E_c/Z72`$ EeV (while the simulation parameters lead to $`E_c/Z85`$ EeV). For the mixed composition cases (Figures 3 and 4) the corresponding values are $`E_c/Z85`$ EeV and $`E_c/Z87`$ EeV respectively.

The angular distribution of the events can also be used to estimate $`L_c/L`$ ($`L`$ here is the distance traversed along the turbulent component, which needs not be the same as that across the regular one). There are two independent ways to do this. One is based on the fact that at energies $`E_c`$ the typical angular separation in the network of caustics of the turbulent magnetic field is $`4L_c/L`$ . Since the effect of the regular field induced deflections is to move these caustics, relative to the source direction, the typical separation in $`\alpha _{}`$ between consecutive lensing peaks is also equal to $`4L_c/L`$. If at least two lensing peaks are observed in the data from one cluster, then $`L_c/L`$ can be estimated. This is the case in the simulations presented above. In the case with just one wavelength depicted in Figure 1 (which has $`L_c=50`$ pc), the separation between peaks is of order $`4^{}`$, as it can be seen in the third panel, leading to $`L_c/L0.017`$. In the Kolmogorov spectrum cases (for which $`L_c25`$ pc) depicted in Figures 2, 3 and 4, the separation among peaks is between $`2^{}`$ and $`3^{}`$, leading to $`L_c/L`$ values in the range 0.008–0.013. All these values are within 50% of the actual values used in the simulations.

Notice that it is likely to observe the first lensing peak rather close to $`2E_c`$, and to observe more than one lensing peak in the range of energies between $`2E_c`$ and $`E_c/2`$, if the deflection caused by the regular magnetic field component is strong enough. Indeed, the network of caustics is already well defined at energies around $`2E_c`$. Although only a small percentage of all potential source directions are sufficiently close to a caustic at energies around $`2E_c`$ to have their flux noticeably magnified, if the deflection caused by the regular magnetic field component is significant, a caustic will reach the source position at energies not far below $`2E_c`$. For instance, if the regular field deflection is at least 3 times larger than the rms deflection caused by the random field, then the source apparent position is displaced, as the energy decreases from $`2E_c`$ down to $`E_c/2`$, by more than $`4L_c/L`$. Since this is the typical separation between caustics, it is likely to produce more than one lensing peak in this energy range. Notice however that if the deflection by the regular field were too strong, the caustics may be crossed too fast, leading to a small integrated effect and thus a less noticeable peak.

Alternatively, $`L_c/L`$ can be determined from the fact that it is just the typical angular scale characterizing the separation between the images present in the first lensing peak. In this peak there are in principle one image corresponding to the principal image which is also present at higher energies, and two new images which are both in the same location of the sky and have appeared separated from the previous one by an angular scale of order $`L_c/L`$. Since this angular separation may be quite small, this approach can be used only when the angular uncertainty in the events is smaller than this scale. For the example shown in Figure 1 this separation is about $`1^{}`$, leading to an $`L_c/L`$ estimate similar to the one obtained with the previous method. For the Kolmogorov spectrum cases, that leads to smaller $`L_c/L`$, we see from the first panel in Figures 2, 3 and 4 that to apply this method in cases like these would require a quite precise angular resolution, probably better than the one achievable with AUGER ($`0.3^{}`$), but possibly within the aim of some proposed detectors .

Once $`E_c/Z`$ and $`L_c/L`$ have been estimated, we can combine them to obtain an estimate of $`B_{rms}L`$ using Eq. (4). In our examples, the recovered values fall in the range $`B_{rms}L10`$$`13\mu `$G kpc, in good agreement with the value $`12\mu `$G kpc used in the simulations.

If the deflection caused by the regular magnetic field is smaller than the rms deflection caused by the turbulent field at the same energy, the overall motion of the images does not provide reliable information about the regular magnetic field component, other than an upper bound on its integrated effect. However, it might still be possible to better reconstruct the source location by looking at the highest energy events, i.e. at energies above those at which the multiple images appear. Indeed, as long as the deflections are small ($`\delta L_c/L`$), the random field component still produces a coherent deflection of the images proportional to $`(E/Z)^1`$, and hence unless there is a significant mixture of different CR compositions there will be a clear signal allowing to reconstruct the source location. One could also extract information about the turbulent field through the spectral features imprinted by the lensing peaks, and by the angular clustering properties of the events. The critical energy $`E_c`$ can be estimated as the highest energy at which a peak is observed in the spectrum. With sufficient angular resolution, $`L_c/L`$ can be estimated from the angular separation between the images present in the first lensing peak, while information about the composition can be extracted from the observed energy ratios of the events clustered in the peaks .

These general features are illustrated in Figure 5, which corresponds to a ray shooting simulation with a small orthogonal component of the regular magnetic field (leading to $`B_{}^{reg}ds=0.3\mu `$G kpc). The random component has $`B_{rms}=4\mu `$G (and spectrum with $`L_{min}=L_{max}=100`$ pc) and $`L=3`$ kpc. Events with energies larger than 40 EeV and injection spectrum $`E^{2.7}`$ are displayed. The left panel is for CR protons alone. The highest energy events move away from the source location (here the origin) with a deflection $`E^1`$, and at energies near 60 EeV two pairs of secondary images appear. They are displaced by $`1.5^{}`$ from the principal image, hence leading to an estimated $`L_c/L0.026`$, in good agreement with the input value 0.017. For the right panel, which has a mixed composition of protons and 40% Be at injection, other secondary images appear for $`E>40`$ EeV, dispersed over a region of $`3^{}`$, and they contain exclusively the heavier Be nuclei, since the protons in them have energies below 40 EeV. In general each of the secondary images has a narrow rigidity distribution due to the lensing origin of the peaks (and hence if two different elements were contributing to an image in the considered energy range, the events would split into groups with energy ratios given by the ratios of their charges).

## 4 Conclusions

One of the most exciting perspectives for the Auger Observatory currently under construction, as well as for other future experiments, is that they may inaugurate the era of UHECR astronomy. A hint that this may be the case is the small angular scale clustering observed by AGASA, which may be the first indication that a few UHECR sources are significantly more powerful than the otherwise apparently isotropic background.

UHECR astronomy is inevitably tied to the magnetic fields along the line of sight towards the most powerful sources, if the bulk of their emission is in the form of charged particles, such as protons or heavier nuclei. Conservative estimates of galactic and intergalactic magnetic fields indicate that the deflection of protons should be small at the highest energies, around and above $`10^{20}`$ eV, so that they should truly point to their birthplace. However, even protons with energies around $`10^{20}`$ eV may have their flux significantly (de)magnified by lensing effects, for instance by the turbulent component of the magnetic field in the Milky Way. A heavier component of UHECRs would undergo deflections larger than protons at comparable energies, and would be subject to similar flux (de)magnification effects at higher energies.

Consequently, the behaviour of the transition towards UHECR astronomy at the highest energies is a very rich source of information both on the source properties as well as on the characteristics of the intervening magnetic fields. This has been exploited, for instance in , as a tool for the reconstruction of intergalactic magnetic field parameters, particularly in the case of bursting sources of UHECRs , where the energy-dependent time delays provide additional handles to perform the task.

Here we developed a strategy to reconstruct the parameters of intervening magnetic fields that profits from the strong lensing effects that take place at energies of the order of $`E_c`$, given by Eq. (4), around which it is likely to observe strongly magnified multiple images of single UHECR sources. We have shown that it is likely that a sort of ‘magnetic spectroscopy’ of UHECR sources may allow to measure the strength of the magnetic field along the line of sight, its coherence properties, and the relative amounts of nuclei with different electric charges among the observed CRs. Our analysis assumed the sources to be steady, what means here essentially that the typical time delays between different images are smaller than the timescales of emission of the sources.

The method works if the number, strength and location of sources, the UHECR composition, and the intervening magnetic field parameters are such that clusters with a large number (several dozens) of events can be identified, without excessive background from events with origin in nearby independent sources. The angular and energy resolution required is also dependent on the actual values of the parameters involved. We have shown, through illustrative examples based in simulated data, that the method is likely to be applicable for realistic values of the parameters of the galactic magnetic field. The overall angular motion of the events as a function of energy allows to reconstruct the integral along the line of sight of the perpendicular component of the regular magnetic field, modulo the CR electric charge. If this can be done for several sources, one may then even map the strength and extent of the magnetic field along different directions, and this could help to establish its overall distribution.444Notice that the information reconstructed with UHECRs is complementary to that provided by Faraday rotation measurements, which are sensitive to the magnetic field component parallel to the line of sight. The reconstructed magnetic field also provides a handle to measure the relative number of CRs with different electric charges. A ‘spectroscopic’ analysis of the peaks in the number of events as a function of this overall motion (or as a function of energy when the effect of the regular component is not strong enough) provides a measure of the critical energy $`E_c/Z`$ at which strong lensing phenomena occur. The enhanced angular clustering of the events at the energies where strong flux magnification takes place can also be used to estimate $`L_c/L`$, the ratio between the coherence length of the turbulent magnetic field component and the path length traversed by the CRs across it. $`E_c/Z`$ and $`L_c/L`$ provide a measure of $`B_{rms}L`$, independent of the CR electric charge. The fact that the magnification peaks due to magnetic lensing occur at energies which are in direct proportion to the electric charge $`Z`$ provides a second handle to grasp the CR composition. It induces a strong correlation of arrival directions at energy ratios fixed by the ratios of the CR electric charges.

Potentially, this method could also extract information about the different scales involved in the turbulence of the intervening magnetic field, if the experiments had sufficient energy and angular resolution (and, of course, statistics). If the spectrum is broad band, long and short wavelengths lead to amplification peaks at somewhat different values of $`E/Z`$ . The density of lensing peaks could thus be used to obtain information on the magnetic field spectral properties. However, in realistic experiments the very narrow lensing peaks due to short wavelengths are likely to go unnoticed.

The method may also provide information on extragalactic magnetic fields. Notice that for an extragalactic turbulent component one has an associated critical energy

$$E_cZ2\times 10^{19}\mathrm{EeV}\frac{B_{rms}}{10^9\mathrm{G}}\left(\frac{L}{10\mathrm{Mpc}}\right)^{3/2}\sqrt{\frac{\mathrm{Mpc}}{L_c}}.$$
(8)
Hence, if these magnetic fields have a large extent ($`L\mathrm{Mpc}`$) and a significant strength ($`10^9`$ G), their effects could produce noticeable signals at energies above $`10^{20}`$ eV even for protons, while for more modest fields, it would be the galactic magnetic fields the ones producing the dominant effects.

###### Acknowledgments.

Work partially supported by ANPCYT, CONICET, and Fundación Antorchas, Argentina.


text:
# 1 Motivation

## 1 Motivation

A main question in QCD is the existence of special classes of gauge fields which are responsible for confinement. In the framework of the dual superconductor, color magnetic monopoles arise as the appropriate candidates. Their investigation necessitates a projection of the non-Abelian gauge theory on the compact U(1) theory. In lattice computations it was demonstrated that monopoles account for more than 90 percent of the string tension in the confinement and they drop toward zero across the phase transition.

The study of chaotic dynamics of classical field configurations in field theory is motivated by phenomenological applications as well as by the understanding of basic principles. The role of chaotic field dynamics for the confinement of quarks is a longstanding question. Here, we analyze the leading Lyapunov exponents of compact U(1) configurations on the lattice. The real-time evolution of the classical field equations was initialized from Euclidean equilibrium configurations created by quantum Monte Carlo simulations. The U(1) gauge field is decomposed into a monopole and photon part. This way we expect to learn some details between the appearance of monopoles and the strength of chaotic behavior in lattice simulations.

## 2 Monopoles in compact quantum field theories

We have investigated $`4d`$ U(1) gauge theory described by the action

$$S\{U_l\}=\underset{p}{}(1\mathrm{cos}\theta _p),$$
(1)
with $`U_l=U_{x,\mu }=\mathrm{exp}(i\theta _{x,\mu })`$ and $`\theta _p=\theta _{x,\mu }+\theta _{x+\widehat{\mu },\nu }\theta _{x+\widehat{\nu },\mu }\theta _{x,\nu }(\nu \mu ).`$ At critical coupling $`\beta _c1.01`$ U(1) gauge theory undergoes a phase transition between a confinement phase with mass gap and monopole excitations for $`\beta <\beta _c`$ and the Coulomb phase with a massless photon for $`\beta >\beta _c`$, see Fig. 1.

We are interested in the relationship between monopoles of U(1) gauge theory and classical chaos across this phase transition. Following Refs. , we have factorized our gauge configurations into monopole and photon fields. The U(1) plaquette angles $`\theta _{x,\mu \nu }`$ are decomposed into the “physical” electromagnetic flux through the plaquette $`\overline{\theta }_{x,\mu \nu }`$ and a number $`m_{x,\mu \nu }`$ of Dirac strings passing through the plaquette

$$\theta _{x,\mu \nu }=\overline{\theta }_{x,\mu \nu }+2\pi m_{x,\mu \nu },$$
(2)
where $`\overline{\theta }_{x,\mu \nu }(\pi ,+\pi ]`$ and $`m_{x,\mu \nu }0`$ is called a Dirac plaquette. Monopole and photon fields are then defined in the following way

$$\theta _{x,\mu }^{\mathrm{mon}}=2\pi \underset{x^{}}{}G_{x,x^{}}_\nu ^{}m_{x^{},\nu \mu }$$
(3)
$$\theta _{x,\mu }^{\mathrm{phot}}=\underset{x^{}}{}G_{x,x^{}}_\nu ^{}\overline{\theta }_{x^{},\nu \mu }.$$
(4)
Here $`^{}`$ acts on $`x^{}`$, the quantities $`m_{x,\mu \nu }`$ and $`\overline{\theta }_{x,\mu \nu }`$ are defined in Eq. (2) and $`G_{x,x^{}}`$ is the lattice Coulomb propagator. One can show that $`\stackrel{~}{\theta }_{x,\mu }\theta _{x,\mu }^{\mathrm{mon}}+\theta _{x,\mu }^{\mathrm{phot}}`$ is up to a gauge transformation identical with the original $`\theta _{x,\mu }`$ defined by $`U_{x,\mu }=\mathrm{exp}(i\theta _{x,\mu })`$.

## 3 Classical chaotic dynamics from quantum Monte Carlo initial states

Chaotic dynamics in general is characterized by the spectrum of Lyapunov exponents. These exponents, if they are positive, reflect an exponential divergence of initially adjacent configurations. In case of symmetries inherent in the Hamiltonian of the system there are corresponding zero values of these exponents. Finally negative exponents belong to irrelevant directions in the phase space: perturbation components in these directions die out exponentially. Pure gauge fields on the lattice show a characteristic Lyapunov spectrum consisting of one third of each kind of exponents . Assuming this general structure of the Lyapunov spectrum we investigate presently its magnitude only, namely the maximal value of the Lyapunov exponent, $`L_{\mathrm{max}}`$.

The general definition of the Lyapunov exponent is based on a distance measure $`d(t)`$ in phase space,

$$L:=\underset{t\mathrm{}}{lim}\underset{d(0)0}{lim}\frac{1}{t}\mathrm{ln}\frac{d(t)}{d(0)}.$$
(5)
In case of conservative dynamics the sum of all Lyapunov exponents is zero according to Liouville’s theorem, $`L_i=0`$. We utilize the gauge invariant distance measure consisting of the local differences of energy densities between two $`3d`$ field configurations on the lattice:

$$d:=\frac{1}{N_P}_P\left|\mathrm{tr}U_P\mathrm{tr}U_P^{}\right|.$$
(6)
Here the symbol $`_P`$ stands for the sum over all $`N_P`$ plaquettes, so this distance is bound in the interval $`(0,2N)`$ for the group SU(N). $`U_P`$ and $`U_P^{}`$ are the familiar plaquette variables, constructed from the basic link variables $`U_{x,i}`$,

$$U_{x,i}=\mathrm{exp}\left(aA_{x,i}^cT^c\right),$$
(7)
located on lattice links pointing from the position $`x=(x_1,x_2,x_3)`$ to $`x+ae_i`$. The generators of the group are $`T^c=ig\tau ^c/2`$ with $`\tau ^c`$ being the Pauli matrices in case of SU(2) and $`A_{x,i}^c`$ is the vector potential. The elementary plaquette variable is constructed for a plaquette with a corner at $`x`$ and lying in the $`ij`$-plane as $`U_{x,ij}=U_{x,i}U_{x+i,j}U_{x+j,i}^{}U_{x,j}^{}`$. It is related to the magnetic field strength $`B_{x,k}^c`$:

$$U_{x,ij}=\mathrm{exp}\left(\epsilon _{ijk}aB_{x,k}^cT^c\right).$$
(8)
The electric field strength $`E_{x,i}^c`$ is related to the canonically conjugate momentum $`P_{x,i}=\dot{U}_{x,i}`$ via

$$E_{x,i}^c=\frac{2a}{g^3}\mathrm{tr}\left(T^c\dot{U}_{x,i}U_{x,i}^{}\right).$$
(9)
The Hamiltonian of the lattice gauge field system can be casted into the form

$$H=\left[\frac{1}{2}P,P+\mathrm{\hspace{0.17em}1}\frac{1}{4}U,V\right].$$
(10)
Here the scalar product stands for $`A,B=\frac{1}{2}\mathrm{tr}(AB^{})`$. The staple variable $`V`$ is a sum of triple products of elementary link variables closing a plaquette with the chosen link $`U`$. This way the Hamiltonian is formally written as a sum over link contributions and $`V`$ plays the role of the classical force acting on the link variable $`U`$.

Initial conditions chosen randomly with a given average magnetic energy per plaquette have been investigated in past years . In the present study we prepare the initial field configurations from a standard four dimensional Euclidean Monte Carlo program on a $`12^3\times 4`$ lattice varying the inverse gauge coupling $`\beta `$ . We relate such four dimensional Euclidean lattice field configurations to Minkowskian momenta and fields for the three dimensional Hamiltonian simulation by identifying a fixed time slice of the four dimensional lattice.

## 4 Chaos, confinement and continuum limit

We start the presentation of our results with a characteristic example of the time evolution of the distance between initially adjacent configurations. An initial state prepared by a standard four dimensional Monte Carlo simulation is evolved according to the classical Hamiltonian dynamics in real time. Afterwards this initial state is rotated locally by group elements which are chosen randomly near to the unity. The time evolution of this slightly rotated configuration is then pursued and finally the distance between these two evolutions is calculated at the corresponding times. A typical exponential rise of this distance followed by a saturation can be inspected in Fig. 2 from an example of U(1) gauge theory in the confinement phase and in the Coulomb phase. While the saturation is an artifact of the compact distance measure of the lattice, the exponential rise (the linear rise of the logarithm) can be used for the determination of the leading Lyapunov exponent. The left plot exhibits that in the confinement phase the original field and its monopole part have similar Lyapunov exponents whereas the photon part has a smaller $`L_{max}`$. The right plot in the Coulomb phase suggests that all slopes and consequently the Lyapunov exponents of all fields decrease substantially.

We now turn to a comparison with the monopole density in U(1) quantum field configurations,

$$Q_{mon}=\frac{1}{4V_4}\underset{x,\mu \nu }{}|m_{x,\mu \nu }|.$$
(11)
The left plot of Fig. 3 exhibits for a statistics of 100 independent configurations that $`Q_{mon}`$ decreases sharply between the strong and the weak coupling regime. It can be seen that the photon part from the decomposition carries also a few monopoles.

The main result of the present study is the dependence of the leading Lyapunov exponent $`L_{\mathrm{max}}`$ on the inverse coupling strength $`\beta `$, displayed in the right plot of Fig. 3. As expected the strong coupling phase, where confinement of static sources has been established many years ago by proving the area law behavior for large Wilson loops, is more chaotic. The transition reflects the critical coupling to the Coulomb phase. The right plot of Fig. 3 shows that the monopole fields carry Lyapunov exponents of nearly the same size as the full U(1) fields. The photon fields yield a non-vanishing value in the confinement ascending toward $`\beta =0`$ for randomized fields which indicates that the decomposition (1) works perfectly for ideal configurations only.

An interesting result concerning the continuum limit can be viewed from Fig. 4 which shows the energy dependence of the Lyapunov exponents for the U(1) theory and its components. One observes an approximately linear relation for the monopole part while a quadratic relation is suggested for the photon part in the weak coupling regime. From scaling arguments one expects a functional relationship between the Lyapunov exponent and the energy

$$L(a)a^{k1}E^k(a),$$
(12)
with the exponent $`k`$ being crucial for the continuum limit of the classical field theory. A value of $`k<1`$ leads to a divergent Lyapunov exponent, while $`k>1`$ yields a vanishing $`L`$ in the continuum. The case $`k=1`$ is special leading to a finite non-zero Lyapunov exponent. Our analysis of the scaling relation (12) gives evidence, that the classical compact U(1) lattice gauge theory and especially the photon field have $`k2`$ and with $`L(a)0`$ a regular continuum theory. The monopole field signals $`k1`$ and stays chaotic approaching the continuum.

## 5 Conclusions

We investigated the classical chaotic dynamics of U(1) lattice gauge field configurations prepared by quantum Monte Carlo simulation. The fields were decomposed into a photon and monopole part. The maximal Lyapunov exponent shows a pronounced transition as a function of the coupling strength indicating that on a finite lattice configurations in the strong coupling phase are substantially more chaotic than in the weak coupling regime. The computations give evidence that the Lyapunov exponents in the original U(1) field and in its monopole part are very similar. The situation for the monopole density is analogous and serves as a consistency check of the decomposition. We conclude that classical chaos in field configurations and the existence of monopoles are intrinsically connected to the confinement of a theory.

We found evidence that the monopole fields stay chaotic in the continuum while the photon fields and the full U(1) fields possess a regular continuum theory. So far the monopoles are extracted from the Euclidean gauge fields. It will be interesting to compute their counterparts in Minkowski space and perform an analysis concerning monopole annihilation.

Acknowledgments.  This work has been supported by the Austrian National Scientific Fund under the project FWF P14435-TPH. We thank Bernd A. Berg and Urs M. Heller as well as Tamás S. Biró and Natascha Hörmann for previous cooperation concerning topological objects and classical chaos in U(1) theory, respectively.


text:
# Pomeranchuk and other instabilities in the 𝑡-𝑡' Hubbard model at the Van Hove filling

(May 10, 2002)

## Abstract

We present a stability analysis of the two-dimensional $`tt^{}`$ Hubbard model for various values of the next-nearest-neighbor hopping $`t^{}`$, and electron concentrations close to the Van Hove filling by means of the flow equation method. For $`t^{}t/3`$ a $`d_{x^2y^2}`$-wave Pomeranchuk instability dominates (apart from antiferromagnetism at small $`t^{}`$). At $`t^{}<t/3`$ the leading instabilities are a $`g`$-wave Pomeranchuk instability and $`p`$-wave particle-hole instability in the triplet channel at temperatures $`T<0.15t`$, and an $`s^{}`$-magnetic phase for $`T>0.15t`$; upon increasing the electron concentration the triplet analog of the flux phase occurs at low temperatures. Other weaker instabilities are found also.

In recent years the two-dimensional (2D) Hubbard model has been used iz ; sc  as the simplest model which maps the electron correlations in the copper-oxide planes of high-temperature superconductors since experimental data suggest that superconductivity in cuprates basically originates from the CuO2 layers tk . Although in the high-temperature cuprate superconductors electron-electron interactions are strong some important features of these systems (in particular, antiferromagnetic and $`d`$-wave superconducting instabilities) are captured already by the 2D Hubbard model at weak to moderate Coulomb coupling.

Apart from the antiferromagnetism and $`d_{x^2y^2}`$-wave superconductivity mentioned above (for review see iz ; sc ; ha  and references therein), a few other instabilities related to symmetry-broken states hm ; hsg ; foh ; ikk ; hs ; hsr ; gkw  in the 2D $`tt^{}`$ Hubbard model with next-nearest-neighbor hopping $`t^{}`$ have been reported recently. Specially, much interest of researchers has been attracted by the case when the Fermi surface passes through the saddle points of the single particle dispersion (Van Hove filling). One of the instabilities found in such a case is a $`d`$-wave Pomeranchuk instability breaking the tetragonal symmetry of the Fermi surface, i.e. a spontaneous deformation of the Fermi surface reducing its symmetry to orthorhombic. It has been recently observed for small values of $`t^{}`$ from renormalization group calculations by Halboth and Metzner hm . They argued that the Pomeranchuk instability occurs more easily if the Fermi surface is close to the saddle points with a sizable $`t^{}`$ (reducing nesting which leads to antiferromagnetism). However, within their technique it is difficult to compare the strength of the Fermi surface deformation with other instabilities and to conclude which one dominates. The authors of Ref. hsr  have investigated the interplay of $`d`$-density wave na ; clmn  and Fermi surface deformation tendencies with those towards $`d`$-wave pairing and antiferromagnetism by means of a similar temperature-flow renormalization group approach. They have found that the $`d`$-wave Pomeranchuk instability never dominates in the 2D $`tt^{}`$ Hubbard model (even under the conditions mentioned above).

On the other hand, Vollhardt et al. vbh  showed that the $`t^{}`$-hopping term destroys the antiferromagnetic nesting instability at weak interactions in two and three dimensions, and supports the stabilization of metallic ferromagnetism in infinite dimensions away from half-filling. Therefore, one could expect also the stabilization of ferromagnetism by a sizable $`t^{}`$ in two dimensions. Indeed, in the $`tt^{}`$ Hubbard model on a 2D square lattice at week to moderate Coulomb coupling, a projection quantum Monte Carlo calculation with $`20\times 20`$ sites and the $`T`$-matrix technique hsg , a generalized random phase approximation including particle-particle scattering foh  point towards a ferromagnetic ground state for large negative values of $`t^{}/t`$ in a density range around the Van Hove filling. Similar tendencies have been found by the authors of Ref. ikk  within the renormalization group and parquet approaches. Honerkamp and Salmhofer recently studied hs  the stability of this ferromagnetic region at finite temperatures by means of the temperature-flow renormalization group technique. They have found that ferromagnetic instability is the leading one at $`t^{}<0.33t`$ and Van Hove filling with critical temperatures depending on the value of $`t^{}`$. When the electron concentration is increased slightly above the Van Hove filling, the ferromagnetic tendencies get cut off at low temperatures and a triplet $`p`$-wave superconducting phase dominates. However, they did not consider the Pomeranchuk instability (which could have the most favorable conditions to occur) and other ones apart from antiferromagnetism, $`d`$\- and $`p`$-wave superconductivity and ferromagnetism.

Therefore, the investigation of interplay and rivalry between the Pomeranchuk instability and ferromagnetism, and other phases in the 2D $`tt^{}`$ Hubbard model at the Van Hove filling is a considerable task. We will consider the leading instabilities depending on the ratio $`U/t`$ (in all papers cited above it was fixed). The main goal of this paper is such a study. We report also a few new instabilities in a range of electron concentration around the Van Hove filling.

We start from the Hamiltonian of the $`tt^{}`$ Hubbard model

$`H={\displaystyle \underset{𝐤\sigma }{}}\epsilon _𝐤c_{𝐤\sigma }^{}c_{𝐤\sigma }+{\displaystyle \frac{U}{N}}{\displaystyle \underset{\genfrac{}{}{0pt}{}{𝐤_\mathrm{𝟏}𝐤_\mathrm{𝟏}^{}}{𝐤_\mathrm{𝟐}𝐤_\mathrm{𝟐}^{}}}{}}c_{𝐤_\mathrm{𝟏}}^{}c_{𝐤_\mathrm{𝟏}^{}}c_{𝐤_\mathrm{𝟐}}^{}c_{𝐤_\mathrm{𝟐}^{}}\delta _{𝐤_\mathrm{𝟏}+𝐤_\mathrm{𝟐},𝐤_\mathrm{𝟏}^{}+𝐤_\mathrm{𝟐}^{}},`$ (1)
where $`\epsilon _𝐤`$ is the Bloch electron energy with the momentum $`𝐤`$, $`c_{𝐤\sigma }^{}(c_{𝐤\sigma })`$ is the creation (annihilation) operator for the electrons with spin projection $`\sigma \{,\}`$, $`U`$ is the local Coulomb repulsion of two electrons of opposite spins, $`N`$ is the number of lattice points, lattice spacing equals unity.

By means of the flow equation method we  the Hamiltonian is transformed into one of molecular-field type. This Hamiltonian is calculated in second order in the coupling $`U`$ gkw . Adopting the notations of Ref. gkw , the expression for the free energy has the form:

$`\beta F={\displaystyle \frac{1}{N}}{\displaystyle \underset{\mathrm{𝐤𝐪}}{}}\beta U\left(1+{\displaystyle \frac{U}{t}}V_{𝐤,𝐪}\right)\mathrm{\Delta }_𝐤^{}\mathrm{\Delta }_𝐪+{\displaystyle \underset{𝐤}{}}f_𝐤\mathrm{\Delta }_𝐤^{}\mathrm{\Delta }_𝐤,`$ (2)
where the first term is the energy contribution and the second term is the entropy contribution, $`\beta =1/(k_BT)`$, $`T`$ is the temperature, $`t`$ is the hopping integral of electrons between nearest neighbors of the lattice, $`V_{𝐤,𝐪}`$ is effective second-order interaction (the factor $`U^2/t`$ has been extracted from it), $`f_𝐤`$ is an entropy coefficient, and $`\mathrm{\Delta }_𝐤`$ are the order parameters. For example, $`\mathrm{\Delta }_{𝐤\sigma ,𝐤\sigma ^{}}=c_{𝐤\sigma }c_{𝐤\sigma ^{}}=(\sigma _y)_{\sigma \sigma ^{}}\mathrm{\Delta }_𝐤^s+_\alpha (\sigma _y\sigma _\alpha )_{\sigma \sigma ^{}}\mathrm{\Delta }_𝐤^{t\alpha }`$, where $`\sigma _\alpha `$ is a Pauli spin matrix ($`\alpha =x,y,z`$), and $`\mathrm{\Delta }_𝐤^s(\mathrm{\Delta }_𝐤^{t\alpha })`$ is the singlet (triplet) amplitude. An expression similar to Eq. (2) is obtained for particle-hole channels with the order parameters $`\nu `$ instead of $`\mathrm{\Delta }`$. In this case, for example, we have $`\nu _{𝐤\sigma ,𝐤\sigma ^{}}=c_{𝐤\sigma }^{}c_{𝐤+𝐐\sigma ^{}}=\nu _𝐤^s\delta _{\sigma ,\sigma ^{}}+_\alpha \nu _𝐤^{t\alpha }(\sigma _\alpha )_{\sigma \sigma ^{}}`$ with $`𝐐=(\pi ,\pi )`$. All quantities of Eq. (2) are defined in Ref. gkw .

For a square lattice the single particle dispersion has the form:

$`\epsilon _𝐤=2t(\mathrm{cos}k_x+\mathrm{cos}k_y)4t^{}\mathrm{cos}k_x\mathrm{cos}k_y.`$ (3)
The spectrum (3) contains Van Hove singularities in the density of states at the energy $`\epsilon _{VH}=4t^{}`$ related to the saddle points of the Fermi surface at $`𝐤=(0,\pm \pi )`$ and $`(\pm \pi ,0)`$. For $`t^{}=0`$ and half-filling the Fermi surface is nested $`\epsilon _{𝐤+𝐐}=\epsilon _𝐤`$, which leads to an antiferromagnetic instability for $`U>0`$. The nesting is removed for $`t^{}/t0`$.

We start from the symmetric state and investigate whether this state is stable against fluctuations of the order parameters $`\mathrm{\Delta }`$ and $`\nu `$. As soon as a non-zero $`\mathrm{\Delta }`$ or $`\nu `$ yields a lower free energy in comparison with the symmetric state with all vanishing $`\mathrm{\Delta }`$ and $`\nu `$, then the symmetric state is unstable and the system will approach a symmetry broken state. This indicates a phase transition.

We perform numerical calculation on a square lattice with $`24\times 24`$ points in the Brillouin zone for the various representations under the point group $`C_{4\nu }`$. The representations of the even-parity states are one-dimensional. We denote them by $`s_+=s_1,s_{}=s_{xy(x^2y^2)},d_+=d_{x^2y^2},d_{}=d_{xy}`$. The odd-parity representation is two-dimensional, here simply denoted by $`p`$. Initially, such numerical calculations have been performed in Refs. gkw ; gr , but they were sensitive to the lattice size at low temperatures. Here we use an improved scheme (for details see Ref. hw ).

We start from $`t^{}=0`$ and half-filling ($`n=1`$) (see Fig. 1). As expected in this case the leading instability is the antiferromagnetic one which disappears at the temperature $`T0.1t`$ or doping $`\delta n1=0.06`$. The next instability is a Pomeranchuk instability with $`d_{x^2y^2}`$-wave symmetry in the singlet channel. The corresponding eigenvectors signals a deformation of the Fermi surface which breaks the point group symmetry of the square lattice. For negative $`t^{}t/3`$ the Pomeranchuk instability dominates at the Van Hove filling (see Fig. 2). The $`d_{x^2y^2}`$\- wave Pomeranchuk instability competes with other instabilities at $`t^{}<t/3`$, and it is not the leading one (Fig. 3). In agreement with the ideas of Ref. hm  the instability is mainly driven by a strong attractive interaction between particles on opposite corners of the Fermi surface near the saddle points and a repulsive interaction between particles on neighboring corners. To favor such a behavior we need a sizable $`t^{}`$ reducing antiferromagnetic correlations.

At half-filling and $`t^{}=0`$ the next instability is a particle-hole instability of singlet type with staggered $`p`$-wave symmetry. It yields gkw  a splitting into two bands and may lead to an energy gap in the charge excitations spectrum. Another mechanism for a charge gap formation has been proposed mj ; klpt  recently in the 2D Hubbard model with $`t^{}=0`$ at weak coupling. The band splitting phase is developed in the region of electron concentration around half-filling, and is one of the strongest in that region. Then the superconducting $`d_{x^2y^2}`$ instability follows which coincides with the $`d_{x^2y^2}`$-wave staggered flux phase (the flux phase has been proposed by the authors of Ref. am ; ko  and discussed recently in Refs. na ; clmn ; le ). Away from half-filling the degeneration disappears, and $`d`$-wave superconductivity dominates at low temperatures in certain regions of electron concentration around half-filling which depends on the value of $`t^{}0`$. Even large values of $`|t^{}|`$ do not destroy the dominant low-temperature behavior of $`d_{x^2y^2}`$-wave superconductivity at doping hw . One phase may suppress another phase. To which extend two order parameters can coexist with each other is a question, which has to be investigated in the future.

For $`t^{}=0`$ the singlet and triplet $`T_c`$ of the particle-hole instabilities with staggered symmetry of $`d_+`$ wave character (that is the flux phase) are degenerate. If $`t^{}0`$ they are different, and the triplet one is higher. Moreover, the triplet analog of flux phase dominates at low temperatures and $`t^{}=5t/12`$ when the electron concentration is slightly above the Van Hove filling (see Fig. 4) in contrast to the results of Ref. hs  which point out the occurrence of triplet superconductivity with $`p`$-wave symmetry in this region. The triplet flux phase is also one of the leading instabilities for $`t^{}t/3`$ and certain region of electron concentrations (see Fig. 2). It has been considered by Nayak na  as a density wave order parameter potentially relevant to the cuprates, but to our knowledge a triplet version of the flux phase has not yet been observed in numerical solutions of the 2D $`tt^{}`$ Hubbard model. We shall discuss this state in more details elsewhere hw .

At $`t^{}=5t/12`$ a few other new instabilities appear to compete at the Van Hove filling and low temperatures (Fig. 3) in disagreement with the conclusions of Ref. hs  on the occurrence of ferromagnetism. The leading one is a Pomeranchuk instability in the $`s_+`$ channel with $`g_+=g_{x^4+y^46x^2y^2}`$ wave character (4 node lines in $`𝐤`$-space). This phase occurs more easily if the electron concentration is close to or slightly smaller than the Van Hove filling (Fig. 5). It requests also sufficiently large absolute values of $`t^{}`$. When the electron concentration is decreased below the Van Hove density, a particle-hole instability of $`p`$-wave symmetry in triplet channel dominates at low temperatures (see Fig. 5), which gives rise to a phase of magnetic currents. In the $`d_{}`$ channel an $`i`$-wave (6 node lines in $`𝐤`$-space) Pomeranchuk instability appears when electron concentration $`n`$ is smaller or close to the Van Hove filling (Figs. 3-5). It is a leading one at small values of the electron concentration hw . We observe (Fig. 5) in the $`s_+`$ channel a $`g_+`$ wave superconductivity below the Van Hove filling, but it requires strong coupling.

Another situation occurs at the temperature region $`T>0.15t`$. Here a particle-hole instability with $`s^{}`$-wave character (its order parameter changes sign close to the Fermi-edge) in the triplet channel dominates at the density range around the Van Hove filling (see Figs. 3-5). It is likely that the order parameter contributions do not compensate exactly, so that a weak ferromagnetism appears. When the electron concentration is increased above the Van Hove filling this instability does not become weaker, but the $`d_{x^2y^2}`$ wave Pomeranchuk and the triplet flux phase instabilities are manifested stronger (they dominate at low temperatures). Then, this $`s^{}`$-magnetic phase disappears at sufficiently large values of electron concentration in comparison with the Van Hove filling, or smaller $`|t^{}|`$.

From Figs. 3-5 one can see a reentrant behavior of the $`s^{}`$-magnetic phase in some region of the values $`U/t`$: approaching $`T_c`$ we get $`T_c^l`$ from low temperatures and $`T_c^u`$ from high temperatures at the same value of coupling $`U/t`$ ($`T_c^lT_c^u`$). This is a result of different behavior of $`T_c(U/t)`$ in two regimes. First regime occurs in the situation where the $`s^{}`$-magnetic instability dominates and the transition from a paramagnetic state to the $`s^{}`$-magnetic phase occurs directly without any intermediate phase, it corresponds to the temperatures $`T>0.15t`$ on Figs. 35. In this case the critical temperature increases with the increase of correlation strength $`U/t`$. Another regime occurs at the temperatures $`T<0.15t`$. In this situation the critical temperature exhibits an anomalous behavior, it decreases with increasing the coupling $`U/t`$. The $`s^{}`$-magnetic phase is reduced. Since at lower temperatures only a smaller region in $`𝐤`$-space around the Fermi-edge contributes, the sign-change of the order parameter reduces the effective interaction.

In conclusion, we have presented a stability analysis of the 2D $`tt^{}`$ Hubbard model on a square lattice for various values of the next-nearest-neighbor hopping $`t^{}`$ and electron concentrations close to the Van Hove filling. A surprising large number of phases has been observed. Some of them have an order parameter with many nodes in $`𝐤`$-space. For $`t^{}t/3`$ the $`d_{x^2y^2}`$-wave Pomeranchuk instability dominates. At $`t^{}<t/3`$ the leading instabilities are a $`g_+`$ wave Pomeranchuk instability and $`p`$-wave particle-hole instability in triplet channel at temperatures $`T<0.15t`$, and $`s^{}`$-magnetic phase for $`T>0.15t`$; upon increasing the electron concentration the triplet flux phase occurs at low temperatures. The $`s^{}`$-magnetic phase is reduced strongly at low temperatures. We have found other weaker instabilities also.