text:
# Elliptic flow orientation, saturation and low 𝑝_𝑡 behavior

## 1 Introduction

Azimuthal momentum space distributions is very interesting effect widely studying at RHIC. The strong elliptic flow generated in non-central collisions was founded in the set of experiments. Basis method is study Fourier expansion of azimuthal particle distributions . For study elliptic flow we interesting $`v_2`$ which is calculated from azimuthal average $`\mathrm{cos}(2\varphi )`$, where angle $`\varphi `$ defined relatively to reaction plain which is must be determined for each event and defined as plain based on the beam axis and line joining center of the nucleous. Next, we must perform average over all set of events.

There are many phenomenological explanation of the elliptic flow phenomenon. The main goal of the theoretical studying is explain how transverse plain asymmetry of matter in the reaction zone is transformed into asymmetry of the momentum space distribution. Existence of the azimuthal asymmetry is striking evidence of reinteractions in the reaction zone which enforce matter into asymmetric radial expansion. The modern most desired explanation of low $`p_t`$ elliptic flow is existence of the Quark Gluon Plasma (QGP) and thermalization stage during nucleus collision. This allow us to use hydrodynamic with relativistic equation of motion $`_\mu T^{\mu \nu }`$. But it is necessary to know pressure as function of energy density and chemical potential. Pressure can be founded from QGP equation of state which unknown, especially near phase transition point. Many numeric calculations performed in the longitudinal boost invariant ideal fluid model with equation of state of ideal gas: $`PT^4`$. In additional dissipative effects, non-equilibrium transport equations and other improvements can be used with hydrodynamic model. Unfortunaly, we can not reach full physical understanding because too much uncontrolled phenomenological assumptions but so small experimental observable variables.

Recent experimental data from RHIC for Au+Au collisions for identified particles  show saturation of $`v_2`$ in the region $`p_t>2GeV`$, Fig. 1. Hydrodynamic calculation can not successful reproduce $`v_2`$ behavior in the saturation region. The origin saturation is still unknown. At high $`p_t`$ was proposed a mechanism of jet quenching – high $`p_t`$ pertrubative gluons propagated through asymmetric reaction zone and lost its momentum interacting with dense matter . Also proposer the contribution to $`v_2`$ from nucleus classical fields in the saturation regime . In our calculation we show that the saturation phenomenon have mostly kinematic origin, but saturation values depend on details of matter motion due to reinteractions.

In this paper we try to separate contributions to $`v_2`$ from macroscopic flow of matter $`u(\stackrel{}{\beta })`$ and from effect of the local matter decay. And try to describe low $`p_t`$ and high $`p_t`$ of $`v_2`$ behavior in unified language. In our model for final particles spectrum we assume factorization between macroscopic flow and matter decay distribution $`f`$. In general case, $`u(\stackrel{}{\beta })`$ can be hydrodynamic flow and $`f`$ be a thermal distribution of specified hardron at the freeze-out temperature. So we can say that we add local decay matter effects to the hydrodynamic model. On the other hand, $`u(\stackrel{}{\beta })`$ can be color string transverse motion due its repulsion and $`f`$ be standard color string fusion spectrum. In same sense we generalize Blast Wave model  and we not concerning on real sense of flow $`u(\stackrel{}{\beta })`$ and decay distribution $`f`$. The model does not cover jet quenching because pertrubative jets produced at the early stage of collision where macroscopic flow is still not formed.

## 2 Elliptic flow at low and high $`p_t`$

Consider elliptic flow in c.m. system at at midrapidity, $`y=0`$. Assuming boost invariance as usual we may reduce problem to the two dimension transverse plain calculations. Let $`u(\stackrel{}{\beta })=u(\beta ,\psi )`$ be a density of expanding matter depend on transverse flow velocity $`\stackrel{}{\beta }`$. We always assume that exists some $`\beta _{max}<1`$ and $`u(\stackrel{}{\beta })=0`$ if $`|\stackrel{}{\beta }|>\beta _{max}`$. For final transverse particle distribution we can write

$$F(\stackrel{}{p})=F(p,\varphi )=\frac{1}{Z(p)}\frac{d^2N}{dp^2}=\frac{1}{Z(p)}u(\stackrel{}{\beta })f_\beta (\stackrel{}{p})d^2\beta $$
(1)
where $`f(\stackrel{}{p})=f(|p|^2)`$ – transverse fusion spectrum of small piece of expanding matter, $`f_\beta (\stackrel{}{p})`$ – modified fusion spectrum of fixed particle type according to it transverse motion with $`\stackrel{}{\beta }`$ velocity (one can speak about ”blue shift”). With aim to exclude final spectrum mixing between particles of different species we consider only elliptic flow for identified particles. Factor $`Z(p)`$ defined to satisfy azimuthal normalization

$$\underset{0}{\overset{2\pi }{}}F(p,\varphi )𝑑\varphi =1$$
(2)
$$Z(p)=\left(u(\stackrel{}{\beta })f_\beta (\stackrel{}{p})d^2\beta \right)𝑑\varphi $$
(3)
We consider only collisions of equal nucleus such as Au+Au. Then the distribution $`u(\stackrel{}{\beta })`$ have very important symmetry properties which follow from the basic properties of the nucleus collision geometry.

$$u(\beta ,\psi )=u(\beta ,\psi +\pi )=u(\beta ,\psi )$$
(4)
This automatic gives

$`{\displaystyle _0^{2\pi }}u(\beta ,\psi )\mathrm{sin}(n\psi )𝑑\psi =0`$ (5)
$`{\displaystyle _0^{2\pi }}u(\beta ,\psi )\mathrm{cos}((2n+1)\psi )𝑑\psi =0`$ (6)
So for $`u(\stackrel{}{\beta })`$ we have only even $`\mathrm{cos}`$ harmonics in the azimuthal Fourier decomposition.

Elliptic flow can be calculated from Fourier decomposition of $`F(p,\varphi )`$

$$v_2=\underset{0}{\overset{2\pi }{}}F(p,\varphi )\mathrm{cos}(2\varphi )𝑑\varphi =\mathrm{cos}(2\varphi )$$
(7)
$$F(p,\varphi )=\frac{1}{2\pi }\left(1+2v_2\mathrm{cos}(2\varphi )+\mathrm{}\right)$$
(8)
Now our task is to calculate $`v_2`$ from (7) using (1).

The integrated $`v_4`$ is about a factor of 10 smaller than $`v_2`$ . In cases when $`v_1`$(directed flow) $`v_4`$ and other higher harmonics enough small relatively to $`v_2`$ (one can say that the azimuthal distribution have shape like $`\mathrm{cos}2\varphi `$) we can simple estimate $`v_2`$ with sufficient accuracy

$$v_2=\frac{\pi }{2}\left(F(p,0)F(p,\frac{\pi }{2})\right)$$
(9)
For find $`f_\beta (\stackrel{}{p})`$ in relativistic form we must boost lab $`\stackrel{}{p}`$ to the comoving frame. Decompose $`\stackrel{}{p}`$ into $`\stackrel{}{\beta }`$ and $`\stackrel{}{\beta }`$ orthogonal directions

$$\stackrel{}{p}=p\mathrm{cos}(\varphi \psi )\frac{\stackrel{}{\beta }}{|\beta |}+p\mathrm{sin}(\varphi \psi )\frac{\stackrel{}{\beta }_{ort}}{|\beta |}$$
(10)
Now we can easy perform boost in $`\stackrel{}{\beta }`$ direction. In relativistic form for $`f_\beta (\stackrel{}{p})`$ we have

$$f_\beta (\stackrel{}{p})=f\left(p^2\mathrm{sin}^2(\varphi \psi )+\frac{\left(p\mathrm{cos}(\varphi \psi )\beta \sqrt{M^2+p^2}\right)^2}{1\beta ^2}\right)$$
(11)
We also need to relativistic transform of phase space $`d^{\mathrm{\hspace{0.25em}2}}p`$ in (1)

$$\frac{d^{\mathrm{\hspace{0.25em}2}}p}{\sqrt{M^2+p^2}}=inv$$
(12)
Performing calculation we obtain

$$F(p,\varphi )=\frac{1}{Z(p)}u(\stackrel{}{\beta })f_\beta (\stackrel{}{p})\left(1\frac{\beta p\mathrm{cos}(\varphi \psi )}{\sqrt{M^2+p^2}}\right)\frac{\beta d\beta d\psi }{\sqrt{1\beta ^2}}$$
(13)
where $`M`$ – identified particle mass and distribution $`f`$ must be considered as function on $`p^2`$.

Consider $`v_2`$ behavior at small $`p`$. In this paper for simplicity we show calculation with $`f_\beta (\stackrel{}{p})`$ in non-relativistic form

$$f_\beta (\stackrel{}{p})=f(|\stackrel{}{p}M\stackrel{}{\beta }|^2)=f(p^2+M^2\beta ^22M\beta p\mathrm{cos}(\varphi \psi ))$$
(14)
Relativistic calculation based on (11) also can be performed. The resulting $`p^2`$ behavior will not changed. Only constant factor definition will changed.

Due to existence factor $`Z^1(p)`$ in (1), any $`\varphi `$-constant multiplier in $`u(\stackrel{}{\beta })f_\beta (\stackrel{}{p})d^2\beta `$ does not contribute to $`F(p,\varphi )`$ and can be dropped. Distribution $`f_\beta (p)`$ at small $`p`$ can be expanded in Taylor series. Zero expansion order does not contribute to $`v_2`$ because there is no dependence on $`\varphi `$. In zero order $`Z(p)`$ does not vanish

$$Z(p)=Z_0=u(\beta ,\psi )f(M^2\beta ^2)\beta 𝑑\beta 𝑑\psi $$
(15)
In the first order we have $`u(\psi )\mathrm{cos}(\varphi \psi )𝑑\psi =0`$ due to (5) and (6). Second order have contribution to $`v_2`$

$$F(p,\varphi )=\frac{1}{Z(p)}u(\beta ,\psi )\frac{1}{2}f^{\prime \prime }(M^2\beta ^2)(2\beta Mp)^2\mathrm{cos}^2(\varphi \psi )\beta 𝑑\beta 𝑑\psi $$
(16)
where we omit all terms which give zero contribution to $`v_2`$.

If we need only $`v_2`$ then

$$\mathrm{cos}^2(\varphi \psi )=\mathrm{cos}^2\varphi \mathrm{cos}^2\psi +\mathrm{sin}^2\varphi \mathrm{sin}^2\psi =\frac{1}{2}\mathrm{cos}2\varphi \mathrm{cos}2\psi $$
(17)
Performing all calculations we find

$$v_2=C(\{u\},\{f\},M)p^2$$
(18)
where

$$C=\frac{2\pi M^2u(\beta ,\psi )\mathrm{cos}2\psi f^{\prime \prime }(M^2\beta ^2)\beta ^3𝑑\beta 𝑑\psi }{u(\beta ,\psi )f(M^2\beta ^2)\beta 𝑑\beta 𝑑\psi }$$
(19)
We found that the $`v_2`$ grow as $`p^2`$ at small p and vanish at $`p=0`$. This fact is not depend on detail form of $`f(\stackrel{}{p})`$ distribution. In the Fig. 1 we can see $`p^2`$-grow of $`v_2`$ for identified particles at small $`p`$.

Consider sign of $`v_2`$. In most cases $`u(\beta ,\psi )\mathrm{cos}2\psi d\psi >0`$ because we set zero azimuthal angle to flow maximum. For any known distribution $`f^{\prime \prime }(p^2)>0`$ for $`p>p_0`$, where $`p_0`$ – same specific sufficient small momentum. So we can conclude that in our model $`v_2>0`$ as observed111In the low energy collision mean $`\beta `$ can be small. So $`f^{\prime \prime }(M^2\beta ^2)`$ and $`v_2`$ may be negative..

Now consider case when $`p`$ is large. For reason of explicit calculation we assume that the $`f(p)`$ at large $`p`$ have power like tail and have asymptotic $`|p|^{2\alpha }`$. This is consistent with any known fits of particle distributions. At first we consider non-relativistic case.

$$F=\frac{1}{Z(p)}u(\stackrel{}{\beta })p^{2\alpha }\left(1+\left(\frac{M\beta }{p}\right)^22\left(\frac{M\beta }{p}\right)\mathrm{cos}(\varphi \psi )\right)^\alpha d^2\beta $$
(20)
We can use $`M\beta /p`$ as small value and expand (20) in Taylor series.

$$F=\frac{1}{Z(p)}u(\stackrel{}{\beta })p^{2\alpha }\left(1+2\alpha ^2\left(\frac{M\beta }{p}\right)^2\mathrm{cos}^2(\varphi \psi )\right)d^2\beta $$
(21)
where we omit first expansion order because it does not contribute to $`v_2`$.

In leading order we can calculate $`Z(p)`$

$$Z(p)=2\pi p^{2\alpha }u(\stackrel{}{\beta })d^2\beta $$
(22)
For $`\mathrm{cos}^2(\varphi \psi )`$ in (21) we can apply rule (17).Then we have

$$v_2=\left(\frac{M\beta \alpha }{p}\right)^2\frac{u(\stackrel{}{\beta })\mathrm{cos}2\psi d^2\beta }{2u(\stackrel{}{\beta })d^2\beta }$$
(23)
So we see $`1/p^2`$ decreasing of $`v_2`$ in the non-relativistic model.

In the relativistic model we must determine asymptotic properties of (11) where $`f|p|^{2\alpha }`$.

$$f_\beta (\stackrel{}{p})=p^{2\alpha }\left|\frac{\beta ^2}{1\beta ^2}\mathrm{cos}^2(\varphi \psi )\frac{2\beta }{1\beta ^2}\mathrm{cos}(\varphi \psi )\right|^\alpha $$
(24)
Inserting it into (1) we see that the $`|p|^{2\alpha }`$ factor absorbed into factor $`Z(p)`$ and does not influence on the $`v_2`$. So $`v_2`$ tend to some constant at high $`p_t`$

$$v_2=\frac{u(\stackrel{}{\beta })R(\mathrm{cos}(\varphi \psi ),\beta )cos(2\varphi )d^2\beta 𝑑\varphi }{u(\stackrel{}{\beta })R(\mathrm{cos}(\varphi \psi ),\beta )d^2\beta 𝑑\varphi }$$
(25)
where $`R(\mathrm{cos}(\varphi \psi ),\beta )=f_\beta (\stackrel{}{p})/p^{2\alpha }`$.

It is interesting to know real $`v_2`$ behavior at high $`p_t`$ to distinguish relativistic and non-relativistic model. Unfortunaly, the current experimental data have large error bars at $`p_t>4Gev`$ and high $`p_t`$ behavior in not clear.

## 3 Elliptic flow orientation and directed flow

Asymmetric nucleus collision is collision between to different nucleus at midrapidity or collision of equivalent nucleus at rapidity near $`\pm y_{beam}`$. If collision is asymmetric then (6) for $`n=0`$ does not holds. And we can expect existence of the directed flow $`v_10`$. Non-zero $`v_1`$ is usually linked with barion stopping, limited fragmentation and projectile nucleon bounce.

Important property of directed flow is orientation relatively to reaction plain. Directed flow can lie only in the reaction plain due to symmetry. Of course, we assume that there is no symmetry breaking in each nucleus collision. If $`v_1`$ directions was determined, the elliptic flow $`v_2`$ must be calculated relatively this direction. If founded $`v_2`$ positive than we can conclude that the elliptic flow is in reaction plain and maximum matter flow in direction where almond shape reaction zone is thin. And vice versa, if $`v_2`$ negative than elliptic flow in direction where reaction zone is thick.

Consider potential expansion model. In this model we assume that there are only potential forces acting on produced matter after collision. Each two near pieces of matter after its production have positive potential energy of reinteraction. Instead abstract ”matter” we can think about strong color fields or about color strings. The origin of repulsive forces ab initio is not known. In pure phenomenological language we can address problem to non-pertrubative QCD vacuum properties.

For simplicity we consider color strings placed in one dimensional line with equal spacing. Let point $`x=0`$ be a center of this string chain, $`\rho `$ – string tension, $`u`$ \- repulsive potential energy of neighboring string per unit length. After time evolution of this system all potential energy transferred into kinetic energy of transverse string motion. Due to symmetry we can consider only strings at $`x>0`$. Let index $`i`$ enumerate strings from $`x=0`$ to positive direction. String having $`i=0`$ after expansion will not moving due to symmetry. Transverse rapidity $`\theta _1`$ of $`i=1`$ string

$$\theta _1=arcch\frac{u}{\rho }1$$
(26)
Other strings rapidity is simple additive

$$\theta _i=i\theta _1$$
(27)
So longest string chain have larger transverse velocity of leftmost and rightmost strings and in thick direction of reaction zone we have flow maximum.

Determination of $`v_2`$ sign relatively to $`v_1`$ orientation is very important for testing of validity of different models describing reinteractions during collision. For example, in the hydrodynamic-inspired models $`v_2`$ positive due to maximum pressure gradient in thin direction. In the potential expansion model where role of matter play color strings we expect negative $`v_2`$ due to maximum energy released in thick direction. At different energies different mechanism can play a main role. First signal of $`v_2`$ positivity at RHIC energies was found in . But this is not full story. In our model we can see in (19) that $`v_2`$ sign can changed not only due preferred expansion direction but due to negative sign of $`f^{\prime \prime }`$.

Near midrapidity at small $`y`$ we can expand $`u(\stackrel{}{\beta },y)`$ in Taylor series

$$u(\stackrel{}{\beta },y)=u_0(\stackrel{}{\beta })+u_1(\stackrel{}{\beta })y+\mathrm{}$$
(28)
where $`u_0(\stackrel{}{\beta })`$ have zero directed flow and $`u_1(\stackrel{}{\beta })`$ – non-zero.

In our model directed flow have form

$$v_1(p)=\underset{0}{\overset{2\pi }{}}\frac{1}{Z(p)}\left(u(\stackrel{}{\beta })f_\beta (p)d^2\beta \right)\mathrm{cos}\varphi d\varphi $$
(29)
Consider $`v_1`$ at low $`y`$ and low $`p`$. Similarly (18) we get

$$v_1=Cyp$$
(30)
where $`C`$ – some constant not depended on $`y`$ and $`p`$. Note that the sign of $`v_1`$ have no sense due to impossibility of its measurement. Such $`y`$ scaling we was observed in . Dependence $`v_1(p)`$ is not measured yet.

## 4 Conclusion

In this paper we have studied azimuthal anisotropy at large and small $`p_t`$ in generalized model with factorization between macroscopic flow and matter decay distribution. We showed universality of $`p^2`$ behavior of elliptic flow at small $`p`$ for specified final particles – $`v_2`$ grow as $`p^2`$ at small p and vanish at $`p=0`$. This fact is not depend on detail form of $`f(\stackrel{}{p})`$ distribution and matter flow $`u(\stackrel{}{\beta })`$. So we can say that it is model independent.

At high $`p_t`$ we compared non-relativistic and relativistic models of Lorentz boosted distribution $`f`$. In the relativistic models with distribution $`f`$ having power-like tail we show existence of elliptic flow saturation. In the non-relativistic model $`v_2`$ decrease as $`1/p^2`$. More precise measurements of $`v_2`$ at region $`p_t>4Gev`$ are needed to make selection between models.

Determination the sign of $`v_2`$ relatively to directed flow is important criteria to for testing of validity of different models describing reinteractions during collision. To make conclusion about model selection we need experimental determination of $`v_2`$ sign in wide range of collision energy. For example, we can compare hydro and potential expansion models.

We also shown that at small $`p_t`$ sign of elliptic flow can change due to sign of second derivative of decay distribution $`f`$. This fact must be considered if we want to make model selection based on elliptic flow orientation. For example, for color string at sufficient small $`p_t`$ we have $`f^{\prime \prime }>0`$.

## Acknowledgments

We thank A.V. Dmitriev and N.V. Prikhod’ko for useful discussions. This work was supported by RFBR Grant RFBR-03-02-16157a and grant of Ministry for Education E02-3.1-282


format: text

text:
# Standard Model Masses and Models of Nuclei

## I Introduction

The measurement of the highest masses in the standard model, namely those of the top quark and the particles W and Z, shows an intriguing coincidence: the mass of the top occurs slightly before than the mass of the doubly magic nucleus of 208 Pb, while the mass of Z and W occur slightly before of the mass of 100 Sn.

These two magic numbers need of an abnormally high spin-orbit splitting, and there is only a third doubly magic number sharing this need: the one of 132 Sn. Then, the ALEPH preliminary measurement ellis  of a mass of 115 GeV for the Higgs boson, just slightly before the mass of this Tin isotope, makes imperative to report this triple coincidence and to state the problem, which I do in this letter. Due to the continuous action of strong nuclear forces, the weak bonding approximation is not valid in this context, so it is plausible that we are seeing some kind of general mesoscopic effect coming from many-body quantum mechanics, and not just a specific effect of nuclear physics.

It is known that magic numbers in nuclei are generated by separating levels via spin-orbit coupling. Low numbers, as the 28 or 20 in doubly magic Ca48 or Ni56, can be fit using traditional potentials with a relatively weak spin-orbit splitting, if any. Numbers 50, 82 and 126 are got because the respective subshells 5g9/2, 6h11/2 and 7i13/2 get large gaps and fall into the lower shell. Traditionally this was obtained via a purely phenomenological Hamiltonian, consisting of an empirical spin-orbit coupling and a Nilsson term $`+k_nj^2`$ depending of the shell. Modernly a relativistic approach based in meson exchange lets one to derive the spin orbit term from first principles. Still, the models are restricted to low nucleon numbers, because of computational effort. Some high (N,Z) models have been essayed, but they fail to reproduce the phenomena of double magicity and its associated quenching.

Here we see that some empirical coincidences are related to the availability of virtual J=0 and J=1 particles from the Standard Model. They seem to imply a need to correct meson exchange by adding the whole set of standard model interacting particles. Such particles should be noticeable during the process of filling the subshells responsable of the magicity, thus they occur slightly before of the double magic numbers.

Let me to organise the letter in small sections, so I will first comment on global issues, then to look into each level in detail, and finally to do some remarks about their possible significance.

## II mass landscape

There are four highly massive particles in the standard model: the vector bosons W,Z, the scalar Higgs boson, and the Top quark. This one can not appear free but in composites, the simplest being the bosons $`(t\overline{u}),(t\overline{d})`$, etc, all of them having a mass near to the Top mass. Thus we will refer to all these particles as ”the standard model bosons”.

If we draw the traditional plot of energy per nucleon, for instance for the Audi-Wapstra experimental tables, we notice three peaks111about or less than 1 percent in magnitude. in the descending part. Converting from GeV to atomic mass units, we can draw over the plot four lines indicating the respective masses of W, Z,Higgs (conjectured by LEP-2), and Top. Then we note that boson masses happen slightly before the peaks.

Now, the peaks in such plot are mostly due to the neutron magic numbers, but if one remembers that inside each neutron line 50, 82, 126 there is a doubly magic number, then it is more sensible to examine the proximity of the SM masses to the doubly magic numbers. This is done in figure 1, by using the traditional N-P plot of nuclides. Here we plot diagonal isobaric lines for each massive particle, and we see that effectively the lines happen near the crosses corresponding to doubly magic numbers.

In order to determine how near is ”near”, we need to add some information to the plot. Namely, we will consider the subshell structure before each double magic number.

Also, as our masses are incorporated directly from the standard model, without direct relation to a nuclear model, we could expect to find some signal of them in mass systematics. Indeed the droplet model FRDM 1992 droplet , which is the state-of-the-art mic-mac model, shows strongly this signal for the W,Z, and a weaker error signal for higgs and top, where a quadrupole correction dominates the adjust.

a)P=82, N=126

The mass of the top quark was measured by the Tevatron in 1994, and it amounts, by direct observation, to 174.3 GeV, about 187.1 amu.

The proton magic number 82 is due to the 6h11/2 subshell, which actually is under two other levels, the 3s1/2 and the 4d3/2. The neutron magic number 126 is due to subshell 7i13/2. If we use Klinkenberg 1952 filling scheme k  to draw the rectangle of nuclei with partial fillings of (6h11/2,7i13/2), the diagonal approaches closely to the $`t`$ quark isobaric line. The situation is shown in figure 2.

We should alert the reader that the filling scheme is not well defined because the energy difference between odd and even number of nucleons is enough to alter the energy levels. Thus if we get the shell scheme for 208Pb from Bohr-Mottelson 1969222indirectly quoted in goriely  this year 2003, so it should be still considered state-of-the-art!, then the levels 5f5/2 and 4p3/2 seem to raise above the 7i13/2, and then the isobar line just cuts a corner of the rectangle. This kind of complications are usual and we will find the same issue in the N=82 range.

One should consider also that a 5 GeV range is available starting from the top mass, if we want to consider here all the ”mesons” from $`t\overline{d}`$ to $`t\overline{b}`$. Of course, the highly unstable $`t\overline{t}`$ ”meson” has a mass far away from the nuclear data.

b)P=50, N=50

Here we still have the same subshell for protons and neutrons, namely the 5g9/2. Thus the relevant nuclei are in a square.

The masses of W and Z were measured by CERN collaborations in 1982-3. They amount to 80.423 and 91.1876 GeV, resp. 86.3 nd 97.9 amu.

While both masses are inside the square bracketed by the subshells, some models could also be interested in its average, 92.1 amu, which roughly closes the diagonal of the square.

If one draws the mass lines over the plot of errors in the microscopic-macroscopic mass formula FRDM-1992droplet , it can be noticed that outside from the square the lines coincide with an huge error area in the model. In figure 3 we have included this data.

c)P=50, N=82

Even without the higgs, the two previous coincidences should be enough to consider an extension of nuclear models. Now, a Higgs-like event was reported by LEP-2 collaborations in their final year 2000 run, with a mass of 115 GeV; roughly 123.46 amu. We can plot it and see how near it is of the extant doubly magic number.

The magic number P=50 is due to the subshell 5g9/2 as in the previous case. The magic number N=82 is due to subshell 6h11/2. This should do a $`10\times 12`$ rectangle.

An unexpected complication is the overlapping of neutron subshells 6h11/2, 3s1/2 and 4d3/2. Single nucleons prefer one subshell, while paired nucleons prefer another. It is difficult to establish the right order, if it exists. Even when the 6h11/2 subshell is supposed to be under the other two (or three!), the first excited level is always provided by it.

If we ignore the overlap and we we accept the old criteria from Klinkenberg, putting the 6h subshell as the lower one, then the LEP-2 event simply cuts a corner of the rectangle. The same discrepancy, at a smaller scale, happens if we take Bohr-Mattelson. And if we draw a whole elongated rectangle accounting for all the possible occupations, then it is not so clear where should we draw the diagonal. In any case, we have reflected all these possibilities in figure 4.

On other hand, we should point out that all the modern models based on HFB mass formulae fail to reproduce the narrow overlap and the priority of 4d3/2, and they assign the highest energy level to the 6h11/2 subshell, which then closes the whole shell (check the plots and data in goriely ). This should be the best approach, as then the $`10\times 12`$ rectangle is directly attached to the doubly magic corner.

The area near this mass isobar has been studied also beyond Z=50, in the context of Chiral Doublet Structures ch

## III Analysis and Remarks

Our empirical approach is partly muddled because of the proximity between different subshell splittings. We refer the reader again to k  to get an idea of the complexity of these determinations. On the other hand, it could be that all the levels are actually trying to compete to fit in the isobaric lines.

While it is possible to believe that this effect is just a random numerological coincidence, it covers all the cases where the naive shells plus spin-orbit do not suffice to justify the observed magicity, and specially it suggests that a motivation can exist for the effect of double magicity. There are not more highly massive Standard Model particles, and there are no more double magic nuclei with high spin-orbit splitting. It could be argued that the 28 shell is not included in the list, but this shell needs only a mild spin-orbit split, at the reach of classical models. It is very unusual to find such level of agreement in numerology with congruent units: most numeric coincidences appear in adimensional considerations, or when comparing numbers coming from different systems of units.

The enhancement happens in the depressed subshells, which agrees with the phenomenological approach that enhances the spin-orbit separation by incorporating a dependence on angular momenta. Still, We can not, at this moment, provide an explanation of the effect. The simultaneous fit of both proton and neutron subshells suggest that the effect includes an isospin symmetry.

For the top quark mesons, both J=1 and J=0 particles are available. More specifically333thanks to P. Walters for this remark the pseudoscalar 0-+ should be expected with more probability that the scalar 0++, while the vector 1– is always ”the backbone of any meson spectroscopy”. For the W and Z, we should look for axial effects influencing angular momenta, and also for other combinations explaining the role of W and Z in the error data of macroscopic models. If the LEP-2 event is confirmed to be the Higgs scalar particle, then we should look for couplings to J=0 scalar. It should be remembered that $`J=0`$ coupling is the trademark of the scalar $`\sigma `$-meson in relativistic models of the nucleus. The status of this $`\sigma `$ meson is experimentally troubled. It is usually assimilated to a $`f_0(600)`$ boson (see (pdg, , pg 450)) by the HEP community, but nuclear authors prefer to remark that ”there is no experimental evidence for a free $`\sigma `$ meson, although the $`\sigma `$ field is a crucial ingredient of relativistic mean-field models” scmf . Usually the $`\sigma `$ scalar is seen as a two-pion particle, ie its scalar nature comes from two identical pseudoscalar couplings.

Our particles directly from the Standard Model could interfere with the effect of the $`\sigma `$. It could happen either as an small modification of the $`1/m_\sigma `$ mass term, or as a direct interaction with the nucleon density. It seems from our data that the meson exchange interaction is enhanced when the mass of the whole ensemble $`N+Z`$ of nucleons equals the mass of some available SM boson. This could be more explicit in the second case.

Also, it could be possible to develop some influence in the effective nucleon mass, which, at the end, is the main responsible for the spin-orbit coupling.

On other hand we could be seeing a quantum many-body effect. The interaction propagated via a massive particle will depend both on the mass of the nucleon, $`M`$, and the mass $`m`$ of the boson. But there is also a multiplicative factor depending of the number of nucleons, i.e. of $`M_n/M`$, being $`M_n`$ the mass of the nucleus. So the total interaction can be written in terms of $`(M_n,M,m)`$ and we could expect some surprise444An intermediate suggestion -to be tested- is to calculate the interaction between the whole neutron and proton balls, with respective masses $`M_N,M_P`$, when the interacting particle has $`mM_N+M_P`$. This is suggested because observation of quenching hints that both subshells must be present for the double magicity phenomena to appear. when $`M_nm`$. Besides, $`M_n`$ could be directly the relevant quantity for kinematics.

The strong nuclear force, as driven by the usual system of meson exchanges, provides a kinematical glue because our forces are small compared to the coupling constant from pions. We are working with bound states; any momenta exchanged between a bound particle and the nucleus is shared between all the particles of the nucleus before the ”orbiting” particle gets able to turn around and interact from the opposite direction. So the acceleration from this momentum exchange depends on the mass of the whole nucleus, and this whole mass is more relevant that the mass of the nucleon -or quark- actually responsible of the exchange. It happens also, in some classical inelastic collisions, that the conversion of energies $`\text{collision}\text{internal}+\text{mass center}`$ depends on the coupling constant of the internal force. The perturbation is to be measurable because in a yukawian potential it depends on the product of mass (inverse range) times coupling.

In any case, if the recoil depends of the whole mass, it is difficult to justify why should one to cut-off the calculation at the scale of the single nucleon. Regretly, most meson theories are effective QFT with a cutoff $`\mathrm{\Lambda }<`$ 1 GeV, so they are not directly suitable to incorporate the highest particles.

With a clear background, and in the absence of serious Effective Field Theory motivations, the ratio $`\mathrm{\Gamma }/m`$ between pole width and pole mass of a given meson should be a consideration more fundamental than simply to look if such mass is above or below the one of the proton. It should be of some help if similar effects were searched in few nucleon systems. The scale of J/psi, B and Upsilon could be an extra contribution respectively related to the variations of mass at 4He, 8B, 12C, or to contribute to the valleys between, or simply to be hidden under other effects. Furthermore, perhaps here the nucleon number is too small or the background of other particles too big. So the most we can say is that the low mass range does not contradict our observation, but that it is not a definitive support. Ab-initio models could be able to fit the data, in the future, without requiring the mildly massive mesons.

Beyond theoretical frameworks, another clue to this effect can arrive from an analysis of the errors in current models of masses. Anomalies are to be expected near mass numbers 86, 98, 123, 187 (W, Z, Higgs, tops). We have remarked this point in figure 3. Indeed the error analysis of mexico  for the same model hints the same results, and also the variations between theoretical and real error in fig 7 of droplet . Besides, the recent analysis of Huertas for FST model shows slope changes in the error at this mass numbers, see figure 4 of huertas1  and fig. 13 of huertas2 .

Of course it shall not escape to any reader that a Mic-Mac model including high SM masses as free parameters could be used to predict these SM values, via least squares or any other variational fit. For instance one can calculate the model parameters dependence on the masses of the nuclei used in the fit; say $`g_>(m)`$ fitting to all nuclei greater than $`m`$, or $`g_<(m)`$ fitting for all the nuclei with mass smaller than $`m`$, or even some intermediate fit $`g_{+10}(m)`$. Then any sudden jump in $`g(m)`$ could be used to signal a unexpected effect around mass $`m`$.

Our work is in some sense falsifiable at mid term: If either no particle appears near 115 GeV or QHD proves able to reproduce the spin couplings for double magicity, the coincidence will be invalid in the first case, accidental in the later. A ten-years lapse will do the final judgement, as computing power increases and the LHC does its hunt.

The author wants to acknowledge L.J.Boya and other members of the Theoretical Physics Department at Zaragoza for discussions, as well to P. Walters at the internet physicsforums pforums .

## IV appendix

As of today, the previous comment has been rejected from practically all the letter journals of its thematic555see my webpage for details on this mundane thing. The rejection has been mostly at editorial level, so I have not got any useful input from referees nor pairs.

During this time I have mostly tried some theoretical preparation, but I can not report of serious advance in this side. On the other side, the empirical, we are luckier, perhaps. I have mentioned in the paper that we should expect huge error peaks in mass models that are unable to make place for the particles we are looking for. Indeed we have shown this effect for W and Z in an area of the FRDM model. But this model does not show huge discrepancies for 115 and top. Do our conjecture fails or, more probably, the model is able to adjust itself to nullify the discrepancy? Last week I did some plottings that hint towards the second possibility.

I took upon myself the task of plotting error graphs for some tables available in the net. I started with the FRDM table, and a mistake when manipulating the data in gnuplot produced the result you can see in figure 5.

I was mixing different column data in the same plot, but still getting an intriguing result. By examining the data, it happened that I was plotting a conditional: for some nuclei the model uses a kind of shape correction, for some others it chooses another. I went to the original paper droplet  to check the selection of nuclei, and I found the coincidences shown in figure 6

The apparition of the electroweak vacuum (the expected value of the higgs field) is surprising because it does not imply a standard model particle. The discarded $`ϵ_6`$ in this area is qualitatively different from the others (fig. 7). Besides, the 246 Gev range is also noticeable in the -very noisy- plot of calculated microscopic energies for the twin model FRLDM. (figure 8).

The emerging role of the neutron drip line is consistent with out thinking of interpreting the rest of the nucleus as a single particle: at the dripline, the neutron skin orbits the rest of the nucleus in a dense version of the halo of light drip line nuclei.

We see that this concrete model uses shape corrections to fill the expected discrepancies. But what about other models? Results vary. Generically, when the error plot has noticeable peaks, the pattern of electroweak scales appears. By looking at otros , about two thirds of the models present very noisy error patterns, but still we have another third showing curious peaks at the appropriate ranges.


text:
# Simulations of Weak Gravitational Lensing

## 1 Introduction

Weak gravitational lensing by large-scale structure has become an indispensable tool for modern cosmology, used already to set constraints on the mass density ($`\mathrm{\Omega }_{\mathrm{mat}}`$, in units of the critical density) and the fluctuation amplitude $`(\sigma _8)`$ (see e.g. van Waerbeke & Mellier, 2003; Hoekstra, Yee, & Gladders, 2002, for the current status) and touted for its potential to constrain cluster scaling relations (Huterer & White, 2002) and dark energy (Benabed & Bernardeau, 2001; Huterer, 2002; Hu, 2002; Heavens, 2003; Abazajian & Dodelson, 2003; Refregier, 2003; Jain & Taylor, 2003; Bernstein & Jain, 2004; Takada & Jain, 2004; Takada & White, 2004). Like the anisotropies in the cosmic microwave background (CMB), the theory of weak gravitational lensing is well understood (Mellier, 1999; Bartelmann & Schneider, 2001). Unlike the CMB however the calculation involves modeling the non-linear evolution of the mass in the universe. This makes the predictions of weak lensing very rich, but also means that an accurate treatment requires numerical simulations.

In this paper we give some details of a pilot program designed to provide weak lensing ray-tracing simulations of a small grid of cosmological models in the currently favored region of parameter space. These models can be used by observers wishing to test their algorithms or fit to their data and by theorists wishing to test or calibrate fast and flexible, but approximate, methods of calculation. The initial grid is small, but the process of its creation is almost entirely automatic allowing it to be expanded simply as the need arises.

This paper focuses on the creation of the model grid. We describe the choice of cosmological parameters in §2.1, the N-body simulations in §2.2 and the ray tracing and numerical issues in §2.3 and §2.4. We give some representative results in §3 before concluding in §4.

## 2 Cosmological models

### 2.1 Choosing parameters

We choose a small number of cosmological models which provide good fits to the current CMB and large-scale structure data. For simplicity all of the models in the pilot program are variants of the cold dark matter model with a dark energy component, and their parameters are shown in Table 1. The parameter variations around a base model are designed to keep the CMB fluctuations almost invariant, in anticipation of increasingly precise data from WMAP and Planck, and the entire process is automated using Perl scripts and C code.

We begin with a relatively small number of models to demonstrate the cost and feasibility of making such grids. As we gain experience in using these data products and understand the drivers better we can extend the model grid and/or increase the fidelity of the simulations. Are our priorities to have larger maps? higher resolution? more redshift range? more ‘sky’? etc.

We first pick the physical matter density $`\omega _{\mathrm{mat}}\mathrm{\Omega }_{\mathrm{mat}}h^2`$, the (comoving) distance to last scattering, $`d_{\mathrm{ls}}`$, and a (constant) equation of state of the dark energy, $`w`$. We approximate the redshift of last scattering as $`z=1080`$, ignoring the slight matter and baryon density dependence. This then allows us to solve for the Hubble constant, $`h`$, and thus $`\mathrm{\Omega }_{\mathrm{mat}}`$ and $`\mathrm{\Omega }_{\mathrm{de}}`$ assuming spatial flatness. Because they are reasonably well known, compared to some other parameters, we fix $`\omega _{\mathrm{mat}}=0.145`$ and $`d_{\mathrm{ls}}=13.7`$Gpc for all of the models in our grid. These numbers are close to the best fit for a recent analysis of WMAP and SDSS data (Tegmark, 2003).

The next step is to specify the other model parameters, for example the optical depth to Thomson scattering, $`\tau `$, the spectral index, $`n`$ and the baryon density $`\omega _b`$. (The optical depth enters primarily through its effect on the normalization of the power spectrum when we fit to WMAP.) We currently hold the number of light neutrinos fixed, and include no massive neutrinos. We deal with pure power-law spectra with no running spectral index. Our primary variations are in the optical depth, the spectral index and the equation of state of the dark energy (see Table 1). The first two affect the amplitude and shape of the density or potential fluctuations at late times while $`w`$ is a parameter of great interest to the cosmology community. We use a fixed equation of state in this initial survey, though it is easy to include any known functional form in the future.

For a specified model we compute the CMB anisotropy spectrum using v4.5 of CMBFAST (Seljak & Zaldarriaga, 1996). The likelihood software provided by the WMAP team (Verde et al., 2003; Hinshaw et al., 2003; Kogut et al., 2003) and the known anisotropy spectrum are then used to find the best fitting normalization of the power spectrum, $`\sigma _8`$, which is then converted into internal code units for the N-body simulations. In this way a file of model parameters, shown in Table 1, is built up. The lensing convergence angular power spectra for our 4 models, assuming the Born and Limber approximations, are shown in Fig. 1. Note that the spectra for each pair of models have similar shapes, differing primarily in amplitude. In each case the model with the most distant sources has the highest amplitude, even though it has a slightly lower matter power spectrum amplitude, $`\sigma _8`$. There will also be differences in the higher order moments and the source redshift dependence of the power spectra.

### 2.2 N-body simulations

Another Perl program reads this cosmological parameter file and generates a directory structure containing the requisite auxiliary files for the N-body simulation, plus scripts for generating the initial conditions and running the N-body code on the IBM-SP Seaborg at NERSC.

All simulations require choices to be made about box size, resolution and sampling. To minimize numerical artifacts we want to run the largest box with the highest dynamic range consistent with the computational resources. Since this is a pilot program we have chosen simulations which use $`256^3`$ particles in a $`200h^1`$Mpc box. As we shall discuss in §2.4, this is expected to be sufficient for our purposes now, but will need to be improved upon in the future. As further resources become available we can increase the volume simulated and/or decrease the force softening.

The initial conditions are generated by displacing equal mass particles from a “fuzzy” Cartesian grid at $`z=50`$ using the Zel’dovich approximation. The initial fluctuations are Gaussian, however we choose from a larger number of initial conditions those whose low few $`k`$-modes are close to the mean. This avoids running any simulations which are significant outliers and may skew statistics for a small number of runs. For each cosmological model we choose two independent initial conditions. Ideally we would run more realizations of each model to increase the volume simulated, but with two runs we will be able to simulate several hundred square degrees which is enough for now.

The TreePM code, described in White (2002), is used to evolve the particles to $`z=0`$ with the phase space information output at selected times after $`z=3`$. The outputs are spaced equally in conformal time, with a spacing equal to the time taken for light to travel $`50h^1`$Mpc (comoving). This spacing is small enough to obtain good lensing maps – a halo with “typical” velocity $`300`$km/s will move only $`50h^1`$kpc between outputs, and its departure from straight line motion will be small. Each simulation gives between 80 and 90 phase space dumps, depending on the cosmology, each of $`500`$MB. The gravitational force softening is of a spline form, with a “Plummer equivalent” softening length of about $`30h^1`$kpc, again small enough for our purposes.

Since the box size is only $`200h^1`$Mpc, leading to a relatively coarse sampling in $`k`$, we generate the initial conditions from a fit to the transfer function which has the baryon oscillations “smoothed out”. For now we use the fit published by Eisenstein & Hu (1999). This fit was made to an older version of CMBFAST, and since that time improvements in the code and the physics (Seljak et al., 2003) have changed the predictions for $`T(k)`$ slightly. The differences are small, and for the purposes of testing and calibrating algorithms it is enough to know the input $`T(k)`$ precisely, but in the future we will improve this aspect.

### 2.3 Making lens maps

We make two sets of maps. The first, which we shall call the Born series, assumes that the convergence field, $`\kappa `$, is simply the integral along the line of sight of the density field weighted by a simple kernel

$$\kappa \frac{3}{2}H_0^2\mathrm{\Omega }_{\mathrm{mat}}𝑑\chi g(\chi )\frac{\delta }{a}$$
(1)
where $`\delta `$ is the overdensity, $`a`$ is the scale-factor, $`\chi `$ is the comoving distance and $`g(\chi )`$ is the lensing weight

$$g(\chi )_\chi ^{\mathrm{}}𝑑\chi _sp(\chi _s)\frac{\chi (\chi _s\chi )}{\chi _s}$$
(2)
for sources with distribution $`p(\chi _s)`$ normalized to $`𝑑p=1`$. This set of maps has the advantage of not requiring a Fourier transform in its construction, providing higher resolution for a given pixelization. Its principal disadvantage is the approximate nature of the calculation.

We make several sets of 16 quasi-independent $`3^{}\times 3^{}`$ maps with $`1024^2`$ pixels. The first set of maps has sources fixed at $`z_s1`$, while the second set uses a source distribution of the form (Brainerd, Blandford & Smail, 1996)

$$\frac{dp}{dz_s}z_s^2\mathrm{exp}\left[(z_s/z_0)^{3/2}\right]$$
(3)
for $`z_0=2/3`$ and $`1`$. For this distribution $`z=\mathrm{\Gamma }(\frac{8}{3})z_01.5z_0`$, but because there is more of the lensing weight at higher $`z`$ the lensing signal is smaller for this distribution than a $`\delta `$-function distribution with $`z_s=1.5z_0`$.

The second set of maps comes from a full multi-plane ray tracing algorithm as described in Vale & White (2003). We make maps of the shear components, $`\gamma _i`$, and the convergence, $`\kappa `$, at a range of source redshifts from $`z0`$ to $`3`$ in steps of $`50h^1`$Mpc. Again we make 16 maps per simulation, choosing different positions and orientations for the ‘observer’. In each case a $`2048^2`$ grid of rays subtending a field of view of $`3^{}`$ is traced through the simulation, with the required Fourier transforms being done on a $`2048^2`$ grid. All assignments to and from the Fourier grid are done with cloud-in-cell (CIC) assignment (Hockney & Eastwood, 1980).

The two shear components and the convergence are output at each source plane. From this a map of the (measurable) reduced shear, $`\gamma _i/(1\kappa )`$, for any source distribution can be computed. For a distribution $`dp/dz_s`$ the weight given to source plane $`j`$ is

$$w_j=\frac{dp}{dz_s}|_jH(z_j)\mathrm{\Delta }\chi $$
(4)
where $`\mathrm{\Delta }\chi =50h^1`$Mpc is the spacing between outputs. The final shear at each point is then $`\gamma _i=_jw_j\gamma _i^{(j)}`$ with $`\gamma _i^{(j)}`$ the $`i`$th component of the shear measured on plane $`j`$. We show an example of the weighting factors for $`z_0=1`$ in Fig. 2, along with the (trivial) weight for a $`\delta `$-function source distribution at $`z_s1.5`$. The integral over $`\chi _s`$ is extremely well approximated by a sum over the output planes. These maps can thus be combined to produce lensing maps for a wide range of source distributions, allowing source tomography to be used. It is straightforward to introduce mock ‘galaxies’, with an appropriate distribution of intrinsic ellipticities, at this stage and produce mock catalogues. On each $`\chi _s`$ slice, the angular positions of the galaxies can be made to trace the dark matter with some fidelity, giving rise to source clustering effects. Reasonable choices produce angular correlations and redshift distributions for the ‘galaxies’ which match current data well. However the parameter space of such catalogues is large, so we have chosen instead to make available simple co-added shear maps which can be post-processed with a wide variety of simulated galaxy properties to make shear maps. To save space we provide the maps downsampled from $`2048^2`$ to $`1024^2`$.

### 2.4 Numerics

To understand the fidelity of the map making procedure we need to account for the finite volume, force and mass resolution of the N-body simulation and the effect of the finite sampling during ray tracing. The comoving distance to our furthest sources is $`3`$Gpc, at which distance our box subtends $`3.8^{}`$. The volume of space probed by a lensing field is $`10^7\mathrm{Mpc}^3`$, comparable to the volume simulated in each run.

In Vale & White (2003), we developed an analytic model to estimate the effects of finite force and mass resolution and the Fourier grid used in the ray tracing. This model was based on modifying the 3D dark matter power spectrum that appears in the Limber integral for $`C_{\mathrm{}}`$ in a way designed to fit the effects seen in a series of N-body simulations. This allows us to compare a given numerical configuration to a “perfect” model. Fig. 3 illustrates these effects individually for model 1 of Table 1 with sources at $`z_s1`$. For simplicity we use the power spectrum fit of Smith et al. (2003) in computing these ratios. Note that the suppression of power by the finite force resolution and FT grid is partly countered by the artificial increase in power due to shot-noise. Based on Fig. 3 we expect our maps to be accurate to several percent, in the 2-point function, for $`\mathrm{}3000`$ and even higher in some situations where the cancellation is accurate.

While the 16 maps per model are not fully independent, they do sample different regions of the box in different stages of evolution and with different projection effects. The fraction of the volume traced by the rays in the field of view, weighted by the contribution to the variance of the convergence from each box, is between 10-20% (for $`100<\mathrm{}<3000`$). This forms a rough upper limit to the degree of correlation between the maps. Including the two independent boxes per model we have simulated several hundred square degrees of sky per cosmology, which is close to the amount of observational data currently existing. Extending this to more sky is a straightforward exercise that simply requires more computing resources.

## 3 Results

These simulations and lensing maps can be used for a variety of purposes. Here we simply show some preliminary analysis to orient the reader. The real uses of these maps will be in testing analysis algorithms and new ideas.

Before turning to the lensing maps let us consider the underlying mass distribution. Fig. 4 shows the 3D matter power spectrum from our 8 boxes at $`z=1`$ and $`z=0`$ compared to the semi-analytic model of Smith et al. (2003). The agreement is good for both the $`\mathrm{\Lambda }`$CDM models and the models with equation of state $`w=0.8`$, as expected. On the scales simulated the dark energy only enters through the Hubble parameter in the drag term, slowing the growth of large-scale fluctuations when the dark energy comes to dominate the expansion. Since the fit of Smith et al. (2003) has been tested for open and $`\mathrm{\Lambda }`$CDM models it is no surprise that it works equally well for $`w=0.8`$. The disagreement at the highest $`k`$ might be due to numerical artifacts in the N-body simulation or in the computation of the power spectrum, however these points are significantly above both the shot-noise and force smoothing limits. We find similar behavior in a number of higher resolution simulations we have analyzed, and a similar level of disagreement can be found in some of the figures of Smith et al. (2003) for CDM models and may reflect inherent inaccuracy in the fitting function. A qualitatively similar level of agreement was also found between particle-mesh simulations and the semi-analytic model by Ishak et al. (2003), though they were unable to probe very high $`k`$-modes due to limited force resolution.

As a further check on the N-body code we investigated the scaling of the power spectrum of an $`n=1`$ self-similar model with the same numerical parameters as the models described. This suggests that the power spectrum should be accurate to a few percent on the range of scales of interest.

The $`\kappa `$ angular power spectrum produced from 16 maps each for models 1 & 2 is shown in Fig. 5. For this comparison we show maps made assuming the Born approximation and with full ray-tracing, both for sources fixed at $`z_s1`$. The agreement is generally good, with the variance at low-$`\mathrm{}`$ expected from sampling 16 $`3^{}\times 3^{}`$ fields. The open circles, corresponding to model 2, are slightly higher than the semi-analytic model at intermediate $`\mathrm{}`$, as expected from the excess power seen in the top panel of Fig. 4 (the triangles) at intermediate $`k`$. At the highest $`\mathrm{}`$ the ‘Born’ maps show a slight excess power which could be due to shot-noise in the simulations or to a breakdown in some of the approximations we have made. Alternatively the semi-analytic model and the ray tracing could be underestimating the power on these scales, as suggested by Figs. 3 and 4. In measurements these scales would be significantly affected by finite galaxy ellipticities and densities, so this disagreement is likely unimportant. Fig. 6 shows the same comparison for models 3 and 4. The power spectra for the remaining four models are very similar and we do not show them explicitly here.

The maps are non-Gaussian, as expected, and we provide the skewness

$$S_3^\kappa \frac{\kappa ^3}{\kappa ^2^2}$$
(5)
for each of the models in Fig. 7. The $`\kappa `$ field is smoothed by a boxcar of side length $`\theta `$ before the moments are computed. The plot shows the average and the variance across the 16 maps as before. As noted previously (White & Hu, 2000; Vale & White, 2003) the skewness suffers from a large sample variance.

We also show the shear 3-point function for isoceles triangles of side length $`\theta _1=\theta _2=2^{}`$ in Fig. 8. Here we plot $`10^6\zeta _{+++}`$ as a function of the cosine of the included angle for the 32 fields, each $`3^{}\times 3^{}`$, simulated for each cosmology assuming the sources are at $`z_s1`$. The points are the mean of the 32 fields, and the error bars represent the variance. The highest 3-point function is for models 1 & 2, which have the highest normalization ($`\sigma _8`$) of all of the models. Models 3 & 4 have the next largest amplitude, and the second highest normalization. For models 5 & 6 the effect of a higher normalization isn’t enough to overcome the less negative $`w`$ of models 7 & 8, which have a slightly higher amplitude. Notice the large scatter from field to field, but the familiar shape (Takada & Jain, 2002). Apart from an amplitude change the shapes of the curves are extremely similar, as shown in Fig. 9. We find a similar behavior for other configurations we have checked, which has interesting implications for the signal encoded in the non-Gaussianity of the maps and how best to extract it. We also note that the peak of $`\zeta _{+++}`$ appears to be shifted slightly from the equilateral configuration, due to the asphericity and substructure of the dark matter halos (Ho & White, 2003; Dolney, Jain & Takada, 2004).

## 4 Conclusions

Gravitational lensing has made rapid advances in recent years, and observations are now pushing the limits of existing theoretical models. If lensing is to fulfill its promise as a precision tool, theorists need to improve the predictions of cosmological models for the $`n`$-point correlation functions. An important step in this program is the creation of grids of N-body models of sufficient resolution and sampling to enable accurate simulation of weak lensing. This paper describes a first step in this direction.

Based on analytic arguments developed in Vale & White (2003) we believe that simulations with $`256^3`$ particles in boxes of side $`200h^1`$Mpc are sufficient to produce lensing maps accurate to $`\mathrm{}20003000`$ or scales of a few arcminutes. To achieve good convergence in the multi-plane ray-tracing algorithm that we use requires time dumps spaced equally in conformal time with a spacing close to $`50h^1`$Mpc.

The total computational cost of this project was very modest. The major human cost is the time spent managing the disk space, archiving and retrieving the data to ensure that the total usage stays below quota. The simulations each took 2000-3000 time steps, requiring slightly under one day each on 32 processors of the IBM-SP Seaborg at NERSC, i.e. between 640 and 720 CPU hours. The phase space data required $`40`$GB per model while the ray tracing simulations required around $`100`$GB of intermediate storage and 100 CPU hours per model. The disk usage is driven by the need for all the phase space data and the $`8090`$ high resolution source planes for each of 16 maps. The final stacked and downsampled maps are $`<1`$GB per model and source distribution. The entire process was monitored and controlled by scripts, making it easy to extend the grid as more resources become available.

We have made the raw maps, along with some auxiliary data products, freely available to the community at http://mwhite.berkeley.edu/ in the hope that they will be useful in taking the next step.

M.W. thanks Tzu-Ching Chang, Dragan Huterer, Bhuvnesh Jain, Jason Rhodes and Masahiro Takada for helpful comments on an earlier draft. The simulations used here were performed on the IBM-SP at the National Energy Research Scientific Computing Center. This research was supported by the NSF and NASA.


text:
# Linear polarization on Gamma-Ray Bursts: from the prompt to the late afterglow

## 1 Introduction

Linear polarization has revealed to be a characteristics of GRBs throughout their entire evolution. The recent claim by Coburn & BoggsCB03  that the prompt emission of GRB 021206 was polarized at a very high level has stimulated a thorough analysis of the polarizing properties of the jet geometries in GRB outflows, drawing attention on the possibility that magnetic fields may be advected from the central source rather than generated by the internal shock. On the other hand, observational and theoretical studies of afterglow polarization has revealed a much more complicated picture than previously thought, emphasising the importance of the jet structure, its dynamics and the properties of the ISM in shaping the polarization curves. Even though direct observations lack, also the optical flash may be highly polarized, especially if due to a reverse shock in the burst ejecta rather than to the pair enrichment of the nearby ISM.

In this paper I review the status of theories and observations of linear polarization in GRBs. The three phases are analysed initially separately and then their relation discussed. The importance and insight of polarization is emphasised in all phases.

## 2 The prompt phase

Analysing the scattering geometry of photons in the RHESSI detector, Coburn & Boggs CB03  were able to measure the average linear polarization of the prompt emission of GRB 021206 in the \[150 keV–2 MeV\] energy range. They find that the prompt emission of the burst is highly polarized, with $`\mathrm{\Pi }=0.8\pm 0.2`$. This measurement was subsequently heavily criticised by Rutledge & Fox RF03 , who performed an independent analysis of the same dataset, obtaining a much smaller number of double-scattered photons and, as a consequence, merely an upper limit on the linear polarization of the event. Despite that, the result CB03  has stimulated a vast theoretical effort in order to understand under which conditions such a large polarization could be obtained. Two classes of models have emerged. In the first class, the origin of polarization is ascribed to the presence of a large scale ordered magnetic field, which ought to be advected from the central engine and may play a role in the launching of the jet itself Ly03 . In the second class of models, the magnetic field is supposed to be shock generated and tangled on small timescales, and the asymmetry required to produce polarization is due to a particular location of the observer with respect to the jet axis W03 . This second class of models can be extended to different emission mechanisms, such as inverse bulk Compton scattering L03a .

### 2.1 Magnetic models

We define here magnetic models those in which polarization is due to the large scale geometry of the magnetic field. The magnetic field is likely to be dominated by a toroidal component, since the radial field decays faster than the tangential one ($`B_rr^2`$ while $`B_{}r^1`$). One important ingredient of these models is that the observer, due to the relativistic aberration of photons, cannot see the whole jet. In fact it is only a small $`1/\mathrm{\Gamma }`$ region of the jet that is observable and therefore the observer does not detect the overall toroidal structure of the field (which would wash out the polarization signal) but a highly ordered patch. This is most important if the magnetic field is not the dominant component of the outflow. Since regions of the jet separated by more than $`1/\mathrm{\Gamma }`$ are causally disconnected, it is difficult to envisage a coherent magnetic field on scales larger than $`1/\mathrm{\Gamma }`$, unless the structure has been created before the acceleration of the jet (when it was still connected) and frozen into it. Such a transport seem easier to attain in a magnetic dominated outflow P03  and even natural in a force-free subsonic bubble Ly03 .

Even if the observer has access only to a fully ordered region of field, the polarization signal cannot be as large as that expected from a non relativistic flow. In the classic case: $`\mathrm{\Pi }=(p+1)/(p+7/3)`$ where $`p`$ is the power-law index of the electron energy distribution. The reduction of the observed polarization in the relativistic case is due to the aberration of photon trajectories. In order to keep the electric and magnetic field of the wave orthogonal to each other and to the wave propagation direction, the position angle of polarization is rotated in different ways as a function of the distance from the line of sight. After integration, the polarization is reduced by a factor that depends on the spectral slope of the radiation, and spans between $`10\%`$ and $`20\%`$. The non-relativistic solution is shown with the relativistic maximum polarization in the left panel of Fig: 1 (see also G03 ).

A characteristic feature of these models is that any observer located within the opening angle of the fireball detects a highly polarized signal, with the exception of those observing the fireball within $`\theta =1/\mathrm{\Gamma }`$ from the symmetry axis. A level of polarization comparable to that detected by Coburn & Boggs cannot be achieved, even though given the large uncertainties a definitive conclusion cannot be drawn.

### 2.2 Geometric models

It has been traditionally assumed that the magnetic field responsible for the synchrotron emission in GRBs is generated at the shock front. This is a robust conclusion in the afterglow phase, where the compression of the interstellar field is far too low to produce the observed radiation. It may hold true also for the prompt phase, in which case a tangled field would be responsible for the observed radiation. If this field is tangled in a plane, but compressed in the direction perpendicular to the plane itself L80 ; M99  it is possible to observe polarized synchrotron radiation since radiation emitted in the plane, which is maximally polarized in the comoving frame, is then aberrated toward the observer with an angle $`\theta =1/\mathrm{\Gamma }`$ GL99 .

In the afterglow phase, this configuration leads to polarization of up to several tens of per cent GL99 ; S99 , but in the prompt phase is usually negligible, unless a narrow jet is observed along the required direction W03 . Polarization in this context has been analyzed in several works with synchrotron as the radiation mechanism G03 ; N03  as well as if the photons are produced by bulk inverse Compton scattering L03a . The difference between the two cases is not of fundamental nature, since the dependence of polarization on angle is the same for the two mechanisms. Inverse Compton, on the other hand, can produce larger polarization since it can be maximally ($`100\%`$) polarized in the comoving frame.

The polarization produced by a narrow jet is shown (for the case of inverse Compton) in the right panel of Fig. 1. In these models a narrow jet is fundamental since the number of observer that see a polarized event is limited to those lying in the region between the edge of the jet and $`1/\mathrm{\Gamma }`$ from it. This region becomes vanishingly small for $`\theta _{\mathrm{jet}}>10/\mathrm{\Gamma }`$. GRB 021206 was exceptionally bright and, assuming it was at cosmological redshift, would have had a narrow jet with opening angle of few degrees at most.

To distinguish between magnetic and geometric models is quite easy once a reasonable number of measurement is available. In the geometric case only a small fraction of the brightest bursts should be polarized, while in the magnetic case most of them should.

## 3 The flash

It has been suggested by Granot & Königl GK03  that the reverse shock emission could be as highly polarized as the prompt GRB, since the plasma responsible for this emission is the same one that produced the gamma-ray photons. This consideration can be included in a more general discussion on optical flashes, i.e. bright optical components that appears at the beginning of the afterglow phase with a fast decay.

Optical flashes can be produced in two ways: by a reverse shock in a baryonic jet or by pair enrichment of the external medium in the vicinity of the GRB B02 . In most cases the optical flash is expected to be polarized at a level comparable to the prompt emission. There is actually only one case in which the optical flash following a polarized GRB can be unpolarized. In a magnetic model, if the flash is due to pair enrichment of the ISM, the flash comes from a shock generated field without the geometric constraints to produce polarization. In all geometric models the geometry of the prompt phase is preserved during the flash emission. Finally, if the flash is due to reverse shock, it should be polarized as discussed by Granot & Königl GK03 .

## 4 Afterglow

The discovery of linear polarization in GRB afterglows dates back to 1999, when a small but highly significant level of polarization was detected in GRB 990510 C99 . The detection took place amid a theoretical effort to predict and/or explain it.

Guzinov & Waxman Gr99  discussed the possibility that the shock generated field organises in coherent patches that expand at a sizable fraction of the speed of light. They calculated that an observer should see approximatively $`N\stackrel{>}{}50`$ patches. If the polarization inside a patch is $`\mathrm{\Pi }_070\%`$, the observed one is $`\mathrm{\Pi }=\mathrm{\Pi }_0/\sqrt{N}\stackrel{<}{}10\%`$. Indeed the degree of linear polarization observed in GRB afterglows is in the per cent range (see C03  for a review). However, given the random nature of the model, the degree and position angle of polarization fluctuates in time. This is not, at least in some cases, detected in polarization curves L03b ; Gr03  and for this reason this model is now not considered particularly promising.

Shortly after the discovery of polarization Ghisellini & Lazzati GL99  and, independently, Sari S99  proposed a model based on the assumption that the fireball is beamed in a come and that the shock generated field is either compressed in the shock plane or elongated in the axial direction. Polarization is observed if the observer is not coincident with the cone axis and has a definite and testable pattern. Polarization at early times is null, increases slowly with time until it reaches a maximum and then starts to decrease again until it vanishes. At this moment, which is roughly coincident with the jet break in the unpolarized lightcurve, the position axis of polarization rotates by $`90^{}`$. Then the polarization curve is characterised by a second peak, of higher intensity, eventually vanishing to an unpolarized flux at long times. The intensity of the polarization signal depends on the off-axis angle: the larger the off-axis angle the larger the polarization.

These models have been further analyzed and extended by Rossi et al. R03 . They studied the effect of different assumptions on the jet sideways expansion showing that the ratio between the peaks of polarization is smaller for faster expansion speed. They also generalised the model to non uniform jets, in which the energy per unit solid angle decreases as $`E_\mathrm{\Omega }\theta ^2`$ where $`\theta `$ is the angle with respect to the flow symmetry axis. In this case the polarization curve is largely different. While it is null for early and late times, as in the uniform case, it has a single peak, correspondent in time with the jet break time, and constant position angle. In Fig 2 we show the comparison between the polarization curves for the homogeneous and the structured models. It is clear that the two curves are different in an easily testable way. This is particularly important if we consider that the light curves of the two models are almost indistinguishable. Further complications to the models have been added by GK03  who consider the presence of a coherent component of the magnetic field in the ISM. The propagation of the polarized light of the OT in the host and Galactic ISM have been instead discussed by L03b .

Comparison of the models with the data has proven difficult. The main limitation of these models is that they assume that the emissivity of the fireball is uniform (or strictly dictated by the $`\theta ^2`$ law). Any deviation from this assumption, or inhomogeneity of the external medium, causes a noise on top of the models in both¡ polarization and position angle. Usually this situation is recognisable in the light curve through the presence of bumps and wiggles on top of the regular power-law decay. Indeed, every time the light curve is complex, the polarization curve has a complex structure, such as in GRB 021004 Rol03 ; L03b  and in GRB 030329 Gr03 . On the other hand, simple polarization curves seem to be associated to power-law afterglows (GRB020813 Gor03 ).

## 5 Comparison of the three phases

The position angle of polarization should be related in the three phases. In the intrinsic models (i.e. neglecting polarization induced by the ISM and by an external magnetic field) the position angle of the polarization can be either contained in or orthogonal to the plane containing the jet axis and the line of sight. It is therefore expected that, should polarization be measured in the future in the three phases of a single GRB event, the position angle should either remain constant throughout the whole evolution or rotate by $`90^{}`$ between the optical flash and afterglow. It may eventually rotate back to the original position. Any difference from this simple behaviour should be considered a sign of an external component in the generation of polarization.

## 6 Summary and conclusions

The study of polarization evolution in GRBs is highly informative, albeit difficult. It carries important informations about the jet structure and dynamics that are hidden in degeneracies of the light curve, but are emphasised in the polarization curve. Observationally, the afterglow phase is the most simple to investigate, even though is may be affected by small-scale inhomogeneities, and constraining the smoothness of the light curve is of fundamental importance in order to model a polarization curve. Polarization in the prompt and optical flash emission is highly informative of the structure of the ejecta, even though further observations are required in order to establish the mere existence of polarization in this phases.