OscProb::PMNS_Transitions Class Reference

#include <PMNS_Transitions.h>

Inheritance diagram for OscProb::PMNS_Transitions:
OscProb::PMNS_Base

List of all members.

Public Member Functions

 PMNS_Transitions ()
virtual ~PMNS_Transitions ()
virtual void SetCMuMu (double cmumu)
virtual void SetCTauTau (double ctautau)
virtual void SetBMu (double bmu)
virtual void SetBTau (double btau)
virtual void InitialiseTransitionParameters ()

Protected Member Functions

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

Protected Attributes

double fcmumu
double fctautau
double fbmu
double fbtau

Detailed Description

Definition at line 19 of file PMNS_Transitions.h.


Constructor & Destructor Documentation

PMNS_Transitions::PMNS_Transitions (  ) 

Definition at line 47 of file PMNS_Transitions.cxx.

References fEvalGSL_Trans, fEvecGSL_Trans, OscProb::PMNS_Base::fNumNus, H_GSL_Trans, InitialiseTransitionParameters(), and W_GSL_Trans.

00047                                    : PMNS_Base(4)
00048 {
00049 
00050   // Allocate memory for the GSL objects
00051   fEvalGSL_Trans = gsl_vector_alloc(fNumNus);
00052   fEvecGSL_Trans = gsl_matrix_complex_alloc(fNumNus, fNumNus);
00053   H_GSL_Trans = gsl_matrix_complex_alloc(fNumNus, fNumNus);
00054   W_GSL_Trans = gsl_eigen_hermv_alloc(fNumNus);
00055   this->InitialiseTransitionParameters();
00056 }

PMNS_Transitions::~PMNS_Transitions (  )  [virtual]

Definition at line 58 of file PMNS_Transitions.cxx.

References fEvalGSL_Trans, fEvecGSL_Trans, H_GSL_Trans, and W_GSL_Trans.

00058                                    {
00059 
00060   // Free memory from GSL objects
00061   gsl_matrix_complex_free(H_GSL_Trans); H_GSL_Trans = 0;
00062   gsl_eigen_hermv_free(W_GSL_Trans); W_GSL_Trans = 0;
00063   gsl_matrix_complex_free(fEvecGSL_Trans); fEvecGSL_Trans = 0;
00064   gsl_vector_free(fEvalGSL_Trans); fEvalGSL_Trans = 0;
00065 
00066 }


Member Function Documentation

void PMNS_Transitions::InitialiseTransitionParameters (  )  [virtual]

Definition at line 69 of file PMNS_Transitions.cxx.

References fbmu, fbtau, fcmumu, and fctautau.

Referenced by PMNS_Transitions().

00070 {
00071   fcmumu=0;
00072   fctautau=0;
00073   fbmu=0;
00074   fbtau=0;
00075 }

virtual void OscProb::PMNS_Transitions::SetBMu ( double  bmu  )  [inline, virtual]

Definition at line 28 of file PMNS_Transitions.h.

References fbmu.

00028 {fbmu=bmu; return;}

virtual void OscProb::PMNS_Transitions::SetBTau ( double  btau  )  [inline, virtual]

Definition at line 29 of file PMNS_Transitions.h.

References fbtau.

00029 {fbtau=btau; return;}

virtual void OscProb::PMNS_Transitions::SetCMuMu ( double  cmumu  )  [inline, virtual]

Definition at line 26 of file PMNS_Transitions.h.

References fcmumu.

00026 {fcmumu=cmumu; return;}

virtual void OscProb::PMNS_Transitions::SetCTauTau ( double  ctautau  )  [inline, virtual]

Definition at line 27 of file PMNS_Transitions.h.

References fctautau.

00027 {fctautau=ctautau; return;}

void PMNS_Transitions::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 84 of file PMNS_Transitions.cxx.

References OscProb::PMNS_Base::BuildHms(), fbmu, fbtau, OscProb::PMNS_Base::fBuiltHms, fcmumu, fctautau, OscProb::PMNS_Base::fEnergy, OscProb::PMNS_Base::fEval, fEvalGSL_Trans, OscProb::PMNS_Base::fEvec, fEvecGSL_Trans, OscProb::PMNS_Base::fGotES, OscProb::PMNS_Base::fHms, OscProb::PMNS_Base::fIsNuBar, OscProb::PMNS_Base::fNumNus, H_GSL_Trans, OscProb::PMNS_Base::kGeV2eV, and W_GSL_Trans.

00085 {
00086 
00087   // Check if anything has changed before recalculating
00088   if(!fBuiltHms){
00089     //This bit sets the dm2 elements, so how do I make two of them equal?
00090     //I think this is also where the magic with the mixing angles happens.
00091     BuildHms();
00092   }
00093 
00094   if (fIsNuBar){cout <<"You're doing this wrong." << endl; assert(false);}
00095 
00096   double lv = 2 * kGeV2eV*fEnergy;     // 2E in eV 
00097 
00098   // Finish building Hamiltonian in matter with dimension of eV
00099   for(size_t i=0;i<size_t(fNumNus);i++){
00100     complex buf = fHms[i][i]/lv;
00101     *gsl_matrix_complex_ptr(H_GSL_Trans, i, i) = gsl_complex_rect(real(buf),0);
00102     for(size_t j=i+1;j<size_t(fNumNus);j++){
00103       buf = fHms[i][j]/lv;
00104       *gsl_matrix_complex_ptr(H_GSL_Trans, j, i) = gsl_complex_rect(real(buf),-imag(buf));
00105     }
00106   }
00107 
00108   //Now add all the transitions parameters
00109   *gsl_matrix_complex_ptr(H_GSL_Trans, 0, 0) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,0,0),
00110                                                               -1.*fEnergy*fcmumu/2.);
00111   *gsl_matrix_complex_ptr(H_GSL_Trans, 1, 1) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,1,1),
00112                                                               -1.*fEnergy*fctautau/2.);
00113   *gsl_matrix_complex_ptr(H_GSL_Trans, 2, 2) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,2,2),
00114                                                               -1.*fEnergy*fcmumu/2.);
00115   *gsl_matrix_complex_ptr(H_GSL_Trans, 3, 3) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,3,3),
00116                                                               -1.*fEnergy*fctautau/2.);
00117 
00118 
00119   *gsl_matrix_complex_ptr(H_GSL_Trans, 2, 0) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,2,0),
00120                                                               fEnergy*fbmu/2.);
00121   *gsl_matrix_complex_ptr(H_GSL_Trans, 3, 1) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,3,1),
00122                                                               fEnergy*fbtau/2.);
00123   *gsl_matrix_complex_ptr(H_GSL_Trans, 0, 2) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,0,2),
00124                                                               fEnergy*fbmu/2.);
00125   *gsl_matrix_complex_ptr(H_GSL_Trans, 1, 3) = gsl_complex_add_real(gsl_matrix_complex_get(H_GSL_Trans,1,3),
00126                                                               fEnergy*fbtau/2.);
00127 
00128   // Solve Hamiltonian for eigensystem
00129   gsl_eigen_hermv(H_GSL_Trans, fEvalGSL_Trans, fEvecGSL_Trans, W_GSL_Trans);
00130 
00131   // Fill fEval and fEvec vectors from GSL objects
00132   for(int i=0;i<fNumNus;i++){
00133     fEval[i] = gsl_vector_get(fEvalGSL_Trans,i);
00134     for(int j=0;j<fNumNus;j++){
00135       gsl_complex buf = gsl_matrix_complex_get(fEvecGSL_Trans,i,j);
00136       fEvec[i][j] = complex( GSL_REAL(buf), GSL_IMAG(buf) );
00137     }
00138   }
00139   
00140   fGotES = true;
00141   
00142 }


Member Data Documentation

double OscProb::PMNS_Transitions::fbmu [protected]

Definition at line 40 of file PMNS_Transitions.h.

Referenced by InitialiseTransitionParameters(), SetBMu(), and SolveHam().

Definition at line 41 of file PMNS_Transitions.h.

Referenced by InitialiseTransitionParameters(), SetBTau(), and SolveHam().

Definition at line 38 of file PMNS_Transitions.h.

Referenced by InitialiseTransitionParameters(), SetCMuMu(), and SolveHam().

Definition at line 39 of file PMNS_Transitions.h.

Referenced by InitialiseTransitionParameters(), SetCTauTau(), and SolveHam().


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

Generated on 16 Apr 2018 for loon by  doxygen 1.6.1