Simple model

.gitmodules

@@ -0,0 +1,3 @@
+[submodule "test/ndiff"]
+	path = test/ndiff
+	url = https://github.com/quinoacomputing/ndiff.git

CMakeLists.txt

@@ -10,7 +10,10 @@ endif(NOT CMAKE_BUILD_TYPE)
 
 project (DDESimulatorTemplate CXX C)
 set(DENSE_SOURCE_DIR "${PROJECT_SOURCE_DIR}")
-set(DENSE_SOURCE_DIR "${PROJECT_SOURCE_DIR}" PARENT_SCOPE)
+get_directory_property(hasParent PARENT_DIRECTORY)
+if(hasParent)
+  set(DENSE_SOURCE_DIR "${PROJECT_SOURCE_DIR}" PARENT_SCOPE)
+endif(hasParent)
 #Would prefer to use CMAKE_CXX_STANDARD in CMake 3.1, but not yet
 include(CheckCXXCompilerFlag)
 CHECK_CXX_COMPILER_FLAG("-std=c++11" COMPILER_SUPPORTS_CXX11)
@@ -31,7 +34,7 @@ endif()
 set(Real "double" CACHE STRING "dense::Real")
 set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DDENSE_CONFIGURABLE_Real=\"${Real}\"")
 
-find_package(CUDA 6.0 QUIET)
+#find_package(CUDA 6.0 QUIET)
 if (${CUDA_FOUND})
   SET(CUDA_PROPAGATE_HOST_FLAGS OFF)
   set(CUDA_HOST_COMPILER ${CMAKE_CXX_COMPILER})

source/CMakeLists.txt

@@ -72,16 +72,6 @@ function(STOCH_NR_SIMULATION localname simdir)
   SIMULATION(${localname} ${simdir} ${DENSE_SOURCE_DIR}/source/Next_Reaction_Main.cpp)
 endfunction(STOCH_NR_SIMULATION localname simdir)
 
-
-function(PARAMSEARCH localname simdir)
-  SIMCORELIB(${localname}_lib ${simdir})
-  add_executable(${localname} ${PROJECT_SOURCE_DIR}/source/param_search_main.cpp)
-  target_include_directories(${localname} PUBLIC ${simdir})
-  target_include_directories(${localname} PUBLIC ${PROJECT_SOURCE_DIR}/source)
-  target_link_libraries(${localname} PUBLIC ${localname}_lib)
-endfunction(PARAMSEARCH localname simdir)
-
-
 function(BENCHMARK localname simdir)
   SIMCORELIB(${localname}_lib ${simdir})
   add_executable(${localname} ${PROJECT_SOURCE_DIR}/source/Fast_Gillespie_Main.cpp)
@@ -130,18 +120,6 @@ function(CUDA_SIMULATION localname simdir)
   target_link_libraries(${localname} ${localname}_lib)
 endfunction(CUDA_SIMULATION localname simdir)
 
-function(PARAM_SEARCH_LIB localname simdir)
-  if (NOT TARGET ${localname})
-    add_library(${localname}
-                   ${DENSE_SOURCE_DIR}/source/search/libsres/ESES.cpp
-                   ${DENSE_SOURCE_DIR}/source/search/libsres/ESSRSort.cpp
-                   ${DENSE_SOURCE_DIR}/source/search/libsres/sharefunc.cpp
-                   ${DENSE_SOURCE_DIR}/source/search/sres.cpp)
-    target_include_directories(${localname} PUBLIC ${simdir})
-    target_include_directories(${localname} PUBLIC ${DENSE_SOURCE_DIR}/source)
-  endif (NOT TARGET ${localname})
-endfunction(PARAM_SEARCH_LIB localname simdir)
-
 
 if( ${CMAKE_CURRENT_BINARY_DIR} STREQUAL ${PROJECT_BINARY_DIR}/source)
   set(BUILD_TESTS_ONLY FALSE)
@@ -174,8 +152,6 @@ if( ${CMAKE_CURRENT_BINARY_DIR} STREQUAL ${PROJECT_BINARY_DIR}/source)
     set_target_properties(cell_growth PROPERTIES RUNTIME_OUTPUT_DIRECTORY "${PROJECT_BINARY_DIR}")
     STOCH_NR_SIMULATION(nextreaction ${PROJECT_BINARY_DIR})
     set_target_properties(nextreaction PROPERTIES RUNTIME_OUTPUT_DIRECTORY "${PROJECT_BINARY_DIR}")
-		  PARAMSEARCH(paramsearch ${PROJECT_BINARY_DIR})
-    set_target_properties(paramsearch PROPERTIES RUNTIME_OUTPUT_DIRECTORY "${PROJECT_BINARY_DIR}")
     if (${CUDA_FOUND})
       CUDA_SIMULATION(cuda_simulation ${PROJECT_BINARY_DIR})
       set_target_properties(cuda_simulation PROPERTIES RUNTIME_OUTPUT_DIRECTORY "${PROJECT_BINARY_DIR}")

source/Next_Reaction_Main.cpp

@@ -49,6 +49,7 @@ int main(int argc, char* argv[]){
   if(args.help == 2){
     return EXIT_FAILURE;
   }
+  int num_vertices = args.adj_graph.num_vertices();
   using Simulation = Next_Reaction_Simulation;
   auto num_vertices = args.adj_graph.num_vertices();
   Sim_Builder<Simulation> sim = Sim_Builder<Simulation>(args.perturbation_factors, args.gradient_factors, std::move(args.adj_graph), ac, av); 

source/measurement/basic.hpp

@@ -50,7 +50,7 @@ class BasicAnalysis : public Analysis<Simulation> {
 
   private:
 
-    std::vector<Real> mins, maxs, means;
+    std::vector<Real> mins, maxs, means, variance;
 
 		bool finalized;
 

source/measurement/basic.ipp

@@ -1,7 +1,9 @@
 #include <limits>
 #include <algorithm>
 #include <iostream>
+#include <string>
 #include <stdexcept>
+#include <cmath>
 
 template <typename Simulation>
 BasicAnalysis<Simulation>::BasicAnalysis (
@@ -13,6 +15,7 @@ BasicAnalysis<Simulation>::BasicAnalysis (
   mins(observed_species.size(), std::numeric_limits<Real>::infinity()),
   maxs(observed_species.size(), Real{0}),
   means(observed_species.size(), Real{0}),
+  variance(observed_species.size(), Real{0}),
   mins_by_context(Analysis<>::max - Analysis<>::min, mins),
   maxs_by_context(Analysis<>::max - Analysis<>::min, maxs),
   means_by_context(Analysis<>::max - Analysis<>::min, means) {
@@ -21,30 +24,40 @@ BasicAnalysis<Simulation>::BasicAnalysis (
 
 template <typename Simulation>
 void BasicAnalysis<Simulation>::update (Simulation& simulation, std::ostream&) {
-  for (Natural cell_no = this->min; cell_no < this->max; ++cell_no) {
-    for (std::size_t i = 0; i < this->observed_species_.size(); ++i) {
+  Natural n_min1 = Analysis<>::samples;
+  for (std::size_t i = 0; i < this->observed_species_.size(); ++i) {
+    Real concsum = 0.0;
+    for (Natural cell_no = this->min; cell_no < this->max; ++cell_no) {
   		Real concentration = simulation.get_concentration(cell_no, this->observed_species_[i]);
       mins[i] = std::min(concentration, mins[i]);
       maxs[i] = std::max(concentration, maxs[i]);
-      means[i] += concentration;
       mins_by_context[cell_no][i] = std::min(concentration, mins_by_context[cell_no][i]);
       maxs_by_context[cell_no][i] = std::max(concentration, maxs_by_context[cell_no][i]);
-  		means_by_context[cell_no][i] += concentration;
+      means_by_context[cell_no][i] = 
+        (means_by_context[cell_no][i]*n_min1 + concentration)/(n_min1 + 1);
+      concsum += concentration;
   	}
+
+    Real oldmean = means[i];
+    means[i] = (means[i]*n_min1 + (concsum/(this->max - this->min)))/(n_min1 + 1);
+    Real meandiff = oldmean - means[i];
+    Real sumsqdiff = 0.0f;
+    for (Natural cell_no = this->min; cell_no < this->max; ++cell_no) {
+      Real concdiff = simulation.get_concentration(cell_no, this->observed_species_[i]) - means[i];
+      sumsqdiff += concdiff*concdiff;
+    }
+    Real prev_samples = n_min1*(this->max - this->min);
+    variance[i] = ((prev_samples - 1)*variance[i] + prev_samples*meandiff*meandiff + sumsqdiff) 
+                                    / (prev_samples + (this->max - this->min) - 1 );
   }
   ++Analysis<>::samples;
 }
 
 template <typename Simulation>
 void BasicAnalysis<Simulation>::finalize () {
 
-if(!finalized){
-	finalized = true;
-}
-  for (auto & cell_means : means_by_context) {
-    for (auto & mean : cell_means) {
-      mean /= Analysis<>::samples;
-    }
+  if(!finalized){
+    finalized = true;
   }
 	detail.concs = means;
 	detail.other_details.emplace_back(mins);
@@ -63,13 +76,14 @@ void BasicAnalysis<Simulation>::show (csvw * csv_out) {
   Analysis<>::show(csv_out);
   if (csv_out) {
     auto & out = *csv_out;
-    out << "Species,Minimum Concentration,Mean Concentration,Maximum Concentration\n";
+    //out << "Species,Minimum Concentration,Mean Concentration,Maximum Concentration,Standard Deviation \n";
     //out << std::scientific << std::setprecision(5);
     for (specie_id species : Analysis<>::observed_species_) {
-      out << specie_str[species] << "," <<
-        mins[species] << "," <<
-        means[species] / (Analysis<>::samples * (Analysis<>::max - Analysis<>::min)) << "," <<
-        maxs[species] << "\n";
+      out << specie_str[species] << "\n" <<
+        "Minimum Concentration:\n" << mins[species] << "\n" <<
+        "Mean Concentration:\n" << means[species] << "\n" <<
+        "Maximum Concentration:\n" << maxs[species] << "\n" <<
+        "Standard Deviation:\n" << std::sqrt(variance[species]) << "\n" ;
     }
   }
-}
+}

test/Belousov-Zhabotinksy/anlys_basic_a.csv

@@ -1,5 +1,28 @@
-# Showing cells 0-99 until 14.5438 min
-Species,Minimum Concentration,Mean Concentration,Maximum Concentration
-FE2,0,0.722222,6
-FE3,0,0.757576,7
-BR,0,0.734848,6
+# Showing cells 0-99 until 14.5482 min
+FE2
+Minimum Concentration:
+0
+Mean Concentration:
+0.777778
+Maximum Concentration:
+8
+Standard Deviation:
+1.49974
+FE3
+Minimum Concentration:
+0
+Mean Concentration:
+0.762626
+Maximum Concentration:
+7
+Standard Deviation:
+1.18609
+BR
+Minimum Concentration:
+0
+Mean Concentration:
+0.825758
+Maximum Concentration:
+9
+Standard Deviation:
+1.56462

test/Belousov-Zhabotinksy/anlys_basic_b.csv

@@ -1,5 +1,28 @@
-# Showing cells 0-99 between 3 min and 14.5438 min
-Species,Minimum Concentration,Mean Concentration,Maximum Concentration
-FE2,0,1.45455,6
-FE3,0,3.03283,13
-BR,0,1.98485,9
+# Showing cells 0-99 between 3 min and 14.5482 min
+FE2
+Minimum Concentration:
+0
+Mean Concentration:
+1.19697
+Maximum Concentration:
+6
+Standard Deviation:
+1.62
+FE3
+Minimum Concentration:
+0
+Mean Concentration:
+3.34343
+Maximum Concentration:
+11
+Standard Deviation:
+2.17191
+BR
+Minimum Concentration:
+0
+Mean Concentration:
+2.05303
+Maximum Concentration:
+12
+Standard Deviation:
+2.45974

test/Belousov-Zhabotinksy/anlys_basic_c.csv

@@ -1,5 +1,28 @@
-# Showing cells 0-99 between 6 min and 14.5438 min
-Species,Minimum Concentration,Mean Concentration,Maximum Concentration
-FE2,0,1.48485,8
-FE3,0,4.7601,17
-BR,0,2.56566,14
+# Showing cells 0-99 between 6 min and 14.5482 min
+FE2
+Minimum Concentration:
+0
+Mean Concentration:
+1.42424
+Maximum Concentration:
+8
+Standard Deviation:
+1.79784
+FE3
+Minimum Concentration:
+0
+Mean Concentration:
+5.60859
+Maximum Concentration:
+21
+Standard Deviation:
+4.03225
+BR
+Minimum Concentration:
+0
+Mean Concentration:
+3
+Maximum Concentration:
+22
+Standard Deviation:
+3.52639

test/Belousov-Zhabotinksy/anlys_basic_d.csv

@@ -1,5 +1,28 @@
-# Showing cells 0-99 between 9 min and 14.5438 min
-Species,Minimum Concentration,Mean Concentration,Maximum Concentration
-FE2,0,1.33838,8
-FE3,0,6.30303,32
-BR,0,3.49495,16
+# Showing cells 0-99 between 9 min and 14.5482 min
+FE2
+Minimum Concentration:
+0
+Mean Concentration:
+1.32323
+Maximum Concentration:
+8
+Standard Deviation:
+1.75645
+FE3
+Minimum Concentration:
+0
+Mean Concentration:
+7.25758
+Maximum Concentration:
+35
+Standard Deviation:
+5.84434
+BR
+Minimum Concentration:
+0
+Mean Concentration:
+3.78283
+Maximum Concentration:
+28
+Standard Deviation:
+4.9935

test/Belousov-Zhabotinksy/anlys_basic_e.csv

@@ -1,5 +1,28 @@
-# Showing cells 0-99 between 12 min and 14.5438 min
-Species,Minimum Concentration,Mean Concentration,Maximum Concentration
-FE2,0,1.30808,8
-FE3,0,8.07071,42
-BR,0,4.82071,33
+# Showing cells 0-99 between 12 min and 14.5482 min
+FE2
+Minimum Concentration:
+0
+Mean Concentration:
+1.29293
+Maximum Concentration:
+8
+Standard Deviation:
+1.73863
+FE3
+Minimum Concentration:
+0
+Mean Concentration:
+9.90909
+Maximum Concentration:
+72
+Standard Deviation:
+10.898
+BR
+Minimum Concentration:
+0
+Mean Concentration:
+5.88384
+Maximum Concentration:
+62
+Standard Deviation:
+8.7651