format: text

text:
# Nonlinear Amplitude Maxwell-Dirac Equations. Optical Leptons

## 1 Introduction

In the last four decades there has been considerable interest in the nonlinear generalizations of the quantum field equations  and in the possibility of obtaining exact stationary solitary solutions of the field equations . As a rule, the nonlinearity has been introduced in an ad-hoc fashion in the Klein-Gordon equation and also for all four spinor components of the Dirac equations. For the usual case of a cubic nonlinearity, exact localized solutions are not found. Our present work, reported in this paper, shows that the optical analogy of nonlinear Dirac equations leads to a nonlinear part only the first coupled equation of Nonlinear Dirac Equations. This result allows to solve the NDE by separation of variables and to obtain solutions representing optical vortices with classical momenta one and spin one-half. The initial investigation of optical vortices began with a scalar theory, based on the well-known 2D+1 paraxial Nonlinear Schrodinger equation (NSE) . The existence of optical vortices was predicted in the self-focusing regime, but as it was shown in many papers, the solutions obtained are modulationally unstable. In spite of this, various interactions (attraction, repulsion, fusion) have been observed. The scalar paraxial approximation is valid for slowly varying amplitudes of the electrical field in weakly nonlinear media. As it was pointed out in Refs. this theory is not valid for a very intense narrow pulses. The first generalization of the scalar paraxial theory of optical vortices is based on investigation of so called spatio-temporal scalar evolution equations . For all cited theories no exact solutions have been found, but numerical and energy momentum techniques are used. The existence of exact stable vortex solutions of these types of nonlinear equations was finally discovered with the vector generalization of the 3D+1 Nonlinear Schrodinger equation . It also have been shown numerically that these vortices are stable at distances comparable to those where localized solutions of the one component scalar equation self-focusing rapidly. All of these exact vortices are a combination of linearly-polarized components and have spin $`\mathrm{}=1`$. To extend the theory to vortices with spin $`j=1/2`$ we return to an analogy between the Maxwell and Dirac equations . As it is shown in Ref.  this analogy is possible only if the electrical and magnet components are represented as a sum of linear and circularly polarized components.

## 2 Maxwell’s equations with non-stationary linear and nonlinear polarization

Consider the Maxwell’s equations in the next two cases:

$`1`$. A source-free medium with non-stationary linear and nonlinear electric polarization and non-stationary magnetic polarization (the case with dispersion).

$`2`$. A source-free medium with stationary linear and nonlinear electric polarization and stationary magnetic polarization (the case without dispersion).

For these cases, the Maxwell’s equations can be written:

$`\times \stackrel{}{E}={\displaystyle \frac{1}{c}}{\displaystyle \frac{\stackrel{}{B}}{t}},`$ (1)
$`\times \stackrel{}{H}={\displaystyle \frac{1}{c}}{\displaystyle \frac{\stackrel{}{D}}{t}},`$ (2)
$`\stackrel{}{D}=0,`$ (3)
$`\stackrel{}{B}=\stackrel{}{H}=0,`$ (4)
$`\stackrel{}{D}=\stackrel{}{P}^{lin}+4\pi \stackrel{}{P}_n\mathrm{}.`$ (5)
The linear magnetic polarization corresponds is :

$`\stackrel{}{B}=\stackrel{}{H}+4\pi \stackrel{}{M}_{lin},`$ (6)
where $`\stackrel{}{E}`$ and $`\stackrel{}{H}`$ are the electric and magnetic intensity fields, $`\stackrel{}{D}`$ and $`\stackrel{}{B}`$ are the electric and magnetic induction fields, $`\stackrel{}{P}^{lin}`$, $`\stackrel{}{P}_n\mathrm{}`$ are the linear and nonlinear polarization of the medium respectively and $`\stackrel{}{M}_{lin}`$ is the linear magnetic polarization. The magnetic polarization (magnetization) $`\stackrel{}{M}_{lin}`$ is written as the product of the linear magnetic susceptibility $`\eta ^{(1)}`$ and the magnetic field $`\stackrel{}{H}`$. The nonstationary linear electric polarization can be written as:

$`\stackrel{}{P}^{lin}={\displaystyle \underset{\mathrm{}}{\overset{t}{}}}\left(\delta (t\tau )+4\pi \chi ^{\left(1\right)}\left(t\tau \right)\right)\stackrel{}{E}(\tau ,x,y,z)𝑑\tau `$
$`={\displaystyle \underset{\mathrm{}}{\overset{t}{}}}\epsilon _0\left(t\tau \right)\stackrel{}{E}(\tau ,x,y,z)𝑑\tau ,`$ (7)
where $`\chi ^{(1)}`$ and $`\epsilon _0`$ are the linear electric susceptibility and the dielectric constant respectively. Similar expression describes the dependence of $`\stackrel{}{B}`$ on $`\stackrel{}{H}`$ in the case of nonstationary linear magnetic polarization :

$`\stackrel{}{B}={\displaystyle \underset{\mathrm{}}{\overset{t}{}}}\left(\delta (t\tau )+4\pi \eta ^{\left(1\right)}\left(t\tau \right)\right)\stackrel{}{H}(\tau ,x,y,z)𝑑\tau `$
$`={\displaystyle \underset{\mathrm{}}{\overset{t}{}}}\mu _0\left(t\tau \right)\stackrel{}{H}(\tau ,x,y,z)𝑑\tau ,`$ (8)
where $`\eta ^{(1)}`$ and $`\mu _0`$ are the linear magnetic susceptibility and magnetic permeability respectively. The magnetic susceptibility of the main part of the dielectrics ranges from $`10^610^4`$ and usually decreases with the increasing of the frequency. In the following, we will study such media with nonstationary cubic nonlinear polarization, where the nonlinear polarization in the case of one carrying frequency can be expressed as:

$`\stackrel{}{P}_{nl}^{(3)}={\displaystyle \frac{3}{4}}{\displaystyle \underset{\mathrm{}}{\overset{t}{}}}{\displaystyle \underset{\mathrm{}}{\overset{t}{}}}{\displaystyle \underset{\mathrm{}}{\overset{t}{}}}\chi ^{\left(3\right)}(t\tau _1,t\tau _2,t\tau _3)`$
$`\times \stackrel{}{E}(\tau _1,x,y,z)\stackrel{}{E}^{}(\tau _2,x,y,z)\stackrel{}{E}(\tau _3,x,y,z)d\tau _1d\tau _2d\tau _3.`$ (9)
The causality request the next conditions on the response functions:

$`\epsilon (t\tau )=0;\chi ^{\left(3\right)}(t\tau _1,t\tau _2,t\tau _3)=0,`$
$`t\tau <0;t\tau _i<0;i=1,2,3.`$ (10)
That why we can prolonged the upper integral boundary to infinity and to use standard Fourier presentation :

$`{\displaystyle \underset{\mathrm{}}{\overset{t}{}}}\epsilon _0(\tau t)\mathrm{exp}(i\omega \tau )𝑑\tau ={\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}\epsilon _0(\tau t)\mathrm{exp}(i\omega \tau )𝑑\tau `$ (11)
$`{\displaystyle \underset{\mathrm{}}{\overset{t}{}}}{\displaystyle \underset{\mathrm{}}{\overset{t}{}}}{\displaystyle \underset{\mathrm{}}{\overset{t}{}}}\chi ^{\left(3\right)}(t\tau _1,t\tau _2,t\tau _3)𝑑\tau _1𝑑\tau _2𝑑\tau _3=`$
$`{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}\chi ^{\left(3\right)}(t\tau _1,t\tau _2,t\tau _3)𝑑\tau _1𝑑\tau _2𝑑\tau _3.`$ (12)
The spectral presentation of linear optical susceptibility $`\widehat{\epsilon _0}(\omega )`$ is connected to the nonstationary optical response function by the next Fourier transform:

$`\widehat{\epsilon }_0(\omega )\mathrm{exp}(i\omega t)={\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}\epsilon _0(t\tau )\mathrm{exp}(i\omega \tau )𝑑\tau ,`$ (13)
Similar expression for the spectral presentation of the non-stationary magnetic response $`\widehat{\mu _0}(\omega )`$:

$`\widehat{\mu }_0(\omega )\mathrm{exp}(i\omega t)={\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}\mu _0(t\tau )\mathrm{exp}(i\omega \tau )𝑑\tau .`$ (14)
and nonlinear optical susceptibility $`\chi ^{(3)}`$

$`\widehat{\chi }^{(3)}={\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}\chi ^{\left(3\right)}(t\tau _1,t\tau _2,t\tau _3)`$
$`\times \mathrm{exp}\left(i\left(\omega \left(\tau _1+\tau _2+\tau _3\right)\right)\right)d\tau _1d\tau _2d\tau _3`$ (15)
can be written. It is important to point here the remark of Akhmanov at all in , that such nonstationary representation applied to slowly varying amplitudes of electrical and magnet fields is valid as well when the optical pulse duration of the pulses $`t_0`$ is greater than the characteristic response time of the media $`\tau _0`$ ($`t_0>>\tau _0`$), as when the time duration of the pulses are equal or less than the time response of the media ($`t_0\tau _0`$). We will discuss this possibility in the process of deriving of the amplitude equations.

## 3 Deriving of the amplitude equations

In this section we derive the slowly varying amplitude approximation in the standard way, as it was done in , . For the case of Maxwell equations with linear and nonlinear dispersion (1)-(6) we define the electric and magnetic field amplitudes with the relations:

$`\stackrel{}{E}(x,y,z,t)=\stackrel{}{A}(x,y,z,t)\mathrm{exp}\left(i\left(\omega _0tg(x,y,z)\right)\right),`$ (16)
$`\stackrel{}{H}(x,y,z,t)=\stackrel{}{C}(x,y,z,t)\mathrm{exp}\left(i\left(\omega _0tq(x,y,z)\right)\right),`$ (17)
where $`\stackrel{}{A}`$, $`\stackrel{}{C}`$, $`\omega _0`$,$`g`$ and $`q`$ are the amplitudes of the electric and magnetic fields, the optical frequency and the real spatial phase functions respectively. In case of monochromatic and quasi-monochromatic fields the Stokes parameters can be constructed from transverse components of the wave field . This leads to two component vector fields in a plane, transverse to the direction of propagation. For electromagnetic fields with spectral bandwidth (our case) the two dimensional coherency tensor cannot be used and the Stokes parameters cannot be found directly. As it was shown by T.Carozzi and all in , using high order of symmetry (SU(3)), in this case six independent Stokes parameters can be found. This corresponds to three component vector field. Here is investigated this case. The increasing of the spectral bandwidth of the vector wave, increases also the depolarization term (component, normal to the standard Stokes coherent polarization plane).

We use the Fourier representation of the response functions (13)-(2) and of the amplitude functions $`\stackrel{}{A}`$ (and $`\stackrel{}{C}`$) to obtain next expressions for the first derivatives in time of the linear polarization, nonlinear polarization and magnetic induction fields:

$`{\displaystyle \frac{1}{c}}{\displaystyle \frac{\stackrel{}{P}^{lin}(x,y,z,t)}{t}}=i\mathrm{exp}\left(i\left(\omega _0tg(x,y,z)\right)\right)`$
$`\times {\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}{\displaystyle \frac{\omega \widehat{\epsilon }_0\left(\omega \right)}{c}}\stackrel{}{A}(x,y,z,\omega \omega _0)\mathrm{exp}\left(i(\omega \omega _0)t\right)d\omega ,`$ (18)
$`{\displaystyle \frac{4\pi }{c}}{\displaystyle \frac{\stackrel{}{P}^{nlin}(x,y,z,t)}{t}}=i\mathrm{exp}\left(i\left(\omega _0tg(x,y,z)\right)\right)`$
$`\times {\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}{\displaystyle \frac{3\pi \omega \widehat{\chi }^{(3)}\left(\omega \right)}{c}}|\stackrel{}{A}(x,y,z,\omega \omega _0)|^2\stackrel{}{A}(x,y,z,\omega \omega _0)`$
$`\times \mathrm{exp}\left(i(\omega \omega _0)t\right)d\omega .`$ (19)
$`{\displaystyle \frac{1}{c}}{\displaystyle \frac{\stackrel{}{B}(x,y,z,t)}{t}}=i\mathrm{exp}\left(i\left(\omega _0tq(x,y,z)\right)\right)`$
$`\times {\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}{\displaystyle \frac{\omega \widehat{\mu }_0\left(\omega \right)}{c}}\stackrel{}{C}(x,y,z,\omega \omega _0)\mathrm{exp}(i(\omega \omega _0)t)d\omega .`$ (20)
At this point we restrict the spectrum of the amplitude of electrical and magnet fields by writing the wave vectors $`k_{1,nl,2}`$ in a Taylor series:

$`k_1\left(\omega \right)={\displaystyle \frac{\omega \widehat{\epsilon }_0\left(\omega \right)}{c}}=k_1^0\left(\omega _0\right)+`$
$`+{\displaystyle \frac{\left(k_1\left(\omega _0\right)\right)}{\omega }}\left(\omega \omega _0\right)+{\displaystyle \frac{1}{2}}{\displaystyle \frac{^2\left(k_1\left(\omega _0\right)\right)}{\omega ^2}}\left(\omega \omega _0\right)^2+\mathrm{}.=`$
$`=k_1^0\left(\omega _0\right)+{\displaystyle \frac{1}{v_1}}\left(\omega \omega _0\right)+{\displaystyle \frac{1}{2}}k_1^\mathrm{"}\left(\omega \omega _0\right)^2+\mathrm{}.,`$ (21)
$`k_{nl}\left(\omega \right)={\displaystyle \frac{3\pi \omega \widehat{\chi }^{(3)}\left(\omega \right)}{c}}=`$
$`=k_{nl}^0\left(\omega _0\right)+{\displaystyle \frac{\left(k_{nl}\left(\omega _0\right)\right)}{\omega }}\left(\omega \omega _0\right)+\mathrm{}.=`$
$`=k_{nl}^0\left(\omega _0\right)+{\displaystyle \frac{1}{v_{nl}}}\left(\omega \omega _0\right)+\mathrm{}.,`$ (22)
$`k_2\left(\omega \right)={\displaystyle \frac{\omega \widehat{\mu }_0\left(\omega \right)}{c}}=k_2^0\left(\omega _0\right)`$
$`+{\displaystyle \frac{\left(k_2\left(\omega _0\right)\right)}{\omega }}\left(\omega \omega _0\right)+{\displaystyle \frac{1}{2}}{\displaystyle \frac{^2\left(k_2\left(\omega _0\right)\right)}{\omega ^2}}\left(\omega \omega _0\right)^2+\mathrm{}.=`$
$`k_2^0\left(\omega _0\right)+{\displaystyle \frac{1}{v_2}}\left(\omega \omega _0\right)+{\displaystyle \frac{1}{2}}k_2^\mathrm{"}\left(\omega \omega _0\right)^2+\mathrm{}.,`$ (23)
where $`v_i`$ and $`k_i^\mathrm{"};i=1,nl,2`$ have dimensions of group velocity, and nonlinear addition to the group velocity and dispersion respectively. The nonlinear wave vector is expressed:

$`k_{nl}^0=3\pi \omega _0\widehat{\chi }^{(3)}(\omega _0)/c={\displaystyle \frac{\omega _0\widehat{\epsilon }_0}{c}}{\displaystyle \frac{3\pi \widehat{\chi }^{(3)}(\omega _0)}{\widehat{\epsilon }_0}}=k_1n_2,`$ (24)
where

$`n_2(\omega _0)={\displaystyle \frac{3\pi \widehat{\chi }^{(3)}(\omega _0)}{\widehat{\epsilon }_0}},`$ (25)
is the nonlinear refractive index. It is important to note here that we do not use any approximation of the response function. There is only one requirement of the spectral presentations (2),(13) and (14) of the response functions: to admit first and second order derivatives in respect to frequency ( to be smooth functions). The restriction is only in respect of the relation between the main frequency $`\omega _0`$ and time duration of the envelope functions $`t_0`$ determinate from the relations (21), (22) and (23) (conditions for slowly varying amplitudes). Putting eqn. (21) in (3), (23) in (3) and (22) in (3), and keeping in mind the expressions for time derivatives of the spectral presentation of the amplitude functions, the first derivatives (3), (3) and (3) are presented in form:

$`{\displaystyle \frac{1}{c}}{\displaystyle \frac{\stackrel{}{P}^{lin}(x,y,z,t)}{t}}=\left(ik_1^0\stackrel{}{A}+{\displaystyle \frac{1}{v_1}}{\displaystyle \frac{\stackrel{}{A}}{t}}{\displaystyle \frac{ik_1^\mathrm{"}}{2}}{\displaystyle \frac{^2\stackrel{}{A}}{t^2}}\right)`$
$`\times \mathrm{exp}\left(i\left(\omega _0tg(x,y,z)\right)\right),`$ (26)
$`{\displaystyle \frac{4\pi }{c}}{\displaystyle \frac{\stackrel{}{P}^{nlin}(x,y,z,t)}{t}}=\left(ik_1^0n_2|\stackrel{}{A}|^2\stackrel{}{A}+{\displaystyle \frac{1}{v_{nl}}}{\displaystyle \frac{|\stackrel{}{A}|^2\stackrel{}{A}}{t}}\right)`$
$`\times \mathrm{exp}\left(i\left(\omega _0tg(x,y,z)\right)\right),`$ (27)
$`{\displaystyle \frac{1}{c}}{\displaystyle \frac{\stackrel{}{B}(x,y,z,t)}{t}}=\left(ik_2^0\stackrel{}{C}+{\displaystyle \frac{1}{v_2}}{\displaystyle \frac{\stackrel{}{C}}{t}}{\displaystyle \frac{ik_2^\mathrm{"}}{2}}{\displaystyle \frac{^2\stackrel{}{C}}{t^2}}\right)`$
$`\times \mathrm{exp}\left(i\left(\omega _0tq(x,y,z)\right)\right).`$ (28)
Finally, from the Maxwell equations (1)-(6), using eq. (3), eq. (3) and eq.(3), and using in fact that $`\stackrel{}{D}\stackrel{}{E}0`$ , we obtain the next system of Nonlinear Maxwell vector amplitude Equations (NME):

$`\times \stackrel{}{A}=ik_2^0\stackrel{}{C}{\displaystyle \frac{1}{v_2}}{\displaystyle \frac{\stackrel{}{C}}{t}}{\displaystyle \frac{ik_2^\mathrm{"}}{2}}{\displaystyle \frac{^2\stackrel{}{C}}{t^2}},`$ (29)
$`\times \stackrel{}{C}=ik_1^0\stackrel{}{A}+{\displaystyle \frac{1}{v_1}}{\displaystyle \frac{\stackrel{}{A}}{t}}{\displaystyle \frac{ik_1^\mathrm{"}}{2}}{\displaystyle \frac{^2\stackrel{}{A}}{t^2}}+ik_1^0n_2\left(\stackrel{}{A}\stackrel{}{A}^{}\right)\stackrel{}{A}`$
$`+\left({\displaystyle \frac{n_2}{v_1}}+k_1{\displaystyle \frac{n_2}{\omega }}\right){\displaystyle \frac{\left(\stackrel{}{A}\stackrel{}{A}^{}\right)\stackrel{}{A}}{t}},`$ (30)
$`\stackrel{}{A}=0,`$ (31)
$`\stackrel{}{C}=0,`$ (32)
if the gradient of the spatial phase functions $`g`$ and $`q`$ satisfied the relations:

$`g\times \stackrel{}{A}=0,`$ (33)
and

$`q\times \stackrel{}{C}=0.`$ (34)
The phase functions whichever satisfied (33)-(34) are determinate in Section ($`10`$).

We investigate the case, when our vector fields are presented as a sum of circular and linear polarizing components:

$`\stackrel{}{A}=\stackrel{}{A}_{lin}+\stackrel{}{A}_{cir},`$ (35)
$`\stackrel{}{C}=\stackrel{}{C}_{lin}+\stackrel{}{C}_{cir}.`$ (36)
The nonlinear polarization admit different nonlinear refractive indexes in the case of linear and circular polarization  ( $`n_2^{lin}n_2^{cir}`$). We will include this difference in the our rescaled equations, defining the next rescaled dependant variables:

$`\stackrel{}{A}=A_0^{lin}\stackrel{}{A}_{lin}^{^{}}+A_0^{cir}\stackrel{}{A}_{cir}^{^{}},`$ (37)
$`\stackrel{}{C}=C_0(\stackrel{}{C}_{lin}^{^{}}+\stackrel{}{C}_{cir}^{^{}}),`$ (38)
$`(A_0^{lin})^2={\displaystyle \frac{n_2^{cir}}{n_2^{lin}}}(A_0^{cir})^2,`$ (39)
and independent variables;

$`x=r_0x^{^{}};y=r_0y^{^{}};z=r_0z^{^{}},t=t_0t^{^{}}.`$ (40)
In additions the next constants are determinate:

$`\begin{array}{c}\alpha _i=k_i^0r_0;\beta _i=k_i^\mathrm{"}r_0/2t_0^2;\gamma _1=r_0k_1n_2^{cir}\left|A_0\right|^2;\hfill \\ \gamma _2=n_2^{cir}\left|A_0\right|^2;\gamma _3=v_1n_2^{cir}\left|A_0\right|^2/c;\hfill \\ \delta =v_1/v_21;v_1r_0/t_0;i=1,2.\hfill \end{array}`$ (42)
The NME (29)-(32)in rescaled variables are transformed into the following (the primes have been omitted for clarity):

$`\times \stackrel{}{A}=i\alpha _2\stackrel{}{C}\delta {\displaystyle \frac{\stackrel{}{C}}{t}}i\beta _2{\displaystyle \frac{^2\stackrel{}{C}}{t^2}},`$ (43)
$`\times \stackrel{}{C}=i\alpha _1\stackrel{}{A}+{\displaystyle \frac{\stackrel{}{A}}{t}}i\beta _1{\displaystyle \frac{^2\stackrel{}{A}}{t^2}}+i\gamma _1\left(\stackrel{}{A}\stackrel{}{A}^{}\right)\stackrel{}{A}`$
$`+\left(\gamma _2+\gamma _3\right){\displaystyle \frac{\left(\left(\stackrel{}{A}\stackrel{}{A}^{}\right)\stackrel{}{A}\right)}{t}},`$ (44)
$`\stackrel{}{A}=0,`$ (45)
$`\stackrel{}{C}=0.`$ (46)
We consider the case of slowly varying amplitudes approximation when the nonlinear constant is $`\gamma _1=1`$. Then the constants in fronts of last therm on right hand are $`\gamma _2\gamma _310^210^3`$. While the constants $`\alpha _i`$ have typical values $`\alpha 10^210^3\left(\alpha _ir_0k_i\right)`$, the dispersion constants $`\beta _i`$ have very small values $`\beta _110^510^6`$ for picosecond and sub-picosecond pulses in the transparency region of nonlinear optical media. Neglecting the small dispersion terms $`\beta _i<<1`$ and the last nonlinear term for witch $`\gamma _2\gamma _3<<1`$, the NME system (43)-(46)can be rewritten:

$`\times \stackrel{}{A}=i\alpha _2\stackrel{}{C}\delta {\displaystyle \frac{\stackrel{}{C}}{t}},`$ (47)
$`\times \stackrel{}{C}=i\alpha _1\stackrel{}{A}+{\displaystyle \frac{\stackrel{}{A}}{t}}+i\gamma _1\left(\stackrel{}{A}\stackrel{}{A}^{}\right)\stackrel{}{A},`$ (48)
$`\stackrel{}{A}=0,`$ (49)
$`\stackrel{}{C}=0.`$ (50)
## 4 Dirac representation of NME

To solve the NME (47)-(50), we apply the method of separation of variables. The slowly varying amplitude vector of the electric field $`\stackrel{}{A}`$ and the magnetic field $`\stackrel{}{C}`$ are represented as:

$`\stackrel{}{A}(x,y,z,t)=\stackrel{}{F}(x,y,z)\mathrm{exp}\left(i\mathrm{\Delta }\alpha t\right),`$ (51)
$`\stackrel{}{C}(x,y,z,t)=\stackrel{}{G}(x,y,z)\mathrm{exp}\left(i\mathrm{\Delta }\alpha t\right).`$ (52)
Substituting these forms into the NME (47)-(50) we obtain:

$`\times \stackrel{}{F}=i\nu _2\stackrel{}{G},`$ (53)
$`\times \stackrel{}{G}=i\nu _1\stackrel{}{F}+i\gamma _1\left(\stackrel{}{F}\stackrel{}{F}^{}\right)\stackrel{}{F},`$ (54)
$$\stackrel{}{F}=0,$$
(55)
$$\stackrel{}{G}=0,$$
(56)
where $`\nu _1=\alpha _1+\mathrm{\Delta }\alpha `$; $`\nu _2=\delta \mathrm{\Delta }\alpha \alpha _2>0`$. In a Cartesian coordinate system, the vector equations (53)-(56) are reduced to a scalar system of eight nonlinear wave equations. When the electric and magnetic fields are represented as a sum of a linear polarization component and a circular polarized one it is possible to reduce eqns. (53) - (56) to a system of four nonlinear equations. Substituting :

$`\psi _1=iF_l=iF_z,`$ (57)
$`\psi _2=F_c=iF_xF_y,`$ (58)
$`\psi _3=G_l=G_z,`$ (59)
$`\psi _4=G_c=G_xiG_y,`$ (60)
into the nonlinear system (53)-(56), we obtain a stationary nonlinear Dirac system of equations (NDE):

$`\left({\displaystyle \frac{}{x}}i{\displaystyle \frac{}{y}}\right)\mathrm{\Psi }_4+{\displaystyle \frac{}{z}}\mathrm{\Psi }_3=i\left(\nu _1+\gamma _1{\displaystyle \underset{i=1}{\overset{2}{}}}\left|\mathrm{\Psi }_i\right|^2\right)\mathrm{\Psi }_1,`$ (61)
$`\left({\displaystyle \frac{}{x}}+i{\displaystyle \frac{}{y}}\right)\mathrm{\Psi }_3{\displaystyle \frac{}{z}}\mathrm{\Psi }_4=i\left(\nu _1+\gamma _1{\displaystyle \underset{i=1}{\overset{2}{}}}\left|\mathrm{\Psi }_i\right|^2\right)\mathrm{\Psi }_2,`$ (62)
$`\left({\displaystyle \frac{}{x}}i{\displaystyle \frac{}{y}}\right)\mathrm{\Psi }_2+{\displaystyle \frac{}{z}}\mathrm{\Psi }_1=i\nu _2\mathrm{\Psi }_3,`$ (63)
$`\left({\displaystyle \frac{}{x}}+i{\displaystyle \frac{}{y}}\right)\mathrm{\Psi }_1{\displaystyle \frac{}{z}}\mathrm{\Psi }_2=i\nu _2\mathrm{\Psi }_4.`$ (64)
This substitution allows to reduce the system of eight equations (53)-(56) to a system of four scalar complex equations (61)-(64). The system (61)-(64) is the optical analog of the nonlinear Dirac equations (NDE). Note that the optical NDE are significantly different from the NDE in the field theory. The nonlinear part appears only in the first two coupled equations of the system.

## 5  Vortex solutions with orbital momentum l=1 and spin j=1/2

The symmetries of the NDE (61)-(64) are used to obtain exact vortex solutions. The NDE (61)-(64) have both spherical and spinor symmetry only in the case where the nonlinear part does not manifest the angular dependence on the radial variable $`\underset{i=1}{\overset{2}{}}\left|\mathrm{\Psi }_i(r,\theta ,\phi )\right|^2=F\left(r\right)`$. This type of solution can be found using with the following technique. Using Pauli matrices, we write the NDE system (61)-(64) as:

$`\left(\stackrel{}{\sigma }\stackrel{}{P}\right)\varphi =i\left(\nu _1+\gamma _1{\displaystyle \underset{i=1}{\overset{2}{}}}\left|\eta _i\right|^2\right)\eta ,`$ (65)
$`\left(\stackrel{}{\sigma }\stackrel{}{P}\right)\eta =i\nu _2\varphi ,`$ (66)
where

$`\stackrel{}{\sigma }=[\left(\begin{array}{cccccccccccccccccccc}0& 1& & & & & & & & & & & & & & & & & & \\ 1& 0& & & & & & & & & & & & & & & & & & \end{array}\right),\left(\begin{array}{cccccccccccccccccccc}0& i& & & & & & & & & & & & & & & & & & \\ i& 0& & & & & & & & & & & & & & & & & & \end{array}\right),\left(\begin{array}{cccccccccccccccccccc}1& 0& & & & & & & & & & & & & & & & & & \\ 0& 1& & & & & & & & & & & & & & & & & & \end{array}\right)]`$,

are the Pauli matrices, $`\stackrel{}{P}=(\frac{}{x},\frac{}{y},\frac{}{z})`$

is the differential operator and

$`\eta =\left(\begin{array}{cccccccccccccccccccc}\mathrm{\Psi }_1& & & & & & & & & & & & & & & & & & & \\ \mathrm{\Psi }_2& & & & & & & & & & & & & & & & & & & \end{array}\right)`$; $`\varphi =\left(\begin{array}{cccccccccccccccccccc}\mathrm{\Psi }_3& & & & & & & & & & & & & & & & & & & \\ \mathrm{\Psi }_4& & & & & & & & & & & & & & & & & & & \end{array}\right)`$,

are the corresponding spinors. After substituting eqn.(66) into eqn.(65) we obtain:

$$\left(\stackrel{}{\sigma }\stackrel{}{P}\right)\left(\stackrel{}{\sigma }\stackrel{}{P}\right)\eta =\nu _2\left(\nu _1+\gamma _1\underset{i=1}{\overset{2}{}}\left|\eta _i\right|^2\right)\eta .$$
(67)
When there is no external electric or magnetic field, the operator on the left-hand side of eqn. (67) is is the Laplacian operator $`\mathrm{\Delta }`$:

$$P^2=\left(\stackrel{}{\sigma }\stackrel{}{P}\right)\left(\stackrel{}{\sigma }\stackrel{}{P}\right)=\mathrm{\Delta }.$$
(68)
From (67) and (68) we obtain:

$$\nu _2\nu _1\eta +\nu _2\gamma _1\underset{i=1}{\overset{2}{}}\left|\eta _i\right|^2\eta +\mathrm{\Delta }\eta =0.$$
(69)
In case of spherical representation of the spinor equations (67), there are two possibilities, $`l=0`$ and $`l=1`$, that will permit only a radial dependence of the nonlinear part:

$$\underset{i=1}{\overset{2}{}}\left|\eta _i(r,\theta ,\phi )\right|^2=F\left(r\right).$$
(70)
For the $`l=0`$ case, the equations (69) are transformed to equations with radial parts:

$$\nu _2\nu _1\eta +\nu _2\gamma _1\underset{i=1}{\overset{2}{}}\left|\eta _i\right|^2\eta +\frac{^2\eta }{r^2}+\frac{2}{r}\frac{\eta }{r}=0.$$
(71)
The scalar variant of these equations has been investigated in many papers but exact localized solutions have not been found. In the case $`l=1`$ we look for spinors in the next two forms:

$`\eta =\left(\begin{array}{cccccccccccccccccccc}\stackrel{~}{\eta }\left(r\right)\mathrm{cos}\left(\theta \right)& & & & & & & & & & & & & & & & & & & \\ \stackrel{~}{\eta }\left(r\right)\mathrm{sin}\left(\theta \right)\mathrm{exp}\left(i\phi \right)& & & & & & & & & & & & & & & & & & & \end{array}\right).`$ (74)
and

$`\eta =\left(\begin{array}{cccccccccccccccccccc}\stackrel{~}{\eta }\left(r\right)\mathrm{sin}\left(\theta \right)\mathrm{exp}\left(i\phi \right)& & & & & & & & & & & & & & & & & & & \\ \stackrel{~}{\eta }\left(r\right)\mathrm{cos}\left(\theta \right)& & & & & & & & & & & & & & & & & & & \end{array}\right).`$ (77)
As it will be seen later this corresponds to two opposite directions of the own orbital momentum $`j=\pm 1/2`$ (opposite charge). After substituting solutions (74)-(77) into equations (69) the following equation describing the radial dependence is obtained:

$$\nu _2\nu _1\stackrel{~}{\eta }+\nu _2\gamma _1\left|\stackrel{~}{\eta }\right|^2\stackrel{~}{\eta }+\frac{^2\stackrel{~}{\eta }}{r^2}+\frac{2}{r}\frac{\stackrel{~}{\eta }}{r}\frac{2}{r^2}\stackrel{~}{\eta }=0.$$
(78)
The angular parts are the standard spherical harmonics with $`l=1`$. This system has exact vortex de Broglie soliton solutions  in the form:

$$\stackrel{~}{\eta }\left(r\right)=\frac{\sqrt{2}}{i}\frac{\mathrm{exp}\left(i\sqrt{\nu _1\nu _2}r\right)}{r},$$
(79)
if $`\nu _2\gamma _1=1`$. The complete solutions for these two cases are written:

$`\eta =\left(\begin{array}{cccccccccccccccccccc}\frac{\sqrt{2}}{i}\frac{\mathrm{exp}(i\sqrt{\nu _1\nu _2}r)}{r}\mathrm{cos}\left(\theta \right)& & & & & & & & & & & & & & & & & & & \\ & & & & & & & & & & & & & & & & & & & \\ \frac{\sqrt{2}}{i}\frac{\mathrm{exp}(i\sqrt{\nu _1\nu _2}r)}{r}\mathrm{sin}\left(\theta \right)\mathrm{exp}\left(i\phi \right)& & & & & & & & & & & & & & & & & & & \end{array}\right),`$ (83)
$`\varphi =\left(\begin{array}{cccccccccccccccccccc}\frac{\sqrt{2}}{i\nu _2}\left(\frac{i\sqrt{\nu _1\nu _2}\mathrm{exp}(i\sqrt{\nu _1\nu _2}r)}{r}+\frac{\mathrm{exp}i(\sqrt{\nu _1\nu _2}r)}{r^2}\right)& & & & & & & & & & & & & & & & & & & \\ & & & & & & & & & & & & & & & & & & & \\ 0& & & & & & & & & & & & & & & & & & & \end{array}\right),`$ (87)
and

$`\eta =\left(\begin{array}{cccccccccccccccccccc}\frac{\sqrt{2}}{i}\frac{\mathrm{exp}(i\sqrt{\nu _1\nu _2}r)}{r}\mathrm{sin}\left(\theta \right)\mathrm{exp}\left(i\phi \right)& & & & & & & & & & & & & & & & & & & \\ & & & & & & & & & & & & & & & & & & & \\ \frac{\sqrt{2}}{i}\frac{\mathrm{exp}(i\sqrt{\nu _1\nu _2}r)}{r}\mathrm{cos}\left(\theta \right)& & & & & & & & & & & & & & & & & & & \end{array}\right),`$ (91)
$`\varphi =\left(\begin{array}{cccccccccccccccccccc}0& & & & & & & & & & & & & & & & & & & \\ & & & & & & & & & & & & & & & & & & & \\ \frac{\sqrt{2}}{i\nu _2}\left(\frac{i\sqrt{\nu _1\nu _2}\mathrm{exp}(i\sqrt{\nu _1\nu _2}r)}{r}+\frac{\mathrm{exp}i(\sqrt{\nu _1\nu _2}r)}{r^2}\right)& & & & & & & & & & & & & & & & & & & \end{array}\right).`$ (95)
There is another, more direct way for separating the variables of the spinor equations (65)-(66). To illustrate this, we represent the NDE (61)-(64) in spherical variables:

$`\mathrm{exp}(i\phi )\left(\mathrm{sin}\theta {\displaystyle \frac{}{r}}+{\displaystyle \frac{\mathrm{cos}\theta }{r}}{\displaystyle \frac{}{\theta }}{\displaystyle \frac{i}{r\mathrm{sin}\theta }}{\displaystyle \frac{}{\phi }}\right)\mathrm{\Psi }_4`$
$`+\left(\mathrm{cos}\theta {\displaystyle \frac{}{r}}{\displaystyle \frac{\mathrm{sin}\theta }{r}}{\displaystyle \frac{}{\theta }}\right)\mathrm{\Psi }_3=i\left(\nu _1+\gamma _1{\displaystyle \underset{i=1}{\overset{2}{}}}\left|\mathrm{\Psi }_i\right|^2\right)\mathrm{\Psi }_1,`$ (96)
$`\mathrm{exp}(i\phi )\left(\mathrm{sin}\theta {\displaystyle \frac{}{r}}+{\displaystyle \frac{\mathrm{cos}\theta }{r}}{\displaystyle \frac{}{\theta }}+{\displaystyle \frac{i}{r\mathrm{sin}\theta }}{\displaystyle \frac{}{\phi }}\right)\mathrm{\Psi }_3`$
$`\left(\mathrm{cos}\theta {\displaystyle \frac{}{r}}{\displaystyle \frac{\mathrm{sin}\theta }{r}}{\displaystyle \frac{}{\theta }}\right)\mathrm{\Psi }_4=i\left(\nu _1+\gamma _1{\displaystyle \underset{i=1}{\overset{2}{}}}\left|\mathrm{\Psi }_i\right|^2\right)\mathrm{\Psi }_2,`$ (97)
$`\mathrm{exp}(i\phi )\left(\mathrm{sin}\theta {\displaystyle \frac{}{r}}+{\displaystyle \frac{\mathrm{cos}\theta }{r}}{\displaystyle \frac{}{\theta }}{\displaystyle \frac{i}{r\mathrm{sin}\theta }}{\displaystyle \frac{}{\phi }}\right)\mathrm{\Psi }_2`$
$`+\left(\mathrm{cos}\theta {\displaystyle \frac{}{r}}{\displaystyle \frac{\mathrm{sin}\theta }{r}}{\displaystyle \frac{}{\theta }}\right)\mathrm{\Psi }_1=i\nu _2\mathrm{\Psi }_3,`$ (98)
$`\mathrm{exp}(i\phi )\left(\mathrm{sin}\theta {\displaystyle \frac{}{r}}+{\displaystyle \frac{\mathrm{cos}\theta }{r}}{\displaystyle \frac{}{\theta }}+{\displaystyle \frac{i}{r\mathrm{sin}\theta }}{\displaystyle \frac{}{\phi }}\right)\mathrm{\Psi }_1`$
$`\left(\mathrm{cos}\theta {\displaystyle \frac{}{r}}{\displaystyle \frac{\mathrm{sin}\theta }{r}}{\displaystyle \frac{}{\theta }}\right)\mathrm{\Psi }_2=i\nu _2\mathrm{\Psi }_4.`$ (99)
We make the following two ansatzes for solutions to the system of nonlinear equations (5)-(5):

$`\begin{array}{c}\mathrm{\Psi }_1=a\left(r\right)\mathrm{cos}\left(\theta \right)\hfill \\ \mathrm{\Psi }_2=a\left(r\right)\mathrm{sin}\left(\theta \right)e^{i\phi }\hfill \\ \mathrm{\Psi }_3=ib\left(r\right)\hfill \\ \mathrm{\Psi }_4=0.\hfill \end{array}`$ (104)
and

$`\begin{array}{c}\mathrm{\Psi }_1=a\left(r\right)\mathrm{sin}\left(\theta \right)e^{i\phi }\hfill \\ \mathrm{\Psi }_2=a\left(r\right)\mathrm{cos}\left(\theta \right)\hfill \\ \mathrm{\Psi }_3=0\hfill \\ \mathrm{\Psi }_4=ib\left(r\right).\hfill \end{array}`$ (109)
Put (104) or (109) in equations (5)-(5) we separate the variables. The following system of equations describing the radial dependence of the amplitudes are obtained:

$$\frac{a\left(r\right)}{r}+\frac{2}{r}a\left(r\right)=\nu _2b\left(r\right),$$
(110)
$$\frac{b\left(r\right)}{r}=\nu _1a\left(r\right)+\gamma _1\left|a\left(r\right)\right|^2a\left(r\right).$$
(111)
In the system of equations (110)-(111) the nonlinearity appears only in the radial part, while for the angular part we have standard spherical spinors with spin $`j=\pm 1`$. Solving eqns.(110) and (111) we straightforwardly show that when $`\nu _2\gamma =1`$ these equations admit the localized de Broglie solitons of eqns. (83)-(87).

## 6  Hamiltonian representation of the NDE. First integrals for vortex solutions with spin $`j=\pm 1/2`$

It is not difficult to show that for the NDE system of eqns. (61)- (64) a Hamiltonian of the form:

$$H=\left(\stackrel{}{\sigma }\stackrel{}{P}\right)+\underset{i=1}{\overset{2}{}}\left|\mathrm{\Psi }_i\right|^2,$$
(112)
can be written. Using this, eqns. (61)-(64) can be rewritten in the form:

$$H\mathrm{\Psi }=\epsilon \mathrm{\Psi },$$
(113)
where $`\epsilon =(i\nu _1,i\nu _1,i\nu _2,i\nu _2)`$ is the energy operator.Here we investigate the case where the nonlinear part of the equation is represented as a number of spinors with a scalar sum that depends only on the radial component.

$$\underset{i=1}{\overset{2}{}}\left|\mathrm{\Psi }_i(r,\theta ,\phi )\right|^2=F\left(r\right).$$
(114)
We also introduce here the well known orbital momentum operator $`\stackrel{}{L}`$, own orbital (spin) momentum $`\stackrel{}{S}`$, and the full momentum $`\stackrel{}{J}`$, as well as:

$$\stackrel{}{L}=\stackrel{}{r}\times \stackrel{}{P}$$
(115)
$`{\displaystyle \frac{1}{2}}\stackrel{}{S}={\displaystyle \frac{1}{2}}\left(\begin{array}{cccccccccccccccccccc}\stackrel{}{\sigma }& 0& & & & & & & & & & & & & & & & & & \\ 0& \stackrel{}{\sigma }& & & & & & & & & & & & & & & & & & \end{array}\right)`$ (118)
$$\stackrel{}{J}=\stackrel{}{L}+\frac{1}{2}\stackrel{}{S}$$
(119)
It straightforward to show that the Hamiltonian (112) of eqn. (113) commutes with the operators $`\stackrel{}{J}^2`$ and $`J_z`$ (the z- projections must be x or y). Using these symmetries and the condition that the nonlinearity is of Kerr type, we can solve the NDE equations (113) by a separation of variables technique. We look for solutions in the form:

$`\mathrm{\Psi }_1=a\left(r\right)\mathrm{\Omega }_{jlM}`$ (120)
$`\mathrm{\Psi }_2=a\left(r\right)\mathrm{\Omega }_{jlM}`$ (121)
$`\mathrm{\Psi }_3=ib\left(r\right)\mathrm{\Omega }_{jl^{^{}}M}`$ (122)
$`\mathrm{\Psi }_4=ib\left(r\right)\mathrm{\Omega }_{jl^{^{}}M},`$ (123)
where $`\mathrm{\Omega }_{jlm}`$ is the spherical spinor, $`l+l^{^{}}=1`$, and $`a\left(r\right)`$ and $`b\left(r\right)`$ are arbitrary radial functions. Using the symmetries of (113) and the fact, that the nonlinear parts depend on $`r`$, we separate variables and obtain the following system of equations for the radial part:

$`{\displaystyle \frac{a\left(r\right)}{r}}+{\displaystyle \frac{1+\chi }{r}}a\left(r\right)=\nu _2b\left(r\right)`$ (124)
$`{\displaystyle \frac{b\left(r\right)}{r}}+{\displaystyle \frac{1\chi }{r}}b\left(r\right)=\nu _1a\left(r\right)+\gamma \left|a\left(r\right)\right|^2a\left(r\right),`$ (125)
where

$`\chi =l\left(l+1\right)j\left(j+1\right)1/4.`$ (126)
Exclude $`b(r)`$ from the system (124)-(125), we obtained the next equation for $`a(r)`$:

$`\nu _1\nu _2a(r)+{\displaystyle \frac{^2a}{r^2}}+{\displaystyle \frac{2}{r}}{\displaystyle \frac{a}{r}}{\displaystyle \frac{(1+\chi )\chi }{r^2}}a+\nu _2\gamma \left|a\right|^2a=0.`$ (127)
Formally this equation admit exact ”de Broglie” soliton solutions for arbitrary number of $`\chi `$. But as we remember our solutions are limited by the conditions (114) the nonlinear part to depends only on the radial components. The condition (114) for a number $`\chi 1`$ can be fulfilled also for a higher number of field on different frequencies. This case including also the parametric processes. The case of one carrying frequencies correspond to localized solutions with $`\chi =1`$ and angular components $`l=1`$ and $`j=\pm 1/2`$. In this case the system (124)-(125) becomes:

$`{\displaystyle \frac{a\left(r\right)}{r}}+{\displaystyle \frac{2}{r}}a\left(r\right)=\nu _2b\left(r\right)`$ (128)
$`{\displaystyle \frac{b\left(r\right)}{r}}=\nu _1a\left(r\right)+\gamma \left|a\left(r\right)\right|^2a\left(r\right).`$ (129)
As was shown above, this system has exact radial solutions of the form (79)-(87):

## 7 Experimental conditions

There are some differences between the nonlinear conditions for localized solutions of the Vector Nonlinear Schrodinger Equation (VNLS) and localized conditions for the Nonlinear Maxwell Equations (NME). The nonlinear parameter for the VNLS is written :

$$\gamma _{vnls}=k_0^2r_0^2n_2\left|A_0\right|^2=1.$$
(130)
For localized solutions in optical region the constant $`\alpha ^2`$ ranges from:

$$\alpha ^2=k_0^2r_0^210^410^6,$$
(131)
which corresponds to a required nonlinear refractive index change of order of:

$$n_2\left|A_0\right|^210^410^6.$$
(132)
On the other hand, the nonlinear condition for localized solutions of the normalized NDE (47)-(50) is:

$$\gamma _{NME}=\nu _2\gamma _1(\delta \mathrm{\Delta }\alpha \alpha _2)k_1r_0n_2\left|A_0\right|^2=1.$$
(133)
For the typical value of the constant $`\nu _21`$ the next required nonlinear refractive index change appears:

$$n_2\left|A_0\right|^210^210^3.$$
(134)
Another difference between these two cases is that the solutions of the VNLS are comprised of linearly polarized components and that the dispersion of the nonlinear medium plays an important role. This leads to the fact that vortex solutions of the VNLS may be observed only in special dispersion regions of nonlinear media. The vortices of the NDE are a combination of linear and circular polarization components, and do not have this marked dependence on the dispersion and may be observed in the transparency region of a nonlinear media.

## 8 Vortex solutions in a nonlinear Kerr type media without dispersion

The case without linear dispersion of the electrical and magnetic susceptibility will correspond to $`\chi ^{(1)}=const`$ and $`\eta ^{(1)}=const`$ . We suppose now that the amplitude functions do not depend from time and will look for 3D+1 monochromatic electric and magnet field of kind:

$`\stackrel{}{E}(x,y,z,t)=\stackrel{}{M}(x,y,z)\mathrm{exp}\left(i\omega _0t\right),`$ (135)
$`\stackrel{}{H}(x,y,z,t)=\stackrel{}{N}(x,y,z)\mathrm{exp}\left(i\omega _0t\right),`$ (136)
where $`\stackrel{}{M}`$ $`\stackrel{}{N}`$ and $`\omega _0`$ are the amplitudes of the electric, magnetic fields and the optical frequency respectively.

Substituting the relations (135)-(135) one obtain the next amplitude equations in the case of dispersionless media:

$`\times \stackrel{}{M}=i\alpha _2\stackrel{}{N},`$ (137)
$$\times \stackrel{}{N}=i\alpha _1\stackrel{}{M}+\gamma \left(\stackrel{}{M}\stackrel{}{M}^{}\right)\stackrel{}{M},$$
(138)
$$\stackrel{}{M}=0$$
(139)
$$\stackrel{}{N}=0$$
(140)
Again the electric and magnetic fields are represented as a sum of a linear polarization component and a circular polarized one. Substituting:

$$\psi _1=iM_l=iM_z$$
(141)
$`\psi _2=M_c=iM_xM_y`$ (142)
$`\psi _3=N_l=iN_z`$ (143)
$$\psi _4=N_c=N_xiN_y,$$
(144)
into the nonlinear system (137)-(140), we obtain the same kind of stationary nonlinear Dirac system of equations (NDE) as (61)-(64) :

$`\left({\displaystyle \frac{}{x}}i{\displaystyle \frac{}{y}}\right)\mathrm{\Psi }_4+{\displaystyle \frac{}{z}}\mathrm{\Psi }_3=i\left(\kappa _1+\zeta {\displaystyle \underset{i=1}{\overset{2}{}}}\left|\mathrm{\Psi }_i\right|^2\right)\mathrm{\Psi }_1`$ (145)
$$\left(\frac{}{x}+i\frac{}{y}\right)\mathrm{\Psi }_3\frac{}{z}\mathrm{\Psi }_4=i\left(\kappa _1+\zeta \underset{i=1}{\overset{2}{}}\left|\mathrm{\Psi }_i\right|^2\right)\mathrm{\Psi }_2$$
(146)
$$\left(\frac{}{x}i\frac{}{y}\right)\mathrm{\Psi }_2+\frac{}{z}\mathrm{\Psi }_1=i\kappa _2\mathrm{\Psi }_3$$
(147)
$`\left({\displaystyle \frac{}{x}}+i{\displaystyle \frac{}{y}}\right)\mathrm{\Psi }_1{\displaystyle \frac{}{z}}\mathrm{\Psi }_2=i\kappa _2\mathrm{\Psi }_4,`$ (148)
but with different coefficients: $`\kappa _1=\frac{\omega _0\epsilon _0}{c}`$; $`\kappa _2=\frac{\omega _0\mu _0}{c}`$; $`\zeta =\frac{4\pi \omega _0\chi (3)}{c}`$. Naturally, the NDE (145)-(148) admit the same kind of solutions (83)-(87)

## 9 Conditions for finiteness of the energy of the vortex solutions

In this section we will investigate the localized energy in the both cases; with and without of linear dispersion. This corresponds to two kind of optical vortices; with and without spectral bandwidth. In a dielectric media without dispersion the expression for the linear part of energy density is:

$`W_{lin}={\displaystyle \frac{1}{8\pi }}\left(\epsilon \left|\stackrel{}{E}\right|^2+\mu \left|\stackrel{}{H}\right|^2\right),`$ (149)
where $`\epsilon `$ and $`\mu `$ are constants. Substituting the vortex solutions (83)-(87) in (149) and integrating over 3D space one obtain that quasi-monochromatic vortices in a dielectric media without dispersion admit infinite energy. Now we come to the case of slowly varying amplitude approximation and to calculation of energy of the vortex solutions of a media with linear electric and magnet dispersion. To prove the finiteness of energy of the vortex solutions (83)-(87) we start with the equations for averaged in time balance of energy density of electrical and magnet field :

$`{\displaystyle \frac{W}{t}}={\displaystyle \frac{1}{16\pi }}\left(\stackrel{}{E}{\displaystyle \frac{\stackrel{}{D}^{}}{t}}+\stackrel{}{E}^{}{\displaystyle \frac{\stackrel{}{D}}{t}}+\stackrel{}{B}{\displaystyle \frac{\stackrel{}{B}^{}}{t}}+\stackrel{}{B}^{}{\displaystyle \frac{\stackrel{}{B}}{t}}\right),`$ (150)
where $`\stackrel{}{D}=\stackrel{}{P}_{lin}+4\pi \stackrel{}{P}_{nlin}`$ is a sum of the linear induction and the nonlinear induction of the electrical field. The calculations of the averaged energy of the optical waves in a dispersive media are worked out considering the first order of slowly varying amplitude approximation of electrical induction (the same as in the NDE). The result comes to the old result of Brillouin 1921 for energy density of electrical field:

$`W_{lin}={\displaystyle \frac{1}{8\pi }}\left({\displaystyle \frac{\left(\omega \widehat{\epsilon }_0\right)}{\omega }}\left|\stackrel{}{A}\right|^2+{\displaystyle \frac{\left(\omega \widehat{\mu }_0\right)}{\omega }}\left|\stackrel{}{C}\right|^2\right).`$ (151)
The conditions of electric constant:

$`\widehat{\epsilon }_0>0;{\displaystyle \frac{\left(\omega \widehat{\epsilon }_0\right)}{\omega }}>0`$ (152)
is fulfilled in most of the dielectrics. The conditions of the magnetic constant can be present in two cases. In the first one:

$`\widehat{\mu }_0>0;{\displaystyle \frac{\left(\omega \widehat{\mu }_0\right)}{\omega }}>0.`$ (153)
This condition corresponds to the case, when there are not included any magnet relaxation processes. To include these processes is possible only for ultrashort optical pulses with time duration of dimension of local structure of one paramagnetic media. The real part of magnet susceptibility in this case is symmetric function in respect to the paramagnetic resonance and the derivative in respect to the frequency is always negative.

$`\widehat{\mu }_0>0;{\displaystyle \frac{\widehat{\mu }_0(\omega )}{\omega }}<0.`$ (154)
Under appropriate conditions this leads to:

$`\widehat{\mu }_0>0;{\displaystyle \frac{\left(\omega \widehat{\mu }_0(\omega )\right)}{\omega }}<0.`$ (155)
We can find spectral region and dispersion parameters when the next condition will be satisfied:

$`{\displaystyle \frac{\left(\omega \widehat{\mu }_0\right)}{\omega }}={\displaystyle \frac{\left(\omega \widehat{\epsilon }_0\right)}{\omega }}.`$ (156)
Thus, the linear (infinity) part of energy density is zero. The nonlinear part of averaged energy density is expressed in  and for the vortex solutions (83)-(87) becomes:

$`W_{nlin}=n_2\left|\mathrm{\Psi }\right|^2\mathrm{}(\mathrm{\Psi }^2)+\omega _0{\displaystyle \frac{n_2}{\omega _0}}\left|\mathrm{\Psi }\right|^2\mathrm{}(\mathrm{\Psi }^2).`$ (157)
These results give the conditions for the finiteness of energy of the vortices. Integrating $`W_{nlin}`$ in the 3D space we obtain a finite value proportional to the main frequency $`\omega _0`$.

## 10  Spatial phase functions, Poynting vector and flow of energy

The kind of the phase functions which satisfied (33)-(34) is determinate for vortex solutions (104) with spin $`j=1/2`$. Using again the relations between the spinors of NDE and the amplitude functions (57)-(60) we have:

$`F_x={\displaystyle \frac{\psi _2\psi _2^{}}{2i}};F_y={\displaystyle \frac{\psi _2+\psi _2^{}}{2}};F_z={\displaystyle \frac{\psi _1\psi _1^{}}{2i}},`$ (158)
$`G_x={\displaystyle \frac{\psi _4+\psi _4^{}}{2}};G_y={\displaystyle \frac{\psi _4\psi _4^{}}{2i}};G_z={\displaystyle \frac{\psi _3\psi _3^{}}{2i}}.`$ (159)
Substituting the solutions (104) with spin $`j=1/2`$ in (158)-(159) for arbitrary real $`a(r)`$ and $`b(r)`$ we obtain:

$`F_x=a(r){\displaystyle \frac{x}{r}};F_y=a(r){\displaystyle \frac{y}{r}};F_z=a(r){\displaystyle \frac{z}{r}},`$ (160)
$`G_x=0;G_y=0;G_z=b(r).`$ (161)
We rewrite again the conditions for the spatial phase functions:

$`g\times \stackrel{}{F}=0;q\times \stackrel{}{G}=0.`$ (162)
These relations for solutions of kind (160)-(161) are satisfied only when:

$`g=k_0rorg=k_0f(r),`$ (163)
and

$`q=k_0zorq=k_0f(z),`$ (164)
where $`k_0`$ is the carrying wave number. The spatial phase functions of kind $`g=k_0r`$ and $`q=k_0z`$ correspond to spectral limited pulses which satisfied additional relations:

$`\mathrm{}k\mathrm{}r=const.`$ (165)
The spatial phase functions of kind $`g=k_0f(r)`$ and $`q=k_0f(z)`$ correspond to phase modulated pulses and for them the relations (165) is not satisfied. The Poynting vector can be expressed by the amplitude functions of the electrical and magnet field:

$`\stackrel{}{S}=\stackrel{}{E}(x,y,z,t)\times \stackrel{}{H}(x,y,z,t)=`$
$`\mathrm{exp}\left(i\left(W(t)K(r)\right)\right)\stackrel{}{F}(x,y,z,t)\times \stackrel{}{G}(x,y,z,t)`$ (166)
where $`W`$ and $`K`$ are scalar phase functions. Substituting the solutions with spin $`j=1/2`$ in above expression we find that:

$`\stackrel{}{S}=\mathrm{exp}\left(i\left(W(t)K(r)\right)\right)(a(r)b(r){\displaystyle \frac{y}{r}};a(r)b(r){\displaystyle \frac{x}{r}};0).`$ (167)
We see that the Poynting vector $`\stackrel{}{S}`$ is one circulation vector for solutions with spin $`j=1/2`$ and its divergency is zero:

$`\stackrel{}{S}=0.`$ (168)
The relation (168) determine that the energy flow through arbitrary closed surface around our vortex solutions with spin $`j=1/2`$ is zero. The relation (167) show that flow of energy of our solutions circulate in $`x,y`$ plane. Now we can generalize the above results for solutions with spin $`j=1/2`$: The vortex solutions with spin $`j=1/2`$ without external fields are immovable and electromagnetic energy oscillated in $`x,y`$ plane. The electrical field oscillating spherically, in ’r’ direction, while the magnet field oscillating in z direction. In the same way was calculated the Poynting vector for solutions with spin $`j=1/2`$. For them we obtain that $`\stackrel{}{S}0`$ and we expect that they are not stable. One exact investigation of stability request investigation also the perturbation of the Poynting theorem and will discussed in one next paper.

## 11  Conclusion

In this paper we derive a set of Nonlinear Maxwell amplitude Equations (NME) for nonlinear optical media with and without dispersion of the electric and magnetic susceptibility. We have shown that in cases of linear and circularly polarized components of the electric and magnetic fields, the NME reduces to the Nonlinear Dirac system of equations (NDE). The equations are represented in a spinor form. Using the method of separation of variables, exact vortex solutions for both cases have been obtained. The optical vortex solutions admit classical orbital momentum $`l=1`$ and classical own momentum $`j=\pm 1/2`$. Here I would like to explain more abbout the differences between solutions, obtained by separation of variables of the usual linear Dirac equation with potential depending only on ’r’ (for example hydrogen atom) and solutions of the NDE. In case of linear Dirac equations with potential to higher order spherical spinors correspond higher order of radial spherical Bessel functions. On the other hand for NDE to higher order of spherical spinors ($`\mathrm{}=1,2,..`$) correspond higher number of field components and higher value of localized energy. For all radial solutions of NDE the zero spherical radial Bessel function $`\frac{\mathrm{sin}\alpha r}{r}`$ are valid. The optical vortices in a media without dispersion admit infinite energy integral. The energy integral of the vortex solutions is finite only in some special cases of paramagnetic media with suitable conditions on linear electric and magnet dispersion . Using the Poynting vector for solutions with spin $`j=1/2`$ we find that the energy flow through arbitrary closed surface around our vortex solutions is zero and the localized energy of our solutions circulate in $`x,y`$ plane. Other important result is, that the vortex solutions with spin $`j=1/2`$ without external fields are immovable. The electrical field in the vortices oscillating spherically, while the magnet field oscillating in z direction. The initial investigations on stability of these solutions show the following: While the vortices with spin $`j=1/2`$ are stable, the vortices with opposite spin (charge) $`j=1/2`$ are not.

All of the above results will be discussed later in relation with nonlinear field theory. However one direct approach to the plasma physics can be found. We consider the old idea presented by Gaponov and Miller in : Confinement of charged particles in potentials of high frequency electromagnetic fields. In our article are presented solutions inspired by this idea. The optical vortices are of such type of high frequency fields with zero intensity in the centrum. In a cold plasma it is possible to obtain trapping of ions in the field of vortices. Thus, conditions for ruled nuclear fusion can be reached.


text:
# Intrinsic Variability and Field Statistics for the Vela Pulsar: 3. Two-Component Fits and Detailed Assessment of Stochastic Growth Theory

## 1 Introduction

The variability of pulsar emissions, both from pulse to pulse at a given phase and from phase to phase within a pulse, has long been unexplained \[Hankins 1996\]. These variations include subpulses \[Drake & Craft 1968\], with durations of order $`1050\%`$ the width of the average profile and sometimes a steady drift in phase, and microstructures superposed on the subpulses (Craft, Commella & Drake 1968, Hankins 1996), which are concentrated bursts of emissions that sometimes appear quasi-periodic. The statistics of the variable fields, such as the probability distribution of fields, have not been characterized until recently (Cairns, Johnston & Das 2001, hereafter Paper I) and few constraints have been placed on pulsar emission mechanisms, despite decades of research. This paper is the third in a series that characterizes the field statistics of pulsar variability in detail, advances interpretations in terms of existing theories for wave growth in inhomogeneous plasmas, and places constraints on pulsar emission mechanisms and the source physics. The initial analysis of Paper I, the second in the series \[Paper II\] \[called paper II hereafter\], and this paper address the Vela pulsar (PSR B8033-45), including both microstructure and subpulse effects simultaneously, while the fourth paper addresses PSRs B1641-45 and B0950+08 \[Cairns et al. 2002b\].

Wave growth in inhomogeneous media naturally results in bursty, variable waves. The field statistics are determined by the intrinsic statistics of radiation generated in an individual source region (predicted by theories for wave growth in inhomogeneous plasmas, which themselves depend on the emission mechanism and physics of the source plasma), effects of spatial variations due to possible superposition of emission from multiple (sub)sources, and scattering and other propagation effects between the observer and source, as reviewed in detail in Section 2 of Paper II. Only a very brief summary is given here, with most references provided in Paper II. First, the simplest assumption is made, that the observed field statistics are determined by the intrinsic statistics of the radiation process and propagation effects. Appeals to spatial variations are warranted only if the simplest assumptions fail to match data, which is not the case here. Second, scattering is expected to produce closely-Gaussian intensity statistics \[Ratcliffe 1956, Rickett 1977, Paper II\], defined by equation (7) of Paper II and hereafter referred to as equation (II.7). Third, stochastic growth theory \[SGT\] \[Robinson & Cairns 2001\], corresponding to a linear instability operating near marginal stability, results in the lognormal statistics defined by (II.5). Fourth, self-organized criticality \[SOC\] \[Bak et al. 1987\], which corresponds to a strongly coupled system that is driven away from marginal stability but relaxes back via both local and system-wide events, has power-law statistics defined by (II.6) with small indices. Fifth, for non-relativistic, weakly magnetized, electron-proton plasmas the nonlinear process of wave collapse \[Robinson 1997\], which involves the modulational (self-focusing) instability of a wavepacket, also results in power-law statistics but with higher indices that depend on wavepacket shape and effective dimensionality. Electromagnetic simulations for the relativistic, strongly magnetized electron-positron plasmas relevant to pulsars also show collapse occurring but the field statistics were not investigated (Weatherall 1997, 1998). Since collapse is qualitatively very similar in both sets of simulations, it is presumed below that collapse in pulsar environments also leads to power-law statistics with similar indices to those for electron-proton plasmas. Sixth, nonlinear decay processes, in which a primary wave decays into product waves, cutoff the field distribution with known functional form near and above the processes nonlinear threshold. Note that SGT can coexist at moderate fields with wave collapse or decay active at higher fields.

Paper II demonstrates multiple results for the Vela pulsar. (i) The observed probability distributions $`P(\mathrm{log}E)`$ of the wave field $`E`$ are well-defined as functions of the pulse phase $`\varphi `$. (ii) These distributions evolve with $`\varphi `$, from Gaussian intensity statistics off-pulse, to approximately power-law statistics in a transition region, to lognormal statistics near the pulse center, and thence back through approximately power-law statistics to Gaussian intensity statistics off-pulse again. Figure 1 illustrates this behaviour. (iii) The off-pulse data are quantitatively consistent with Gaussian intensity statistics. (iv) The field statistics near the pulse center are consistent with the lognormal form (II.5) provided analysis is restricted to fields well above the noise background. (v) However, at low $`E`$ the observed $`P(\mathrm{log}E)`$ distributions lie consistently above the lognormal fit. This suggests that a second population of waves is superposed on the lognormal component.

Figure 1 also shows the restricted phase range ($`430434`$) where Vela’s giant micropulses are observed \[Johnston et al. 2001, Kramer et al. 2002\]: defined by analogy with giant pulses \[Lundgren et al. 1995\], which have pulse-integrated fluxes exceeding $`10`$ times the average pulse-integrated flux, giant micropulses have fluxes greater than $`10`$ times the average flux at that phase. Giant micropulses and pulses are discussed further in Section 7.

Superposition of at least two wave populations is expected for the pulsar, since the Gaussian population observed off-pulse should extend on-pulse if it corresponds to measurement noise, sky background, or scattered radiation. Accordingly here we analyze the Vela data in detail in terms of two superposed wave populations, using a new theory for the statistics of two populations of transverse wave fields with known statistics that are vectorially combined together (Cairns, Robinson & Das 2002c). This theory is summarized in Section 2. The dataset is described briefly in Section 3. In Section 4 we demonstrate in detail that the approximately power-law $`P(\mathrm{log}E)`$ distributions observed in the transition regions are consistent with vector convolution of a Gaussian intensity distribution with a lognormal distribution. These fits consider almost all available data and typically have very high statistical significance. Section 5 addresses the lognormal region, showing that the combination of two lognormals provides very good fits to the available data, being demonstrably superior to the single-lognormal fits in paper I and to Gaussian-lognormal combinations. Smooth evolution in the fit parameters across the source is demonstrated in Section 6, providing additional confidence in the fits and fitting procedure. Moreover, this smooth evolution also provides an additional argument that the approximately power-law statistics seen at phases in the transition regions are not intrinsic. It is argued in Section 7 that these analyses demonstrate that the Vela pulsar’s field statistics are consistent with lognormal statistics, and so SGT, being relevant. The analyses are consistent with the pulsar’s emission mechanism being purely linear (either direct or indirect) in these phase ranges, with no evidence for a nonlinear emission mechanism, consistent with our earlier analyses (Papers I & II). The results are also set in a wider context, comparing them with those for the Vela pulsar’s giant micropulses and with solar system radio and plasma wave emissions. The paper’s conclusions are summarized in Section 7.

## 2 Statistics of Two Vectorially-Superposed Wave Populations

Suppose the observed field $`𝐄`$ is formed by vector addition of two fields $`𝐄_1`$ and $`𝐄_2`$ that are transverse to the measurement plane, result from different wave populations or source regions, and have separate probability distributions $`P_i(E_i^2)`$ for $`i=1,2`$. The probability distribution $`P(E^2)`$ is then defined in terms of the magnitudes $`E_1`$ and $`E_2`$ and the angle $`\theta `$ between $`𝐄_1`$ and $`𝐄_2`$ by \[Cairns et al. 2002c\]

$`P(E^2)`$ $`=`$ $`2\pi A{\displaystyle 𝑑E_1^2𝑑E_2^2𝑑\theta P(E_1^2)P(E_2^2)}`$ (1)
$`\times `$ $`P_\theta (\theta )\delta (E^2E_1^2E_2^22E_1E_2\mathrm{cos}\theta ).`$
Here $`P_\theta (\theta )`$ is the probability distribution for $`\theta `$, $`A`$ is a normalization constant, and all probability distributions are normalized by $`𝑑E_i^2P_i(E_i^2)=1=𝑑E^2P(E^2)`$. The delta function enforces the vector addition $`𝐄=𝐄_1+𝐄_2`$.

Assuming that $`\theta `$ is uniformly distributed, meaning that the two signals are produced independently in either the same source or different regions, $`P(\theta )=(2\pi )^1`$). Then, integrating over the $`\theta `$ integral using the delta function leads to \[Cairns et al. 2002c\]

$`P(E^2)`$ $`=`$ $`A{\displaystyle 𝑑E_1^2𝑑E_2^2P_1(E_1^2)P_2(E_2^2)}`$ (2)
$`\times `$ $`|[E^2(E_1+E_2)^2][E^2(E_1E_2)^2]|^{1/2}.`$
The square root in (2) contains integrable singularities corresponding to $`\mathrm{cos}\theta =\pm 1`$, where $`𝐄_1`$ and $`𝐄_2`$ are parallel or antiparallel. Figure 2 shows the integration domain for (2), limited by the constraint $`1\mathrm{cos}\theta 1`$ and the physical portions of the singularities ($`E`$, $`E_1`$, and $`E_2`$ must be positive-definite). Two of the three singularity lines in Figure 2 correspond to $`\mathrm{cos}\theta =1`$ and the factor $`E^2=(E_1E_2)^2`$ in (2), implying that antiparallel fields $`𝐄_1`$ and $`𝐄_2`$ significantly affect the distribution $`P(E^2)`$. The integral (2) is performed numerically.

The $`P(E^2)`$ distribution is related to the corresponding intensity distribution $`P(I)`$ and the logarithmically-binned distribution $`P(\mathrm{log}E)`$ through their differentials and normalization conditions by

$$P(E^2)=P(\mathrm{log}E)/(2E^2\mathrm{ln}10)=aP(I).$$
(3)
Explicit expressions for lognormal field distributions and for Gaussian intensity distributions are given in equations (II.5) and (II.7), respectively.

Detailed analyses and explanations of the $`P(E^2)`$ distributions that result from integrating (2) over the integration domain in Figure 2 are described elsewhere \[Cairns et al. 2002c\]. Only four of their important qualitative results are summarized here and illustrated using the $`P(E^2)`$ distribution predicted for two lognormal distributions in Figure 3.

First, typically, as in Figure 3, the combined $`P(E^2)`$ distribution, at large $`E`$ above its peak, has the functional form of the individual distribution $`P_i`$ that extends to higher $`E_i`$. Expected intuitively when one distribution is much more intense then the other, this result explains why previous single-component analyses were usually viable.

Second, a different functional form can be develop in the transition region between the two individual distributions. Specifically a nearly power-law form over $`12`$ decades can result from combining a Gaussian or lognormal distribution with a lognormal distribution that is centered at lower $`E`$ but extends to higher $`E`$ \[Cairns et al. 2002c\].

Third, the prediction (2) for vectorial convolution of $`𝐄_1`$ and $`𝐄_2`$ differs markedly from the corresponding predictions for convolution of the wave intensities or field amplitudes, as illustrated in Figure 3. Accordingly, despite their familiarity and frequent use \[Romani & Johnston 2001, Johnston & Romani 2002\], intensity and amplitude convolution should not be used to interpret the detailed field statistics.

Fourth, vectorial combination often produces a relatively flat tail at low $`E`$ below the peak of $`P(E^2)`$. In Figure 3 this tail clearly does not come from the individual lognormal distributions, which show very rapid falloffs, but is instead due to overlap between the $`P_i`$ distributions for antiparallel vectors: since combining antiparallel vectors $`𝐄_1`$ and $`𝐄_2`$ with similar magnitudes results in $`E0`$, the $`P(E^2)`$ distribution will be significant at low $`E`$ if the $`P_1`$ and $`P_2`$ distributions overlap significantly. These features are shown below to be directly relevant to the Vela data. Specifically, the power-law distributions discussed in Section 3 are due to the second effect, the first result applies primarily in Section 4, and the tail at low $`E`$ is important observationally for all the Vela data. Note, via (3), the relatively flat tail for $`P(E^2)`$ in Figure 3 becomes an $`E^2`$ power-law for the $`P(\mathrm{log}E)`$ distribution. Only $`P(\mathrm{log}E)`$ distributions are analyzed below.

## 3 Data Set and Previous Fits

The Vela dataset analyzed is described in detail in Paper II and elsewhere (Johnston et al. 2001, Kramer et al. 2002), so only a minimal description is given here. The radiation fields on the Parkes antenna are detected as voltages, processed via a standard backend system, and then converted into calibrated fluxes $`F`$ averaged over the detector’s 20 MHz bandwidth. These flux samples are recast in terms of fields and intensities to permit direct comparisons between the data and theories for wave statistics. The variables $`E^{}`$ and $`I^{}`$ defined by

$`E^{}=`$ $`(F+I_{off}^{})^{1/2}`$ $`E,`$ (4)
$`I^{}=`$ $`F+I_{off}^{}`$ $`I,`$ (5)
are related directly to the calibrated fluxes and are proxies for the incident field $`E`$ and intensity $`I`$ with units of (mJy)1/2 and mJy, respectively. The telescope’s backend system removes the time-steady levels of receiver noise, sky and supernova background, and pulsar emission. Accordingly, an offset $`I_{off}^{}=1250`$ mJy, equal approximately to the rms noise level off-pulse (Johnston et al. 2001, Papers I & II), is added to mitigate removal of the time-steady pulsar emission $`F_{psr}`$, as justified in detail in paper II. Plausible causes for non-zero $`F_{psr}`$ are (i) coherently produced radiation that undergoes scattering and diffusion as it propagates from its source, perhaps changing the phase bin at which it is observed and producing a significant background, as well as (ii) synchrotron emission from the pulsar magnetosphere or jets. Detailed arguments against this procedure and the backend system significantly modifying the true statistics are presented in Sections 3 and 9 of paper II. Note that information on both microstructure and subpulses is retained.

Figure 1 above, Figure 1 of Paper I, and Figures 1, 3–5, and 12 of Paper II show the observational context in detail. As indicated in Figure 1, and shown specifically in Figure 4 of paper II, the $`P(\mathrm{log}E^{})`$ distributions observed for phases $`415460`$ and $`545610`$ are approximately power-law at high $`E^{}`$. These distributions, in the so-called transition region, are fitted for the first time in section 4.

The $`P(\mathrm{log}E^{})`$ distributions observed in the approximate range $`460540`$ appear strongly lognormal by eye and previous single-component lognormal fits are good and have reasonable statistical significance if restricted to bins with $`10^2`$ counts and fields $`\stackrel{>}{}100`$ (mJy)1/2 (Papers I & II). Typically these fits accurately model the observed distribution well at fields close to and above the peak in $`P(\mathrm{log}E^{})`$, but lie significantly below the observations at lower $`E^{}`$. This latter point argues for a second component contributing to the observed distribution \[Paper II\], as shown in Section 4 for the power-law regions. Two-component double-lognormal fits for almost the entire range of $`E^{}`$ are presented for phases $`460540`$ in Section 5.

## 4 Two-Component Fits in the Transition Regions

Figure 4 shows the $`P(\mathrm{log}E^{})`$ distribution observed at phase 440, the best fit to (2) for the combination of a Gaussian intensity distribution and a lognormal distribution, and the two individual distributions. The best fit was obtained by minimizing $`\chi ^2`$ for all bins with $`10`$ counts using the amoeba algorithm \[Press et al. 1986\]. In comparison, previous fits of similar Vela data (Paper I and II) considered only bins with $`100`$ counts and fields $`E^{}\stackrel{>}{}100`$ mJy1/2 ($`I^{}\stackrel{>}{}10^4`$ mJy) and fitted only a single component.

Excellent quantitative agreement is evident in Figure 4 between the data and best-fit curve. The fit is very good statistically: Table 1 lists the fit parameters for the Gaussian component, $`I_0^{}`$ and $`\sigma _I^{}`$ in (II.7), and the lognormal component, $`\mu `$ and $`\sigma `$ in (II.5), as well as $`\chi ^2`$, the number of degrees of freedom $`N_{df}`$, and the significance probability of the fit, $`P(\chi ^2)`$. Note that good fits have $`\chi ^2N_{df}`$ and $`P(\chi ^2)\stackrel{>}{}1\%`$. Comparing the Gaussian fit parameters with the properties of the off-pulse (phases $`391400`$) Gaussian intensity noise fitted in paper II, $`I_0^{}=1215`$ mJy and $`\sigma _I^{}=1420`$ mJy, it is clear that the Gaussian component found here corresponds closely to the off-pulse Gaussian background and to the offset $`I_{off}^{}`$ introduced in (5) – (5).

The best-fit line in Figure 4 is very close to the underlying Gaussian distribution except for the $`4`$ data points with largest $`E^{}`$, where the data exceed the Gaussian predictions. This excess is due to the underlying lognormal distribution and corresponds to the highest $`E`$ signals from the lognormal adding to the Gaussian background. The simplest interpretation is that the lognormal distribution corresponds to the pulsar’s intrinsic emission while the Gaussian component is a background signal. Evidence for evolution of the Gaussian background with phase is presented in Section 6 below. Thus, for phase 440, the approximately power-law form below the peak of the observed $`P(\mathrm{log}E^{})`$ distribution corresponds to the Gaussian background while the upturn at high $`E^{}`$, which develops into a power-law at larger $`\varphi `$, is due to the lognormal component.

Figure 5 shows the observed $`P(\mathrm{log}E^{})`$ distribution and associated Gaussian-lognormal fits at phase 450, where the observed distribution has developed closely power-law forms at both low and high $`E^{}`$. For reference, $`P(\mathrm{log}E^{})E^{}^{4.1\pm 0.5}`$ at high $`E^{}`$ and $`E^{}^{+3.1\pm 0.3}`$ at low $`E^{}`$ and at first sight the distribution appears rather similar to theory and observations for wave collapse discussed elsewhere \[Robinson & Newman 1990, Robinson 1997\]. Theoretical difficulties in interpreting these data in terms of collapse are described in Section 7. In contrast, the best-fit Gaussian-lognormal combination for (2) agrees very well with the observed distribution, with high statistical significance. Comparing Figures 45 and the fit parameters in Table 1, it is clear that evolution in the lognormal component explains simply the development of power-law character in the $`P(\mathrm{log}E^{})`$ distributions at high $`E^{}`$ for these phases. The ensuing examples show that this interpretation is also consistent for phase bins $`<440`$ and in the other transition region (bins $`545610`$).

Fitting the Gaussian-lognormal combination in the phase range $`465540`$, where the $`P(\mathrm{log}E^{})`$ distributions look strongly lognormal as shown in Section 5 and in Papers I & II, typically leads to worse fits and lower statistical significance than a double-lognormal combination. This is shown primarily in the next 2 sections, with Figure 6 shown here mostly for consistency. While Figure 6 shows reasonable agreement at high $`E^{}`$, the fit clearly fails for both the tail and the peak.

Consider next the transition region in phase bins $`545610`$ where the $`P(\mathrm{log}E^{})`$ distributions range from being initially power-law at high $`E^{}`$ towards being Gaussian in $`I^{}`$. Figure 7 shows the observed $`P(\mathrm{log}E^{})`$ distribution and associated Gaussian-lognormal fit at phase 560, where the distribution appears power-law both above and below its peak (indices $`6.8\pm 0.9`$ and $`2.2\pm 0.4`$, respectively). Once again the agreement between the fit and data is very good, both quantitatively and in terms of statistical significance.

As the pulsar’s average intensity decreases for $`\varphi >560`$, and so the value $`\mu `$ of the inferred lognormal component decreases, the power-law feature at high $`E^{}`$ is expected to weaken and retreat into the background distribution by analogy with the results for phases 440 and 450. Figures 8 and 9 illustrate this trend for phase bins 590 and 610, respectively, and show that the observed $`P(\mathrm{log}E^{})`$ distributions remain very well fitted by the Gaussian-lognormal combination. The lognormal component indeed has $`\mu `$ decreasing with increasing phase (Table 1), while $`\sigma `$ remains approximately constant, showing that the lognormal moves towards and below the fields characteristic of the relatively constant Gaussian component. Similar results are obtained for phases 410 to 440, where the pulsar and lognormal component are moving above background.

Phase 590 (Figure 8) corresponds to the so-called “bump” region \[Johnston et al. 2001\] and shows that the bump region is not unusual, since the pulsar’s emission continues to be well described as a lognormal component. Comparison with Figure 8a of Kramer et al. (2002) emphasizes the importance of vectorial convolution of Gaussian and lognormal components via (2) for obtaining a good fit.

In summary, the fits in Figures 4 to 9 and the associated variations in Table 1’s fit parameters are strong observational evidence that the development of power-law-like features in the $`P(\mathrm{log}E^{})`$ distributions at high $`E^{}`$ in phase bins $`415460`$ and $`545610`$ correspond to combination of a Gaussian background component with an underlying lognormal component of the pulsar’s emission. These power-law-like features are thus transitional stages between Gaussian intensity statistics off-pulse and strongly lognormal statistics near the pulse peak, justifying the term “transition regions” adopted earlier \[Paper II\] for these phase domains. These data therefore provide no evidence of intrinsic power-law character for the pulsar’s field statistics. In Section 7 below, it is shown that the power-law exponents for Figures 5 and 7 are inconsistent with collapse theory at low $`E^{}`$ and at best marginally consistent at high $`E^{}`$. These analyses are consistent with the dominant pulsar emission on-pulse being represented well by lognormal statistics.

## 5 Lognormal Region and Two-Component Fits

For most phase bins in the range $`460540`$, fitting the observed $`P(\mathrm{log}E^{})`$ distributions with a Gaussian-lognormal combination leads to increased $`\chi ^2`$ and worse fits, when the fitting algorithm converges at all, than for phases below 460 and above 540. Pursuing alternatives, we found that double-lognormal fits work very well for phases 475 – 510, with monotonically decreasing statistical significance from phases 475 to 450 and from 520 to 550. Adequate double-lognormal fits were not found at phases below $`450`$ or above $`550`$.

Figure 10 shows that the double-lognormal combination works well for phase $`470`$, which is very near the peak of the average profile (Figure 1), fitting the observed distribution well except at very low and high $`E^{}`$. The fit parameters are listed in Table 2. Note that the two individual lognormal distributions have very different parameters, with one centered at high $`E^{}`$ and the other at low $`E^{}`$ corresponding to the centroid $`I_0^{}=E^2`$ in (II.7) of the off-pulse Gaussian component. This is why single-component fits restricted to high $`E^{}10^2`$ (mJy)1/2 worked well in previous analyses (Papers I & II). Similarly, at this phase fitting two components rather than one (ignoring only bins with $`<10`$ counts) decreases $`\chi ^2`$ from $`188`$ to $`119`$ but does not improve the fit’s statistical significance qualitatively.

The results for phase 480 are excellent (Figure 11), with very high statistical significance for the fit and good agreement at all values of $`\mathrm{log}E^{}`$. Comparing Figures 10 and 11, the range of the $`P(\mathrm{log}E^{})`$ distribution narrows significantly ($`\sigma `$ decreases) while its centroid ($`\mu `$) varies relatively little. This necessarily decreases the contribution of the off-pulse (Gaussian) component to the fit. Now $`\chi ^2`$ decreases from $`516`$ to $`5.4`$ for the one- and two-component fits, respectively, corresponding to a qualitative change in statistical significance. A similar result is found at phase 475.

Figure 12 presents the two-component fit for phase 490, showing it to agree very well with the data, both at large and small $`E^{}`$, with reasonable statistical significance. Comparing this with earlier single-component fits (Papers I & II), the two-component fit agrees very well with the data over a much larger domain of $`E^{}`$.

Figure 13 displays the double-lognormal fits for phase 510, which can be compared with the Gaussian-lognormal fit in Figure 6. The double-lognormal fits are clearly superior, correctly fitting the tail at low $`E^{}`$ and the peak of the distribution, which were missed by the Gaussian-lognormal fit. Tables 1 and 2 show the differing statistical significances of the fits. This case demonstrates that the nature of the individual distribution (e.g., Gaussian versus lognormal) is important and an optimum choice can be discerned from the fitting results, even when the difference is primarily important for only a few data points at low $`E^{}`$.

## 6 Evolution of the Fit Parameters With Phase

The Gaussian and lognormal fit parameters for the $`P(\mathrm{log}E^{})`$ distributions vary with the Vela pulsar’s phase (Figure 14), as is evident also from comparing Figures 4 – 14 and Tables 1–2. Values of $`\chi ^2`$ less than or approximately equal to the number of degrees of freedom $`N_{df}`$ correspond to good fits from a statistical point of view. A number of important results are apparent:

1. Either the Gaussian-lognormal or double-lognormal combinations fit the observed $`P(\mathrm{log}E^{})`$ distributions well at essentially all phases. Phases where only one symbol appears indicates that either no fit was found for the other combination or else the other fit had $`\chi ^2\stackrel{>}{}10^4`$.
2. Double-lognormal fits are superior for $`465\stackrel{<}{}\varphi \stackrel{<}{}540`$ where the pulsar’s average profile peaks, while the Gaussian-lognormal combination is superior outside this domain.
3. The fit parameters vary smoothly with phase, except near phase bins $`455`$ and $`555`$ where the best-fit combinations are changing from Gaussian-lognormal to double-lognormal or vice-versa. Variations in the lognormal parameters are very large, consistent with the source plasma’s parameters varying significantly with phase.
4. The values of $`\mu `$ and $`\sigma `$ for the Gaussian-lognormal fits agree very well with those for the primary lognormal ($`\mu _1`$ and $`\sigma _1`$) from the double-lognormal fits. This shows that a high degree of confidence can be attached both to the fitting procedures and to interpretation of the pulsar’s statistics and variability in terms of lognormal field distributions.
5. Since the primary lognormal parameters vary smoothly where the superior fitting function switches between Gaussian-lognormal and double-lognormal, these transitions primarily correspond to the functional form of the low $`E^{}`$ wave distribution. This transition can be interpreted as evidence for evolution in the statistics of the less intense wave distribution. However, while these transitions are statistically significant, as implied by Figures 4 \- 14, the interpretation is complicated by the fits yielding $`I_0^{}e^{2\mu _2}I_{off}^{}=1250`$ mJy (Figure 14). Here $`\mu _2`$ corresponds to the weaker lognormal for the double-lognormal fit. Any errors in $`I_{off}^{}`$ thus limit the quantitative (but not qualitative) significance of fit parameters $`I_0^{}`$ and $`\mu _2`$.
6. The Gaussian component off-pulse and in the outer transition regions is best interpreted as background or measurement noise, which should be present at all phases. This is consistent with single-component Gaussian fits at off-pulse phases \[Paper II\] and the finding that $`I_0^{}I_{off}^{}`$. There is weak evidence, however, that $`I_0^{}`$ increases from phase bins 430 to 455 and from bins 600 to 560, in the transition regions. This can be interpreted in terms of another Gaussian component dominating or adding to the off-pulse Gaussian component in these bins, perhaps being interpretable in terms of scattered radiation or another source of radiation with Gaussian statistics. Similar evidence for evolution of the Gaussian component across the source exists for other pulsars \[Cairns et al. 2002b\].
7. The transitions from Gaussian-lognormal to double-lognormal fits and vice versa can be interpreted as evidence for a third wave component appearing near the peak of Vela’s average profile and dominating the off-pulse Gaussian component and the evolving Gaussian component possibly present in the transition regions.
8. Finally, the bottom panel of Figure 14 shows that $`\sigma `$ is large near the pulse edges and small near the pulse center. This is consistent with earlier statements that pulsar modulation indices behave in the same way \[Taylor et al. 1975, Johnston & Romani 2002\]; however, the new results place the earlier statements on a much stronger physical basis, since they demonstrate that the Vela pulsar’s variability corresponds to lognormal statistics and that the lognormal parameter $`\sigma `$ varies in this way across the source, both significant steps towards developing a detailed physical model for the source plasma and emission physics. Such a model would link the wave parameters with instability physics and medium inhomogeneities, as done in several solar system contexts \[Robinson et al. 1993, Cairns & Menietti 2001\].

## 7 Discussion

The existence of very good Gaussian-lognormal or double-lognormal fits to the observed $`P(\mathrm{log}E^{})`$ distributions at almost all on-pulse phases (excluding the giant micropulses discussed below) and the smooth evolution of the fit parameters with phase constitute very strong evidence that lognormal statistics, consistent with SGT, are relevant at all on-pulse phases for the Vela pulsar. Put another way, the pulsar’s variability at all phases can be interpreted as the result of a pure SGT system coupled with a second group of weaker waves that have either Gaussian intensity statistics or else lognormal statistics. These results thus generalize and strengthen the conclusions of Papers I & II.

As remarked in Section 4, some $`P(\mathrm{log}E^{})`$ distributions in the transition region (Figures 5 and 7, but not Figures 4 and 9), have approximately power-law form and appear strongly reminiscent of wave collapse at first sight. On detailed examination, as shown next, the distributions are inconsistent with existing collapse theory. First, existing theory predicts that below the peak in the $`P(E^{})`$ distribution the index $`\beta `$ (with $`P(E^{})E^\beta `$) should vary with the dimension $`D`$ as $`\beta =2D1`$ for isotropic, $`2D3`$ for prolate, and $`1`$ for oblate collapse \[Robinson & Newman 1990, Robinson 1997\]. Converting from the spectral index for $`P(\mathrm{log}E^{})`$ to that for $`P(E^{})`$ using $`E^{}P(E^{})=P(\mathrm{log}E^{})`$, the observed values of $`\beta `$ are $`2.1\pm 0.3`$ and $`2.2\pm 0.4`$ for phase bins 450 and 560, respectively. Isotropic collapse theory predicts $`\beta =3`$ and $`5`$, for $`D=2`$ and $`3`$, while prolate collapse predicts $`1`$ and $`3`$, respectively. The exponents observed at low $`E^{}`$ are thus inconsistent with collapse. Second, the indices predicted for high $`E^{}`$ above the peak in $`P(E^{})`$, with $`P(E^{})E^\alpha `$, are $`\alpha =D+2`$, $`D+3`$, and $`2D+1`$ for isotropic, prolate, and oblate collapse. For $`D=23`$, these predictions are the integers $`(4,5)`$, $`(5,6`$, and $`(5,7)`$, respectively. The values observed for phase bins $`450`$ and $`560`$ are $`5.1\pm 0.5`$ and $`6.8\pm 0.9`$, respectively. While the predicted and observed ranges thus overlap for $`D2`$, simultaneous consistency requires oblate collapse with $`D=2`$ for phase $`450`$ and $`D=3`$ for phase $`560`$. One theoretical difficulty here is, intuitively, that strong magnetization of pulsar magnetospheres is expected to make $`D`$ less than $`3`$ (cf. Asseo et al. 1990, Asseo 1996, Weatherall 1998). Moreover, the field statistics for collapse in pulsar magnetospheres need to be investigated (cf. Weatherall1997, 1998).

The analyses in Section 4 and the smooth evolution in fit parameters in Figure 14 thus show that the $`P(\mathrm{log}E^{})`$ distributions observed in the transition regions (phases 435 – 460 and 545 – 610) are best interpreted as the vectorial convolution of a lognormal component with a Gaussian intensity component and not as an intrinsic power-law component. SOC, modulational instability, wave collapse, and driven thermal waves are then inconsistent with the data. Moreover, these results show that the existence of a power-law trend in the high-$`E^{}`$ field statistics does not automatically provide evidence for SOC, modulational instability, etc. but may instead be evidence for an SGT wave population convolved vectorially with a second population distributed as a Gaussian in $`I^{}`$ or another lognormal (implying a second SGT wave population).

The results in Sections 4 – 6 therefore provide no evidence for an intrinsic power-law tail or for a nonlinear cutoff in the observed $`P(\mathrm{log}E^{})`$ distributions at high $`E^{}`$, which might have corresponded to nonlinear self-focusing or decay processes, respectively. This confirms, strengthens, and generalizes to almost the entire on-pulse domain earlier results for Vela (Papers I & II): the simplest interpretation is that the field statistics are consistent with the emission processes at these phases being purely linear (either direct or indirect) and inconsistent with nonlinear emission mechanisms like wave collapse \[Asseo et al. 1990, Asseo 1996, Weatherall 1998\]. More generally, since spatial and temporal variations across the source might be important, the observations constrain viable nonlinear mechanisms to yield lognormal statistics when averaged over the ensemble of emitting structures in the source. This constraint is strong. Wave collapse, for instance, yields power-law statistics when ensemble-averaged over a homogeneous source \[Robinson & Newman 1990, Robinson 1997\], and so is implausible.

The origin of the second population contributing to the wave statistics on-pulse is of interest. Note that contributions from the background sky and supernova remnant are expected to be negligible due to removal of the DC offset. Where the second population has Gaussian intensity statistics, it is interpretable in terms of measurement noise (e.g., thermal receiver noise), scattering and diffusion of pulsar radiation, and (less probably) superposition of radiation from multiple subsources. The close similarity of the Gaussian’s parameters in the outer transition regions with those off-pulse suggests that measurement noise is the most probable explanation there. Scattering is not ruled out, however, and appears to be an attractive way to interpret the increases in the Gaussian fit parameter $`I_0^{}`$ in the inner portions of the transition regions.

Where the second population is lognormal, it is interpretable in terms of emission from a second linear SGT system, either via a second emission process if generated in the same source as the primary population or else emission via the same or a different process in a second source region. Since the second lognormal’s fit parameters are so different from those for the primary lognormal component, being restricted to fields $`E^{}`$ close to the off-pulse background, these parameters may be contaminated by background effects. It might be tempting to suppose that only receiver noise need be considered, but the evidence for evolution in the Gaussian component for Vela (Section 6) and other pulsars \[Cairns et al. 2002b\] and the transitions from the best fitting function being Gaussian-lognormal to being double-lognormal (and vice versa) with phase argues against this. Finally, while the present analyses do not distinguish between subpulse and microstructure effects, analyzing both without artificial distinction, it could be that the primary lognormal component is primarily associated with microstructure physics while the second component is associated with subpulse effects, scattering, and measurement noise.

The preceding analyses and comments specifically exclude Vela’s phase bins 430 – 434, where Johnston et al. (2001) observed “giant micropulses” with fluxes greater than 10 times the average flux at those phases. These giant micropulses have a cumulative probability distribution of the flux $`F`$, $`_{\mathrm{}}^F𝑑FP(F)`$, that appears power-law in $`F`$ with index $`2.85`$ with an associated uncertainty (by eye) $`\pm 0.3`$ \[Kramer et al. 2002\]. Converting into a $`P(E^{})`$ distribution using (5) yields $`P(E^{})E^\alpha `$ with $`\alpha =6.7\pm 0.6`$ for Vela \[Cairns 2002\]. Similar results are found for other pulsars with known giant micropulses or giant pulses and, moreover, this index lies within the range $`47`$ predicted for wave collapse but lies outside the range $`0.53`$ for known SOC phenomena \[Cairns 2002\]. Further work on collapse theory is thus recommended for these giant phenomena, along the lines mentioned above. However, SGT also remains viable provided the power-law distributions result from vector convolution of multiple wave populations, as found in Section 4 for the Vela pulsar’s normal emissions.

The foregoing results point to a richness in possible emission mechanisms and source plasma characteristics for pulsars and, more generally, for solar system and other astrophysical sources. Specifically, SGT appears to apply widely (Robinson & Cairns 2001, Papers I, II and references therein) to plasma waves and propagating radio emissions, including both Vela and other pulsars, but does not apply to all emissions. Examples include the giant micropulse and pulse phenomena discussed just above, Jovian “S” bursts \[Queinnec & Zarka 2001\], which have power-law flux statistics with an index of $`2.0\pm 0.5`$ that lies in the range expected for SOC, and solar decimetric spike bursts which likely have exponential statistics \[Isliker & Benz 2001\]. These results show that analyses of field statistics can strongly constrain the source physics and mechanisms of natural emissions.

## 8 Conclusions

The variability and emission processes of the Vela pulsar can be investigated and strongly constrained using the statistics of the observed time-varying flux. The strong evolution of the pulsar’s field statistics with pulse phase \[Paper II\], from Gaussian in the intensity off-pulse to approximately power-law in $`E^{}`$ in the transition regions, where the average pulse profile is increasing or decreasing, and then lognormal near the peak of the average profile, show that field statistics allow probing of the source characteristics. Moreover, together with single-component fits to the pulsar’s field statistics (Papers I & II), these variations are strong evidence that multiple wave populations are often superposed.

This paper presents detailed two-component fits and associated interpretations of the Vela pulsar’s field statistics, analyzing microstructure and subpulse effects simultaneously, using a new prediction for the amplitude statistics of the vector sum of two transverse fields with known statistics \[Cairns et al. 2002c\]. Excluding Vela’s giant micropulses, it is shown that the approximately power-law field distributions observed in the transition regions are very well fitted by vector convolution of a lognormal with a Gaussian distribution in the intensity, the Gaussian distribution has properties very similar to the Gaussian noise observed at off-pulse phases, and the evolution with phase of the high-$`E^{}`$ power-law form occurs as a result of the lognormal component moving relative to the background level. Accordingly, these power-laws are not intrinsic but are instead evidence for the primary pulsar emission having lognormal statistics. In the phase domain where the field statistics are strongly lognormal to the eye, double-lognormal fits are superior to Gaussian-lognormal fits, with significantly improved domains of applicability and statistical significance compared with a single-component lognormal fit. Moreover, the fit parameters for the lognormal and Gaussian distributions vary smoothly, but also significantly, with phase. Thus, at essentially all on-pulse phases the observed field statistics are very well fitted by the vector combination of a lognormal with either a Gaussian or a second lognormal. This second component may be interpreted as measurement noise, scattered radiation, or multiple superposed sources if it is Gaussian, or a second population of waves produced in a linear SGT system (either via a different mechanism in the same source or alternatively by either the same or a different mechanism in a different source region) if it is lognormal. Put another way, the Vela pulsar’s emission above background is always well represented in terms of at least one lognormal. Accordingly, these results demonstrate in detail that the Vela pulsar’s field statistics are consistent with SGT being relevant, and with the emission mechanisms being purely linear, whenever the pulsar is detectable above background. The pulsar’s variability thus corresponds to well-defined field statistics that are consistent with emission via a linear plasma instability in an SGT state. Similar conclusions are reached elsewhere \[Cairns et al. 2002b\] for PSRs B1641-45 and B0950+08.

The finding that Vela’s lognormal parameter $`\sigma `$ is large near the pulse edges and small near the center is consistent with earlier results for pulsar modulation indices \[Taylor et al. 1975, Johnston & Romani 2002\]. These variations in $`\sigma `$ and also $`\mu `$ with phase will allowing future probing of the source plasma, particularly when a theory for the pulsar’s SGT parameters $`\mu `$ and $`\sigma `$ is developed. Finally, the richness observed in the field statistics of pulsars and of solar system emissions, plus the important constraints placed here on pulsar physics, demonstrate the power and utility of analyzing field statistics of astrophysical sources.

## Acknowledgments

The authors thank B. J. Rickett for helpful conversations and the Australian Research Council for financial support.