OscProb::PMNS Class Reference

#include <PMNS.h>

List of all members.

Public Member Functions

 PMNS ()
virtual ~PMNS ()
void PrintMix () const
 Print the oscillation matrix.
void PrintDeltaMsqrs () const
 Print the matrix of mass-squared differneces.
double P (int i, int j) const
void SetMix (double th12, double th23, double th13, double deltacp)
void SetMixBWCP (double th1, double th2, double th3, double deltacp)
void SetDeltaMsqrs (double dm21, double dm32)
void PropVacuum (double L, double E, int anti)
void PropMatter (double L, double E, double Ne, int anti)
void Reset ()

Private Types

typedef std::complex< double > complex

Private Member Functions

void Multi (complex A[][3], const complex B[][3], const complex C[][3])
void EvalEqn2 (complex A[][3], const complex U[][3], const complex Udagg[][3], const double dmsqr[][3], double L, double E)
void EvalEqn5 (complex twoEH[][3], const complex U[][3], const complex Udagg[][3], const double dmsqr[][3], double E, double Gf, double Ne)
void EvalEqn10 (complex A[][3], const complex U[][3], const complex X[][3], const complex Udagg[][3])
void EvalEqn11 (complex X[][3], double L, double E, const complex twoEH[][3], const double Msqr[], const double dMsqr[][3])
void EvalEqn21 (double Msqr[], double alpha, double beta, double gamma)
void EvalEqn22 (double &alpha, double &beta, double &gamma, double E, double Gf, double Ne, const double dmsqr[][3], const complex U[][3])
void SortEigenvalues (double dMsqr[][3], const double dmsqr[][3], const double MsqrVac[], double Msqr[])
void DumpMatrix (const complex M[][3]) const
 Utility to print a genetrix complex matrix M.

Private Attributes

double fdmsqr [3][3]
 m^2_i - m^2_j in vacuum
complex fU [3][3]
 PMNS Matrix in all its forms.
complex fUdagg [3][3]
complex fUstar [3][3]
complex fUtran [3][3]
complex fA [3][3]
 The neutrino transition matrix.

Detailed Description

Definition at line 37 of file PMNS.h.


Member Typedef Documentation

typedef std::complex<double> OscProb::PMNS::complex [private]

Definition at line 93 of file PMNS.h.


Constructor & Destructor Documentation

PMNS::PMNS (  ) 

Definition at line 56 of file PMNS.cxx.

References Reset(), SetDeltaMsqrs(), and SetMix().

00057 {
00058   this->SetMix(0.,0.,0.,0.);
00059   this->SetDeltaMsqrs(0.,0.);
00060   this->Reset();
00061 }

PMNS::~PMNS (  )  [virtual]

Definition at line 63 of file PMNS.cxx.

00063            {
00064 }


Member Function Documentation

void PMNS::DumpMatrix ( const complex  M[][3]  )  const [private]

Utility to print a genetrix complex matrix M.

Definition at line 68 of file PMNS.cxx.

Referenced by PrintMix().

00069 {
00070   std::cout 
00071     <<"| "<<M[0][0]<<"\t"<<M[0][1]<<"\t"<<M[0][2]<<" |\n"
00072     <<"| "<<M[1][0]<<"\t"<<M[1][1]<<"\t"<<M[1][2]<<" |\n"
00073     <<"| "<<M[2][0]<<"\t"<<M[2][1]<<"\t"<<M[2][2]<<" |\n"
00074     <<std::endl;
00075 }

void PMNS::EvalEqn10 ( complex  A[][3],
const complex  U[][3],
const complex  X[][3],
const complex  Udagg[][3] 
) [private]

Find the transition matrix for matter using Eqn.10

Parameters:
A - the output transition matrix
U - the neutrino mixing matrix
X - matrix resulting from Eqn.11
Udagg - the adjoint of the mixing matrix

Compute Equation 10, the transition matrix for propagation across matter slab of constant density, from Barger et al:

A = U * X * U^(dagger) for neutrinos

Definition at line 273 of file PMNS.cxx.

References Multi().

Referenced by PropMatter().

00277 {
00278   complex tmp[3][3];
00279   this->Multi(tmp, X, Udagg);
00280   this->Multi(A,   U, tmp);
00281 }

void PMNS::EvalEqn11 ( complex  X[][3],
double  L,
double  E,
const complex  twoEH[][3],
const double  Msqr[],
const double  dMsqr[][3] 
) [private]

Compute the X matrix from Eqn.11

Parameters:
X - the "X" matrix. Sorry, no better name...
L - the neutrino flight distance in km
E - the neutrino energy in GeV
twoEH - the Hamiltonian in units of eV
Msqr - the matter eigenvalues in units of eV^2
dMsqr - the matrix of eigenvalue differences M_i^2-M_j^2

Evaluate Eqn. 11, the Lagragne formula for the matrix e(-iHt), from Barger et al.

Definition at line 288 of file PMNS.cxx.

References MuELoss::a, kK1, Multi(), and zero().

Referenced by PropMatter().

00294 {
00295   // The identity matrix
00296   static const double One[3][3] = {{1.,0.,0.},
00297                                    {0.,1.,0.},
00298                                    {0.,0.,1.}
00299   };
00300   register int a, b, k;
00301   complex phase;
00302   complex EHM0[3][3];
00303   complex EHM1[3][3];
00304   complex EHM2[3][3];
00305   complex EHM21[3][3];
00306   complex EHM20[3][3];
00307   complex EHM10[3][3];
00308 
00309   // There are three matricies which can apper inside the product on
00310   // j!=k. Calculate them before hand
00311   for (a=0; a<3; ++a) {
00312     for (b=0; b<3; ++b) {
00313       EHM0[a][b] = twoEH[a][b]-Msqr[0]*One[a][b];
00314       EHM1[a][b] = twoEH[a][b]-Msqr[1]*One[a][b];
00315       EHM2[a][b] = twoEH[a][b]-Msqr[2]*One[a][b];
00316     }
00317   }
00318   this->Multi(EHM21,EHM2,EHM1);
00319   this->Multi(EHM20,EHM2,EHM0);
00320   this->Multi(EHM10,EHM1,EHM0);
00321 
00322   // Refer to Msqr_j as dMsqr[j][0] since only mass differences matter
00323   for (a=0; a<3; ++a) {
00324     for (b=0; b<3; ++b) {
00325       X[a][b] = zero;
00326       for (k=0; k<3; ++k) {
00327         phase = exp(complex(0.0,-kK1*dMsqr[k][0]*L/E));
00328         if (k==0) {
00329           X[a][b] += (EHM21[a][b]/(dMsqr[k][2]*dMsqr[k][1]))*phase;
00330         }
00331         else if (k==1) { 
00332           X[a][b] += (EHM20[a][b]/(dMsqr[k][2]*dMsqr[k][0]))*phase;
00333         }
00334         else if (k==2) {
00335           X[a][b] += (EHM10[a][b]/(dMsqr[k][1]*dMsqr[k][0]))*phase;
00336         } // Switch for product on j!=k
00337       } // Sum on k
00338     } // Loop on b
00339   } // Loop on a
00340 }

void PMNS::EvalEqn2 ( complex  A[][3],
const complex  U[][3],
const complex  Udagg[][3],
const double  dmsqr[][3],
double  L,
double  E 
) [private]

Find the transition matrix for vacuum oscillations using Eqn.2

Parameters:
A - output transition matrix
U - neutrino mixing matrix
Udagg - adjoint of neutrino mixing matrix
dmsqr - matrix of mass-splittings m_i^2-m_j^2 in eV^2
L - neutrino flight distance in km
E - neutrino energy in GeV

Compute Equation 2, transition matrix for propagation through vaccum, from Barger et al:

A(nua->nub) = sum_i U_ai exp(-1/2 i m_i^2 L/E) Udagg_ib

Definition at line 224 of file PMNS.cxx.

References MuELoss::a, kK1, and zero().

Referenced by PropVacuum().

00230 {
00231   register int a, b, i;
00232   for (a=0; a<3; ++a) {
00233     for (b=0; b<3; ++b) {
00234       A[a][b] = zero;
00235       for (i=0; i<3; ++i) {
00236         complex phase(0.0,-kK1*dmsqr[i][0]*L/E);
00237         A[a][b] += U[a][i]*exp(phase)*Udagg[i][b];
00238       }
00239     }
00240   }
00241 }

void PMNS::EvalEqn21 ( double  Msqr[],
double  alpha,
double  beta,
double  gamma 
) [private]

Compute the matter eigenvalues according to Eqn.21

Parameters:
Msqr - Output eigenvalues in eV^2
alpha - See Eqn22
beta - See Eqn22
gamma - See Eqn22

Definition at line 347 of file PMNS.cxx.

Referenced by PropMatter().

00351 {
00352   double arg;      // Argument of the acos()
00353   double theta0;   // First of the three roots of acos()
00354   double theta1;   // Second of the three roots of acos()
00355   double theta2;   // Third of the three roots of acos()
00356   double fac;      // Factor in front of cos() terms
00357   double constant; // Offset for all eigenvalues
00358   static const double k2PiOver3 = 2.0*M_PI/3.0;
00359   double alpha2           = alpha*alpha;
00360   double alpha3           = alpha*alpha2;
00361   double alpha2Minus3beta = alpha2-3.0*beta;
00362 
00363   arg =
00364     (2.0*alpha3 - 9.0*alpha*beta + 27.0*gamma)/
00365     (2.0*pow(alpha2Minus3beta,1.5));
00366 
00367   // Occasionally round off errors mean that arg wanders outside of
00368   // its allowed range. If its way off (1 part per billion), stop the
00369   // program. Otherwise, set it to its rail value.
00370 
00371   //  The code had a small bug:
00372   //   if( fabs(arg-1.0)>1.E-9 ) abort();
00373   //  I have changed this to:
00374   //   if( fabs(arg)-1.0>1.E-9 ) abort();
00375   
00376   if (fabs(arg)>1.0) {
00377     // if( fabs(arg-1.0)>1.E-9 ) abort();
00378     if (fabs(arg)-1.0>1.E-9) abort();
00379     if (arg<-1.0) arg = -1.00;
00380     if (arg>+1.0) arg = +1.00;
00381   }
00382   
00383   // The three roots, find the first by computing the acos() the
00384   // others are nearby
00385   theta0 = acos(arg)/3.0;
00386   theta1 = theta0-k2PiOver3;
00387   theta2 = theta0+k2PiOver3;
00388   
00389   // The multiplier and offset
00390   fac      = -2.0*sqrt(alpha2Minus3beta)/3.0;
00391   constant = -alpha/3.0; // The constant offset m1^2 is irrelevant
00392   
00393   // And the eigenvalues themselves
00394   Msqr[0] = fac*cos(theta0) + constant;
00395   Msqr[1] = fac*cos(theta1) + constant;
00396   Msqr[2] = fac*cos(theta2) + constant;
00397 }

void PMNS::EvalEqn22 ( double &  alpha,
double &  beta,
double &  gamma,
double  E,
double  Gf,
double  Ne,
const double  dmsqr[][3],
const complex  U[][3] 
) [private]

Compute the expressions alpha, beta, gamma which appear in Eqns.21 and 22

Parameters:
alpha - output parameter
beta - output parameter
gamma - output parameter
E - neutrino energy in GeV
Gf - Fermi's constant in units of GeV^-2
Ne - electron number denstiy in mole/cm^3
dmsqr - matrix of vacuum mass-splittings m_i^2-m_j^2
U - neutrino mixing matrix

.....................................................................

Compute the values of the simplifying expressions alpha, beta, and gamma. This is Eqn22 from the Barger et al paper

Definition at line 404 of file PMNS.cxx.

References kGeV2eV, and kK2.

Referenced by PropMatter().

00412 {
00413   // 2*sqrt(2)*Gf*Ne*E in units of eV^2
00414   double k2r2EGNe = kK2*2.0*M_SQRT2*Gf*Ne*(kGeV2eV*E);
00415   
00416   alpha = k2r2EGNe + dmsqr[0][1] + dmsqr[0][2];
00417   
00418   beta  =
00419     dmsqr[0][1]*dmsqr[0][2] + 
00420     k2r2EGNe*(dmsqr[0][1]*(1.0-norm(U[0][1])) + 
00421               dmsqr[0][2]*(1.0-norm(U[0][2])));
00422 
00423   gamma = k2r2EGNe*dmsqr[0][1]*dmsqr[0][2]*norm(U[0][0]);
00424 }

void PMNS::EvalEqn5 ( complex  twoEH[][3],
const complex  U[][3],
const complex  Udagg[][3],
const double  dmsqr[][3],
double  E,
double  Gf,
double  Ne 
) [private]

Find the Hamiltonian using Eqn.5

Parameters:
twoEH - Hamiltonian matrix scaled by 2E in units of eV
U - neutrino mixing matrix
Udagg - adjoint of neutrino mixing matrix
dmsqr - Matrix of mass-splittings m_i^2-m_j^2 in eV^2
E - neutrino energy in GeV
Gf - Fermi's constant in units of GeV^-2
Ne - electron number density in mole/cm^3

Compute Eqn. 5, the matter 2*E*Hamiltonian in the mass basis

Definition at line 247 of file PMNS.cxx.

References kGeV2eV, kK2, and zero().

Referenced by PropMatter().

00254 {
00255   register int j, k;
00256   double       k2r2GNeE = kK2*2.0*M_SQRT2*Gf*Ne*(kGeV2eV*E);
00257   for (k=0; k<3; ++k) {
00258     for (j=0; j<3; ++j) {
00259       twoEH[k][j] = zero;
00260       if (k==j) twoEH[k][j] = dmsqr[j][0];
00261       twoEH[k][j] -= k2r2GNeE*U[0][j]*Udagg[k][0];
00262     }
00263   }
00264 }

void PMNS::Multi ( complex  A[][3],
const complex  B[][3],
const complex  C[][3] 
) [private]

Multiply two matricies: A = B*C

Parameters:
A - output matrix
B - input matrix
C - input matrix

Compute matrix multiplication A = B*C

Definition at line 200 of file PMNS.cxx.

References zero().

Referenced by EvalEqn10(), EvalEqn11(), PropMatter(), and PropVacuum().

00201 {
00202   register int i, j, k;
00203   for (i=0; i<3; ++i) {
00204     for (j=0; j<3; ++j) {
00205       A[i][j] = zero;
00206       for (k=0; k<3; ++k) {
00207         complex b = B[i][k];
00208         complex c = C[k][j];
00209         complex bc = b*c;
00210 
00211         A[i][j] += B[i][k]*C[k][j];
00212       }
00213     }
00214   }
00215 }

double PMNS::P ( int  i,
int  j 
) const

Return the oscillation probability from flavor i to flavor j

Parameters:
i - initial flavor (0,1,2) = (nue,numu,nutau)
j - final flavor (0,1,2) = (nue,numu,nutau)

Compute oscillation probability from flavor i to flavor j

0 = nue, 1 = numu, 2 = nutau

Definition at line 568 of file PMNS.cxx.

References fA.

Referenced by OscFit::OscillationCalculator::OscProb(), and AtNuOscillate::OscProb().

00569 {
00570   assert(i>=0 && i<3);
00571   assert(j>=0 && j<3);
00572   return norm(fA[i][j]);
00573 }

void PMNS::PrintDeltaMsqrs (  )  const

Print the matrix of mass-squared differneces.

Definition at line 154 of file PMNS.cxx.

References fdmsqr.

Referenced by OscFit::OscillationCalculator::PrintPMNS().

00155 {
00156   std::cout 
00157     <<"|"<<fdmsqr[0][0]<<"\t"<<fdmsqr[0][1]<<"\t"<<fdmsqr[0][2]<<"|\n" 
00158     <<"|"<<fdmsqr[1][0]<<"\t"<<fdmsqr[1][1]<<"\t"<<fdmsqr[1][2]<<"|\n" 
00159     <<"|"<<fdmsqr[2][0]<<"\t"<<fdmsqr[2][1]<<"\t"<<fdmsqr[2][2]<<"|"
00160     << std::endl;
00161 }

void PMNS::PrintMix (  )  const

Print the oscillation matrix.

Definition at line 79 of file PMNS.cxx.

References DumpMatrix(), and fU.

Referenced by OscFit::OscillationCalculator::PrintPMNS().

00079 { this->DumpMatrix(fU); }

void PMNS::PropMatter ( double  L,
double  E,
double  Ne,
int  anti 
)

Propagate a neutrino through a slab of matter

Parameters:
L - length of slab (km)
E - neutrino energy in GeV
Ne - electron number density of matter in mole/cm^3
anti - +1 = neutrino case, -1 = anti-neutrino case

.....................................................................

Update the transition matrix for a step across a slab of matter

Definition at line 490 of file PMNS.cxx.

References EvalEqn10(), EvalEqn11(), EvalEqn21(), EvalEqn22(), EvalEqn5(), fA, fdmsqr, fU, fUdagg, fUstar, fUtran, OscPar::Gf, Multi(), and SortEigenvalues().

Referenced by AtNuOscillate::PropagateAtNuMatter(), and OscFit::OscillationCalculator::PropMatter().

00491 {
00492   static const double  Gf = 1.166371E-5; // G_F in units of GeV^-2
00493   register int i, j;
00494   complex twoEH[3][3];
00495   complex X[3][3];
00496   double  Msqr[3];
00497   double  MsqrV[3];
00498   double  dMsqr[3][3];
00499   double  alpha,  beta,  gamma;
00500   double  alpha0, beta0, gamma0;
00501   complex A[3][3];
00502   complex Aout[3][3];
00503   
00504   // Find the transition matrix. The series of steps are to:
00505   if (anti>0) {
00506     // [1] Find the matter Hamiltonian in the mass basis...
00507     this->EvalEqn5(twoEH, fU, fUdagg, fdmsqr, E, Gf, Ne);
00508 
00509     // [2] Find the eigenvalues and sort them.
00510     this->EvalEqn22(alpha, beta, gamma, E, Gf, Ne, fdmsqr, fU);
00511     this->EvalEqn21(Msqr,  alpha, beta, gamma);
00512 
00513     // Repeat the above, but for vacuum (Ne=0.0) to sort out the order
00514     // of the eigenvalues
00515     this->EvalEqn22(alpha0, beta0, gamma0, E, 0.0, 0.0, fdmsqr, fU);
00516     this->EvalEqn21(MsqrV, alpha0, beta0, gamma0);
00517     this->SortEigenvalues(dMsqr, fdmsqr, MsqrV, Msqr);
00518 
00519     // [3] Evaluate the transition matrix
00520     this->EvalEqn11(X, L, E, twoEH, Msqr, dMsqr);
00521     this->EvalEqn10(A, fU, X, fUdagg);
00522   }
00523   else if (anti<0) {
00524     // As above, but make required substitutions for anti-neutrinos:
00525     // Gf=>-Gf, U=>Ustar, U^dagger=>U^dagger^*=U^T
00526     this->EvalEqn5(twoEH, fUstar, fUtran, fdmsqr, E, -Gf, Ne);
00527     this->EvalEqn22(alpha, beta, gamma, E, -Gf, Ne, fdmsqr, fUstar);
00528     this->EvalEqn21(Msqr,  alpha, beta, gamma);
00529     this->EvalEqn22(alpha0, beta0, gamma0, E, 0.0, 0.0, fdmsqr, fUstar);
00530     this->EvalEqn21(MsqrV, alpha0, beta0, gamma0);
00531     this->SortEigenvalues(dMsqr, fdmsqr, MsqrV, Msqr);
00532     this->EvalEqn11(X, L, E, twoEH, Msqr, dMsqr);
00533     this->EvalEqn10(A, fUstar, X, fUtran);
00534   }
00535   else abort();
00536   
00537   // [4] Apply the transition matrix to the matrix we started with...
00538   this->Multi(Aout, A, fA);
00539   for (i=0; i<3; ++i) {
00540     for (j=0; j<3; ++j) {
00541       fA[i][j] = Aout[i][j];
00542     }
00543   }
00544 }

void PMNS::PropVacuum ( double  L,
double  E,
int  anti 
)

Propagate a neutrino through vacuum

Parameters:
L - flight distance in km
E - neutrino energy in GeV
anti - +1 = neutrino case, -1 = anti-neutrino case

.....................................................................

Update the transition matrix for a step across the vacuum

Definition at line 469 of file PMNS.cxx.

References EvalEqn2(), fA, fdmsqr, fU, fUdagg, fUstar, fUtran, and Multi().

Referenced by AtNuOscillate::PropagateAtNuVacuum(), and OscFit::OscillationCalculator::PropVacuum().

00470 {
00471   register int i, j;
00472   complex A[3][3];
00473   complex Aout[3][3];
00474 
00475   if      (anti>0) this->EvalEqn2(A, fU,     fUdagg, fdmsqr, L, E);
00476   else if (anti<0) this->EvalEqn2(A, fUstar, fUtran, fdmsqr, L, E);
00477   else abort();
00478   this->Multi(Aout, A, fA);
00479   for (i=0; i<3; ++i) {
00480     for (j=0; j<3; ++j) {
00481       fA[i][j] = Aout[i][j];
00482     }
00483   }
00484 }

void PMNS::Reset (  ) 

Erase memory of neutrino propagate and reset transition matrix to unity. Preserves values of mixing and mass-splittings

Reset the transition matrix back to the identity matrix from where it starts

Definition at line 551 of file PMNS.cxx.

References fA, one(), and zero().

Referenced by PMNS(), AtNuOscillate::PropagateAtNuMatter(), AtNuOscillate::PropagateAtNuVacuum(), OscFit::OscillationCalculator::ResetWithNewFlavour(), AtNuOscillate::SetPMNS(), and OscFit::OscillationCalculator::SetPMNS().

00552 {
00553   register int i, j;
00554   for (i=0; i<3; ++i) {
00555     for (j=0; j<3; ++j) {
00556       if (i==j) fA[i][j] = one;
00557       else      fA[i][j] = zero;
00558     }
00559   }
00560 }

void PMNS::SetDeltaMsqrs ( double  dm21,
double  dm32 
)

Set the mass-splittings

Parameters:
dm21 - m2^2-m1^2 in eV^2
dm32 - m3^2-m2^2 in eV^2

Set the mass-splittings. These are m_2^2-m_1^2 and m_3^2-m_2^2 in eV^2

Definition at line 168 of file PMNS.cxx.

References fdmsqr.

Referenced by PMNS(), AtNuOscillate::SetPMNS(), and OscFit::OscillationCalculator::SetPMNS().

00169 {
00170   double eps = 5.0E-9;
00171   double msqr[3];
00172   
00173   msqr[0] = 0.0;
00174   msqr[1] = dm21;
00175   msqr[2] = dm21+dm32;
00176   
00177   // Degeneracies cause problems with diagonalization, so break them
00178   // ever so slightly
00179   if (dm21==0.0) {msqr[0] -= 0.5*eps; msqr[1] += 0.5*eps; }
00180   if (dm32==0.0) {msqr[1] -= 0.5*eps; msqr[2] += 0.5*eps; }
00181 
00182   // Assign the mass splittings fdmsqr_ij = msqr_i - msqr_j by
00183   // convention
00184   for (register int i=0; i<3; ++i) {
00185     for (register int j=0; j<3; ++j) {
00186       // A somewhat subtle point: Barger et al. refer to the sign of
00187       // m1-m2 being positive which implies dm^2_12>0 and
00188       // dm^2_21<0. The labeling in more common use is to assume m3 is
00189       // heaviest such that dm_12<0 and dm_21>0. Rather than reverse
00190       // all the indicies in all the equations, flip the signs here.
00191       fdmsqr[i][j] = -(msqr[i] - msqr[j]);
00192     }
00193   }
00194 }

void PMNS::SetMix ( double  th12,
double  th23,
double  th13,
double  deltacp 
)

Set the parameters of the PMNS matrix

Parameters:
th12 - The angle theta_12 in radians
th23 - The angle theta_23 in radians
th13 - The angle theta_13 in radians
deltacp - The CPV phase delta_CP in radians

Definition at line 83 of file PMNS.cxx.

References fU, fUdagg, fUstar, and fUtran.

Referenced by PMNS(), AtNuOscillate::SetPMNS(), and OscFit::OscillationCalculator::SetPMNS().

00084 {
00085   register int i, j;
00086   double  s12, s23, s13, c12, c23, c13;
00087   complex idelta(0.0,deltacp);
00088   
00089   s12 = sin(th12);  s23 = sin(th23);  s13 = sin(th13);
00090   c12 = cos(th12);  c23 = cos(th23);  c13 = cos(th13);
00091   
00092   fU[0][0] =  c12*c13;
00093   fU[0][1] =  s12*c13;
00094   fU[0][2] =  s13*exp(-idelta);
00095   
00096   fU[1][0] = -s12*c23-c12*s23*s13*exp(idelta);
00097   fU[1][1] =  c12*c23-s12*s23*s13*exp(idelta);
00098   fU[1][2] =  s23*c13;
00099   
00100   fU[2][0] =  s12*s23-c12*c23*s13*exp(idelta);
00101   fU[2][1] = -c12*s23-s12*c23*s13*exp(idelta);
00102   fU[2][2] =  c23*c13;
00103   
00104   // Compute derived forms of the mixing matrix
00105   for (i=0; i<3; ++i) {
00106     for (j=0; j<3; ++j) {
00107       fUtran[i][j] = fU[j][i];
00108       fUstar[i][j] = conj(fU[i][j]);
00109       fUdagg[i][j] = conj(fU[j][i]);
00110     }
00111   }
00112 }

void PMNS::SetMixBWCP ( double  th1,
double  th2,
double  th3,
double  d 
)

Set the mixing matrix using the old convention used in the original 1980 paper.

Warning:
Useful for testing but should not be used for calculations. Use SetMix above instead
Parameters:
th1,th2,th3,deltacp - parameters of matrix

.....................................................................

Initialize the mixing matrix using the older form referenced by the Barger et al paper

Warning:
This should not be used except for testing. Use SetMix above.

Definition at line 122 of file PMNS.cxx.

References fU, fUdagg, fUstar, fUtran, and id.

00123 {
00124   int i, j;
00125   double s1, s2, s3, c1, c2, c3;
00126   complex id(0.0,d);
00127   s1 = sin(th1);  s2 = sin(th2);  s3 = sin(th3);
00128   c1 = cos(th1);  c2 = cos(th2);  c3 = cos(th3);
00129   
00130   fU[0][0] =  c1;
00131   fU[0][1] =  s1*c3;
00132   fU[0][2] =  s1*s3;
00133 
00134   fU[1][0] = -s1*c2;
00135   fU[1][1] =  c1*c2*c3+s2*s3*exp(id);
00136   fU[1][2] =  c1*c2*s3-s2*c3*exp(id);
00137 
00138   fU[2][0] = -s1*s2;
00139   fU[2][1] =  c1*s2*c3-c2*s3*exp(id);
00140   fU[2][2] =  c1*s2*s3+c2*c3*exp(id);
00141 
00142   // Compute derived forms of the mixing matrix
00143   for (i=0; i<3; ++i) {
00144     for (j=0; j<3; ++j) {
00145       fUtran[i][j] = fU[j][i];
00146       fUstar[i][j] = conj(fU[i][j]);
00147       fUdagg[i][j] = conj(fU[j][i]);
00148     }
00149   }
00150 }

void PMNS::SortEigenvalues ( double  dMsqr[][3],
const double  dmsqr[][3],
const double  MsqrVac[],
double  Msqr[] 
) [private]

Sort the matter eignvalues to that they appear in the same order as the vacuum eigenvalues.

Parameters:
dMsqr - Output matrix of eigenvalue differences M_i^2-M_j^2
dmsqr - the vacuum eigenvalues
MsqrVac - vacuum eigenvalues
Msqr - Eigenvalues in matter. In: unsorted, out: sorted

Sort out the eigenvalues

Definition at line 430 of file PMNS.cxx.

Referenced by PropMatter().

00434 {
00435   int i, j, k;
00436   double best, delta;
00437   double MsqrTmp[3] = {-99.9E9,-99.9E9,-99.9E9};
00438   int flg[3] = {0,0,0};
00439 
00440   // Attempt to figure out which of the eigenvalues match between
00441   // dmsqr and MsqrVac
00442   for (i=0; i<3; ++i) {
00443     best =  1.E30;
00444     k    = -1;
00445     for (j=0; j<3; ++j) {
00446       delta = fabs(MsqrVac[i] - dmsqr[j][0]);
00447       if (delta<best) { best = delta; k = j; }
00448     }
00449     if (best>1.E-9) abort();
00450     if (k<0 || k>2) abort();
00451     flg[k] = 1;
00452     MsqrTmp[i] = Msqr[k];
00453   }
00454   // Check that each eigenvalue is used
00455   for (i=0; i<3; ++i) if (flg[i]!=1) abort();
00456   
00457   for (i=0; i<3; ++i) {
00458     Msqr[i] = MsqrTmp[i];
00459     for (j=0; j<3; ++j) {
00460       dMsqr[i][j] = (MsqrTmp[i] - MsqrTmp[j]);
00461     }
00462   }
00463 }


Member Data Documentation

complex OscProb::PMNS::fA[3][3] [private]

The neutrino transition matrix.

Definition at line 204 of file PMNS.h.

Referenced by P(), PropMatter(), PropVacuum(), and Reset().

double OscProb::PMNS::fdmsqr[3][3] [private]

m^2_i - m^2_j in vacuum

Definition at line 199 of file PMNS.h.

Referenced by PrintDeltaMsqrs(), PropMatter(), PropVacuum(), and SetDeltaMsqrs().

complex OscProb::PMNS::fU[3][3] [private]

PMNS Matrix in all its forms.

Definition at line 200 of file PMNS.h.

Referenced by PrintMix(), PropMatter(), PropVacuum(), SetMix(), and SetMixBWCP().

complex OscProb::PMNS::fUdagg[3][3] [private]

Definition at line 201 of file PMNS.h.

Referenced by PropMatter(), PropVacuum(), SetMix(), and SetMixBWCP().

complex OscProb::PMNS::fUstar[3][3] [private]

Definition at line 202 of file PMNS.h.

Referenced by PropMatter(), PropVacuum(), SetMix(), and SetMixBWCP().

complex OscProb::PMNS::fUtran[3][3] [private]

Definition at line 203 of file PMNS.h.

Referenced by PropMatter(), PropVacuum(), SetMix(), and SetMixBWCP().


The documentation for this class was generated from the following files:

Generated on 3 Oct 2018 for loon by  doxygen 1.6.1