OscProb::PMNS_Sterile Class Reference

#include <PMNS_Sterile.h>

Inheritance diagram for OscProb::PMNS_Sterile:
OscProb::PMNS_Base

List of all members.

Public Member Functions

 PMNS_Sterile (int NumNus)
virtual ~PMNS_Sterile ()

Protected Member Functions

virtual void SolveHam ()
 Solve the full Hamiltonian for eigenvectors and eigenvalues.

Detailed Description

Definition at line 19 of file PMNS_Sterile.h.


Constructor & Destructor Documentation

PMNS_Sterile::PMNS_Sterile ( int  NumNus  ) 

Definition at line 47 of file PMNS_Sterile.cxx.

References fEvalGSL, fEvecGSL, OscProb::PMNS_Base::fNumNus, H_GSL, and W_GSL.

00047                                      : PMNS_Base(numNus)
00048 {
00049 
00050   // Allocate memory for the GSL objects
00051   fEvalGSL = gsl_vector_alloc(fNumNus);
00052   fEvecGSL = gsl_matrix_complex_alloc(fNumNus, fNumNus);
00053   H_GSL = gsl_matrix_complex_alloc(fNumNus, fNumNus);
00054   W_GSL = gsl_eigen_hermv_alloc(fNumNus);
00055 
00056 }

PMNS_Sterile::~PMNS_Sterile (  )  [virtual]

Definition at line 58 of file PMNS_Sterile.cxx.

References fEvalGSL, fEvecGSL, H_GSL, and W_GSL.

00058                            {
00059 
00060   // Free memory from GSL objects
00061   gsl_matrix_complex_free(H_GSL); H_GSL = 0;
00062   gsl_eigen_hermv_free(W_GSL); W_GSL = 0;
00063   gsl_matrix_complex_free(fEvecGSL); fEvecGSL = 0;
00064   gsl_vector_free(fEvalGSL); fEvalGSL = 0;
00065 
00066 }


Member Function Documentation

void PMNS_Sterile::SolveHam (  )  [protected, virtual]

Solve the full Hamiltonian for eigenvectors and eigenvalues.

Solve the full Hamiltonian for eigenvectors and eigenvalues This uses GSL to solve the eigensystem.

The matter effect assumes # of neutrons = # of electrons.

Reimplemented from OscProb::PMNS_Base.

Definition at line 75 of file PMNS_Sterile.cxx.

References OscProb::PMNS_Base::BuildHms(), OscProb::PMNS_Base::fBuiltHms, OscProb::PMNS_Base::fEnergy, OscProb::PMNS_Base::fEval, fEvalGSL, OscProb::PMNS_Base::fEvec, fEvecGSL, OscProb::PMNS_Base::fGotES, OscProb::PMNS_Base::fHms, OscProb::PMNS_Base::fIsNuBar, OscProb::PMNS_Base::fNe, OscProb::PMNS_Base::fNumNus, H_GSL, OscProb::PMNS_Base::kGeV2eV, OscProb::PMNS_Base::kGf, OscProb::PMNS_Base::kK2, and W_GSL.

00076 {
00077 
00078   // Check if anything has changed before recalculating
00079   if(!fBuiltHms){
00080     BuildHms();
00081   }
00082 
00083   double lv = 2 * kGeV2eV*fEnergy;     // 2E in eV 
00084   double kr2GNe = kK2*M_SQRT2*kGf*fNe; // Matter potential in eV
00085 
00086   // Finish building Hamiltonian in matter with dimension of eV
00087   for(size_t i=0;i<size_t(fNumNus);i++){
00088     complex buf = fHms[i][i]/lv;
00089     *gsl_matrix_complex_ptr(H_GSL, i, i) = gsl_complex_rect(real(buf),0);
00090     for(size_t j=i+1;j<size_t(fNumNus);j++){
00091       buf = fHms[i][j]/lv;
00092       if(!fIsNuBar) *gsl_matrix_complex_ptr(H_GSL, j, i) = gsl_complex_rect(real(buf),-imag(buf));
00093       else          *gsl_matrix_complex_ptr(H_GSL, j, i) = gsl_complex_rect(real(buf), imag(buf));
00094     }
00095     if(i>2){
00096       // Subtract NC coherent forward scattering from sterile neutrinos. See arXiv:hep-ph/0606054v3, eq. 3.15, for example. 
00097       if(!fIsNuBar) *gsl_matrix_complex_ptr(H_GSL, i, i) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL,i,i) ,  kr2GNe/2);
00098       else          *gsl_matrix_complex_ptr(H_GSL, i, i) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL,i,i) , -kr2GNe/2);;
00099     }
00100   }
00101   // Add nue CC coherent forward scattering from sterile neutrinos. 
00102   if(!fIsNuBar) *gsl_matrix_complex_ptr(H_GSL, 0, 0) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL,0,0) ,  kr2GNe);
00103   else          *gsl_matrix_complex_ptr(H_GSL, 0, 0) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL,0,0) , -kr2GNe);
00104 
00105   // Solve Hamiltonian for eigensystem
00106   gsl_eigen_hermv(H_GSL, fEvalGSL, fEvecGSL, W_GSL);
00107 
00108   // Fill fEval and fEvec vectors from GSL objects
00109   for(int i=0;i<fNumNus;i++){
00110     fEval[i] = gsl_vector_get(fEvalGSL,i);
00111     for(int j=0;j<fNumNus;j++){
00112       gsl_complex buf = gsl_matrix_complex_get(fEvecGSL,i,j);
00113       fEvec[i][j] = complex( GSL_REAL(buf), GSL_IMAG(buf) );
00114     }
00115   }
00116   
00117   fGotES = true;
00118   
00119 }


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

Generated on 15 Jul 2018 for loon by  doxygen 1.6.1