Allow calculating energy without including epitope bonus.

diff --git a/hamiltonian.cpp b/hamiltonian.cpp
index af10101..bcedac6 100644 (file)
--- a/hamiltonian.cpp
+++ b/hamiltonian.cpp
@@ -212,8 +212,6 @@ double Hamiltonian::get_energy(const std::set<unsigned int> &mutated_sites) cons
                 
     }
     
-    // Check for mismatch between targeted epitope and current sequence
-    
     return ((efield * bh) + (ecoupling * bJ));
     
 }
@@ -283,8 +281,23 @@ bool EpitopeHamiltonian::escaped(const std::set<unsigned int> &mutated_sites, co
 
 
 // Return the energy for a given set of mutated sites
+// default: do include contribution from epitopes
+
+double EpitopeHamiltonian::get_energy(const std::set<unsigned int> &mutated_sites) const
+{
+    return get_energy<true>(mutated_sites);
+}
+
+// tried to have a forward declaration of explicit specializations,
+// which "must occur before instantiation",
+// but then linker complained about undefined reference
+double ensure_get_energy_false_instantiated(const std::set<unsigned int> &x) { EpitopeHamiltonian A((std::string&)""); return A.get_energy<false>(x); }
+
+// Return the energy for a given set of mutated sites
 
-double EpitopeHamiltonian::get_energy(const std::set<unsigned int> &mutated_sites) const {
+template <bool include_epitope>
+double EpitopeHamiltonian::get_energy(const std::set<unsigned int> &mutated_sites) const
+{
 
     double efield = 0.0;
     double ecoupling = 0.0;
@@ -302,15 +315,19 @@ double EpitopeHamiltonian::get_energy(const std::set<unsigned int> &mutated_site
                 
     }
     
-    // Check for mismatch between targeted epitope and current sequence
-    
     double energy = (efield * bh) + (ecoupling * bJ);
 
-    for (unsigned ep=0; ep<penalty.size(); ++ep) {
-        if (!escaped(mutated_sites, epitopeWT[ep], epitopeMut[ep])) {
-            energy += penalty[ep] * bh;
+    if (include_epitope) {
+
+        // Check for mismatch between targeted epitope and current sequence
+        
+        for (unsigned ep=0; ep<penalty.size(); ++ep) {
+            if (!escaped(mutated_sites, epitopeWT[ep], epitopeMut[ep])) {
+                energy += penalty[ep] * bh;
+            }
+            // alternatively, if(escaped): energy -= penalty
         }
-        // alternatively, if(escaped): energy -= penalty
+
     }
 
     return energy;
@@ -318,6 +335,7 @@ double EpitopeHamiltonian::get_energy(const std::set<unsigned int> &mutated_site
 }
 
 
+
 /***********************************************
  * my simple 2-site 2-allele system
  ***********************************************/

diff --git a/hamiltonian.h b/hamiltonian.h
index 32ca26c..3d2faf4 100644 (file)
--- a/hamiltonian.h
+++ b/hamiltonian.h
@@ -56,6 +56,7 @@ public:
     bool escaped(const Virus &v, const std::vector<unsigned int> &eWT, const std::vector<unsigned int> &eMut);
     bool escaped(const std::set<unsigned int> &mutated_sites, const std::vector<unsigned int> &eWT, const std::vector<unsigned int> &eMut) const;
     double get_energy(const std::set<unsigned int> &mutated_sites) const;
+    template<bool include_epitope> double get_energy(const std::set<unsigned int> &mutated_sites) const;
 
 };