Small refactoring

Author: mdimura

CMakeLists.txt

@@ -22,6 +22,7 @@ target_link_libraries(sparsehc-dm ${LIBS})
 add_definitions(-std=c++11)
 
 file(GLOB SHARED_SOURCES src/averagecluster.cpp src/cluster.cpp src/common.cpp src/completecluster.cpp src/dendrogram.cpp src/inmatrix.cpp src/matrix.cpp src/singlecluster.cpp)
-ADD_LIBRARY(sparsehc_dm SHARED sparsehc-dm_python.cpp ${SHARED_SOURCES})
+ADD_LIBRARY(sparsehc_dm SHARED python/sparsehc-dm_python.cpp ${SHARED_SOURCES})
 target_link_libraries(sparsehc_dm ${LIBS})
-set_target_properties(sparsehc_dm PROPERTIES PREFIX "")
+set_target_properties(sparsehc_dm PROPERTIES PREFIX "" SUFFIX ".so")
+

README

@@ -21,7 +21,7 @@ for i in range(0,N-1):
 Z=linkage(D, method='complete',preserve_input=False)
 ```
 ###Instalation
-####Prerequisits: boost graph and stxxl library
+####Prerequisites: boost graph and stxxl library
 ```
 sudo apt-get install libboost-graph-dev libstxxl-dev libstxxl1
 ```

README.md

@@ -1,16 +1,17 @@
 #!/bin/bash
-touch fastcluster_load.txt  sparsehc-dm_load.txt sparsehc-dm_load_pyRMSD.txt
+
 for Nk in 10 20 40 80
 do
-  /usr/bin/time -v python fastcluster_load.py $Nk 2>&1 | grep 'finished\|Maximum' >> fastcluster_load.txt 
+  /usr/bin/time -v python fastcluster_load.py $Nk 2>&1 | grep 'finished\|Maximum' >> fastcluster_perf.txt 
 done
 
-for Nk in 10 20 40 80 100 200 400
+for Nk in 10 20 40
 do
-  /usr/bin/time -v python sparsehc-dm_load.py $Nk 2>&1 | grep 'finished\|Maximum' >> sparsehc-dm_load.txt 
+ ./cpptraj_cluster.sh $Nk 2>&1 >> cpptraj_perf.txt 
 done
 
-for Nk in 10 20 40 80 100 200 400
+for Nk in 10 20 40 80 100 200 #400 700
 do
-  /usr/bin/time -v python sparsehc-dm_load_pyRMSD.py $Nk 2>&1 | grep 'finished\|Maximum' >> sparsehc-dm_load_pyRMSD.txt 
-done
+  /usr/bin/time -v python sparsehc-dm_load.py $Nk 2>&1 | grep 'finished\|Maximum' >> sparsehc-dm_perf.txt 
+done
+

examples/benchmark.sh

@@ -0,0 +1,17 @@
+#!/bin/bash
+Nk=$1
+PARMFILE=aMD-148l-all_1.pdb
+TRAJINFILE=aMD-148l-first${Nk}k.nc
+
+EPSILON=30.0
+
+cat <<EOF | /usr/bin/time -v cpptraj.OMP
+parm ${PARMFILE}
+trajin ${TRAJINFILE}
+
+cluster hieragglo epsilon $EPSILON complete rms :1-9999@CA out frame_clstr.dat summary clstr-summary.dat \
+info info.dat
+
+go
+EOF
+cd ..

examples/cpptraj_cluster.sh

@@ -4,6 +4,7 @@
 import random
 import time
 import mdtraj as md
+import numpy
 
 N=int(sys.argv[1])
 traj_filename='aMD-148l-first{}k.nc'.format(N)
@@ -14,15 +15,21 @@
 
 start=time.time()
 traj=md.load(traj_filename,top=top_filename, atom_indices=atoms_to_keep)
+
+m=sparsehc_dm.InMatrix()
 finishedLoad=time.time()
 print ("finished loading ({}k): {}".format(N,finishedLoad-start))
 
-m=sparsehc_dm.InMatrix()
 N=traj.n_frames
+rmsds=list()
 for i in range(0,N-1):
-  rmsds=md.rmsd(traj, traj, i)[i+1:]
+  rmsds=md.rmsd(traj, traj, i)[i+1:].tolist()
   sparsehc_dm.push(m,rmsds,i)
-    
+  #if( (i+1)%20==0):
+    #sparsehc_dm.push(m,rmsds,i-19,N)
+    #rmsds=list()
+
+
 #for i in range(0,N-1):
   #rmsds=md.rmsd(traj, traj, i)
   #for j in range(i+1,N):

examples/sparsehc-dm_load.py

@@ -0,0 +1,11 @@
+from distutils.core import setup
+from distutils.extension import Extension
+setup(name="sparsehc-dm",    
+      ext_modules=[
+                    Extension("sparsehc-dm",["sparsehc-dm_python.cpp",],
+                    library_dirs=["/usr/local/lib",],
+                    libraries=["",],
+                    include_dirs=[".","../src"],
+                    depends=[]),
+                    ]
+     )

python/setup.py

@@ -1,6 +1,7 @@
 #include <memory>
 #include <boost/python.hpp>
 #include <boost/python/numeric.hpp>
+#include <boost/python/stl_iterator.hpp>
 #include "common.h"
 #include "inmatrix.h"
 #include "averagecluster.h"
@@ -40,23 +41,37 @@ boost::python::list linkage(InMatrix& mat, const std::string& method) {
 
 }
 
-void pushList(InMatrix& mat,const boost::python::list& list, int i)
+/*TODO: improve pushList() performance by introducing a buffer to first copy values
+from Python::list and then push them to InMatrix in a separate thread.
+*/
+template<typename T>
+void pushList(InMatrix& mat,const T& list, unsigned i/*,const unsigned rowLength*/)
 {
-	Element e;
-	for(unsigned j=0; j<boost::python::len(list); ++j) {
-		float v=boost::python::extract<float>(list[j]);
-		e.update(i,j+i+1,v);
-		mat.push(e);
-	}
-}
-void pushNp(InMatrix& mat,const boost::python::numeric::array& list, int i)
-{
-	Element e;
-	for(unsigned j=0; j<boost::python::len(list); ++j) {
-		float v=boost::python::extract<float>(list[j]);
-		e.update(i,j+i+1,v);
-		mat.push(e);
+	/*const unsigned n=boost::python::len(list);
+	std::vector<float> vec;
+	vec.resize(n);*/
+	unsigned j=i+1;
+	std::for_each(boost::python::stl_input_iterator<float>(list),
+		       boost::python::stl_input_iterator<float>(),
+		       [&mat,&i,&j](const float& v) {
+				mat.push(Element(i,j++,v));
+			 });
+	/*std::copy(boost::python::stl_input_iterator<float>(list),
+		   boost::python::stl_input_iterator<float>(),vec.begin());
+
+	auto it=boost::python::stl_input_iterator<float>(list);
+	while(it!=boost::python::stl_input_iterator<float>()) {
+		for(unsigned j=i+1; j<rowLength; j++) {
+			//std::cout<<"i="<<i<<" j="<<j<<" n="<<n<<std::endl;
+			vec.emplace_back(i,j,*it++);
+		}
+		++i;
 	}
+	{ vec.emplace_back(i,i+1+j++,v);}
+	mat.push(Element(i,i+1+j++,v));
+	for(auto& e:vec) {
+		mat.push(std::move(e));
+	}*/
 }
 void (InMatrix::*push)(unsigned, unsigned, float) = &InMatrix::push;
 void (InMatrix::*pushElement)(Element) = &InMatrix::push;
@@ -69,6 +84,8 @@ BOOST_PYTHON_MODULE(sparsehc_dm)
 	class_<InMatrix,boost::noncopyable>("InMatrix")
 			.def("push",push);
 	def("linkage", linkage);
-	def("push", pushNp);
+	//boost::python::list
+	//boost::python::numeric::array
+	def("push", pushList<boost::python::list>);
 
 }

sparsehc-dm_python.cpp renamed to python/sparsehc-dm_python.cpp

@@ -47,9 +47,13 @@ struct InMatrix: public Matrix {
 		push(Element(row,col,value));
 	}
 	void push(Element e) {
-		threshold=threshold<e.value?nextafterf(e.value,e.value+100.0f):threshold;
+		threshold=threshold<e.value?e.value:threshold;
 		sorter.push(std::move(e));
 	}
+	/*void push(float v) {
+		threshold=threshold<v?v:threshold;
+		sorter.push(std::move(e));
+	}*/
 
 	/*void pushvec(std::vector<float> v,unsigned i) {
 		for(unsigned j=0; j<v.size(); ++j) {
@@ -60,6 +64,7 @@ struct InMatrix: public Matrix {
 		}
 	}*/
 	void sort() {
+		threshold=nextafterf(threshold,threshold+100.0f);
 		sorter.sort();
 		numPoints=(1L+sqrt(1L+sorter.size()*8L))/2L;
 	}