From c0b667c019e657f5347a388110c3d619cfbc0300 Mon Sep 17 00:00:00 2001
From: Dariusz Murakowski <murakdar@mit.edu>
Date: Fri, 1 May 2015 00:18:36 -0400
Subject: [PATCH] Copy constructor for Virus and NTVirus object to reduce
 unnecessary energy calculations.

---
 pop_ss.cpp   | 36 +++++++++++++++++++++++++++++-------
 pop_ss.h     |  2 ++
 reaction.cpp | 13 ++++++++-----
 3 files changed, 39 insertions(+), 12 deletions(-)

diff --git a/pop_ss.cpp b/pop_ss.cpp
index 2556ead..69ebd18 100644
--- a/pop_ss.cpp
+++ b/pop_ss.cpp
@@ -11,17 +11,28 @@
 #include "seqTools.h"
 
 
+// copy existing virus
+Virus::Virus(const Virus &v)
+: energy(v.energy)
+, mutated_sites(v.mutated_sites)
+, L(v.L)
+{
+}
+
 //Constuct a virus object of wildtype
-Virus::Virus(const Hamiltonian &H) {
-    this->mutated_sites.clear();
-    L=H.size;
-    energy=0;
+Virus::Virus(const Hamiltonian &H)
+: energy(0.0)
+, mutated_sites()
+, L(H.size)
+{
 }
 
 //Construct a virus and compute its energy.
-Virus::Virus(const Hamiltonian &H, const std::set<unsigned int> &mutated_sites) {
-    this->mutated_sites=mutated_sites;
-    L=H.size;
+Virus::Virus(const Hamiltonian &H, const std::set<unsigned int> &mutated_sites_in)
+: energy(0.0)
+, mutated_sites(mutated_sites_in)
+, L(H.size)
+{
     update_energy(H);
 }
 
@@ -56,6 +67,17 @@ std::string aaseq2str(const std::vector<aa> &aa_seq)
     return out;
 }
 
+// copy existing virus
+NTVirus::NTVirus(const NTVirus &v)
+: energy(v.energy)
+, nt_seq(v.nt_seq)
+, aa_seq(v.aa_seq)
+, config(v.config)
+, L_nt(v.L_nt)
+, L_aa(v.L_aa)
+{
+}
+
 NTVirus::NTVirus(const std::string &NT)
 : energy(0.0)
 , nt_seq(NT)
diff --git a/pop_ss.h b/pop_ss.h
index 221276f..4312b48 100644
--- a/pop_ss.h
+++ b/pop_ss.h
@@ -23,6 +23,7 @@ public:
     double energy;
     std::set<unsigned int> mutated_sites;
     
+    Virus(const Virus &v);
     Virus(const Hamiltonian &H);
     Virus(const Hamiltonian &H, const std::set<unsigned int> &mutated_sites);
     
@@ -57,6 +58,7 @@ public:
     unsigned int L_nt;
     unsigned int L_aa;
 
+    NTVirus(const NTVirus &v);
     NTVirus(const std::string &NT);
     NTVirus(const PottsHamiltonian &H, const std::string &NT);
     void update_energy(const PottsHamiltonian &H);
diff --git a/reaction.cpp b/reaction.cpp
index 4fc010f..f77aecc 100644
--- a/reaction.cpp
+++ b/reaction.cpp
@@ -63,8 +63,8 @@ double RC_fun(double E)
     // sigmoid function
     //double z = exp(-E);
     //return z/(1.0+z);
-    //return (E>0.0) ? 1.0/(1.0+exp(-E)) : exp(E)/(1.0+exp(E));
-    return exp(-E);
+    return (E>0.0) ? 1.0/(1.0+exp(-E)) : exp(E)/(1.0+exp(E));
+    //return exp(-E);
 }
 
 double VirusReaction::recalc()
@@ -99,7 +99,8 @@ void VirusReaction::fire()
 
     // pick how many sites to mutate in the new virus
     // adapted from WF Virus::mutate()
-    Virus v(H,iter->first.mutated_sites);   // TODO: replace with copy constructor
+    //Virus v(H,iter->first.mutated_sites);   // TODO: replace with copy constructor
+    Virus v(iter->first);
     unsigned int n = gsl_ran_binomial(rnd,mu,H.size);
     MutatedSiteSequence sites_to_mutate;
     while (sites_to_mutate.size() < n) {
@@ -120,7 +121,8 @@ void VirusReaction::fire()
         else                             v.mutated_sites.erase(x);
     }
     */
-    v.update_energy(H);
+    if (n > 0)  // sequence changed, so energy might have too
+        v.update_energy(H);
 
     // update copy numbers
     if (V->pop.count(v)==0) V->pop[v]  = 1;  // new variant
@@ -267,7 +269,8 @@ void NTVirusReaction::fire()
 
     // pick how many sites (nucleotides) to mutate in the new virus
     // adapted from WF Virus::mutate()
-    NTVirus v(H,iter->first.nt_seq);    // TODO: replace with copy constructor
+    //NTVirus v(H,iter->first.nt_seq);    // TODO: replace with copy constructor
+    NTVirus v(iter->first);
     unsigned int n = gsl_ran_binomial(rnd,mu,3*H.size);
     MutatedSiteSequence sites_to_mutate;
     while (sites_to_mutate.size() < n) {
-- 
2.7.4