From 1795e1ceed6463403f610ff33725f5608d322e68 Mon Sep 17 00:00:00 2001 From: "Gavin S. Davies" Date: Fri, 31 Jul 2026 17:31:21 +0200 Subject: [PATCH] Remove dead eigenvalue-solve machinery from OscCalcDecayEigen complexsolver() filled fEval with the Hamiltonian's eigenvalues, but fEval was never read anywhere (propagation uses fHam.exp() in PropagatePath instead), so every uncached energy point paid for a ComplexEigenSolver call whose result was discarded. The paired fGotES cache flag was never checked either. Removes both members, all their writers (constructors, SetIsNuBar, BuildHms), the complexsolver method itself, and the now-dead call site in SolveHam. Per @chenel's review on #48. --- OscLib/OscCalcDecayEigen.cxx | 35 ++++++----------------------------- OscLib/OscCalcDecayEigen.h | 13 ++++--------- 2 files changed, 10 insertions(+), 38 deletions(-) diff --git a/OscLib/OscCalcDecayEigen.cxx b/OscLib/OscCalcDecayEigen.cxx index 015758f..785eee9 100644 --- a/OscLib/OscCalcDecayEigen.cxx +++ b/OscLib/OscCalcDecayEigen.cxx @@ -27,7 +27,7 @@ namespace osc //--------------------------------------------------------------------------- - OscCalcDecayEigen::OscCalcDecayEigen() : fNumNus(3), fBuiltHms(false), fIsNuBar(false), fGotES(false) + OscCalcDecayEigen::OscCalcDecayEigen() : fNumNus(3), fBuiltHms(false), fIsNuBar(false) { this->SetStdPars(); this->ResetToFlavour(1); @@ -36,7 +36,7 @@ namespace osc } //--------------------------------------------------------------------------- - OscCalcDecayEigen::OscCalcDecayEigen(const OscCalcDecayEigen& calc) : fNumNus(3), fBuiltHms(false), fGotES(false) + OscCalcDecayEigen::OscCalcDecayEigen(const OscCalcDecayEigen& calc) : fNumNus(3), fBuiltHms(false) { fRho = calc.fRho; @@ -77,8 +77,7 @@ namespace osc fHd.setZero(); fBuffer.setZero(); fAlpha.setZero(); - fEval.setZero(); - + } void OscCalcDecayEigen::SaveTo(TDirectory* dir, const std::string& name) const @@ -128,7 +127,6 @@ namespace osc void OscCalcDecayEigen::SetIsNuBar(bool isNuBar) { - fGotES = fGotES && (fIsNuBar == isNuBar); fBuiltHms = fBuiltHms && (fIsNuBar == isNuBar); fIsNuBar = isNuBar; } @@ -360,11 +358,9 @@ namespace osc void OscCalcDecayEigen::BuildHms() { - // Check if anything changed + // Check if anything changed if(fBuiltHms) return; - - // Tag to recompute eigensystem - fGotES = false; + for(int j=0; j eigensolver; - eigensolver.compute(A); - for (int t = 0; t < w.size(); t++) { - w(t) = eigensolver.eigenvalues()(t).real(); - } - } - // + //--------------------------------------------------------------------------- //............................................................................. void OscCalcDecayEigen::SolveHam(double E, double Ne) @@ -445,11 +427,6 @@ namespace osc fHam(0, 0) += kr2GNe;} else{ fHam(0, 0) -= kr2GNe;} - - // Solve Hamiltonian for eigenvalues using the Eigen library method - complexsolver(fHam, fEval); - - fGotES = true; } diff --git a/OscLib/OscCalcDecayEigen.h b/OscLib/OscCalcDecayEigen.h index 35c0568..146fca9 100644 --- a/OscLib/OscCalcDecayEigen.h +++ b/OscLib/OscCalcDecayEigen.h @@ -92,10 +92,7 @@ namespace osc ///Rotate the Hamiltonian by theta_ij and delta_ij virtual void RotateH(int i, int j, Eigen::Matrix3cd& Ham); - /// Wrapper to solve non-hermitian matrix eigenvalues. - void complexsolver(const Eigen::Matrix3cd& A, Eigen::Vector3d& w); - - /// Build Hms = H*2E, where H is the Hamiltonian in vacuum on flavour basis + /// Build Hms = H*2E, where H is the Hamiltonian in vacuum on flavour basis /// and E is the neutrino energy. This is effectively the matrix of masses squared. virtual void BuildHms(); @@ -120,12 +117,10 @@ namespace osc double fCachedNe; ///Cached electron density double fCachedE; ///Cached neutrino energy bool fBuiltHms; ///Tag to avoid rebuilding Hms - bool fIsNuBar; /// anti-neutrino flag - bool fGotES; ///Tag to avoid recalculating eigensystem + bool fIsNuBar; /// anti-neutrino flag - Eigen::Vector3cd fBuffer; ///Buffer for neutrino state tranformations - Eigen::Vector3d fEval; ///Eigenvalues of the Hamiltonian - Eigen::Matrix3cd fHd; ///Decay hamiltonian + Eigen::Vector3cd fBuffer; ///Buffer for neutrino state tranformations + Eigen::Matrix3cd fHd; ///Decay hamiltonian Eigen::Vector3d fDm; ///m^2_i - m^2_1 in vacuum Eigen::Matrix3d fTheta; ///theta[i][j] mixing angle Eigen::Matrix3d fDelta; ///delta[i][j] CP violating phase