Allow concurrent use of multiple precisions.

Makefile

@@ -1,3 +1,6 @@
+# TODO: -arch doesn't work anymore?
+# TODO: -m32 needed on 64-bit machines, otherwise SWIG won't work?
+
 CXXFLAGS                 = -Wall
 # OUTPUT_OPTION            = -o .obj/$@
 swig_flags               = -Wall -c++ -python -outputtuple
@@ -107,21 +110,28 @@ install: libbct.a
 		mkdir $(install_dir)/include/bct/matlab; \
 	fi
 	if [[ "$(CXXFLAGS)" == *GSL_FLOAT* ]]; then \
-		cat bct_float.h bct.h > _bct.h; \
+		cat bct_float.h bct.h > $(install_dir)/include/bct/bct_float.h; \
 	elif [[ "$(CXXFLAGS)" == *GSL_DOUBLE* ]]; then \
-		cat bct_double.h bct.h > _bct.h; \
+		cat bct_double.h bct.h > $(install_dir)/include/bct/bct.h; \
 	elif [[ "$(CXXFLAGS)" == *GSL_LONG_DOUBLE* ]]; then \
-		cat bct_long_double.h bct.h > _bct.h; \
+		cat bct_long_double.h bct.h > $(install_dir)/include/bct/bct_long_double.h; \
 	fi
-	mv _bct.h $(install_dir)/include/bct/bct.h
 	cp matlab/matlab.h $(install_dir)/include/bct/matlab
 	cp matlab/sort.h $(install_dir)/include/bct/matlab
 	cp precision.h $(install_dir)/include/bct
-	cp libbct.a $(install_dir)/lib
+	if [[ "$(CXXFLAGS)" == *GSL_FLOAT* ]]; then \
+		cp libbct.a $(install_dir)/lib/libbct_float.a; \
+	elif [[ "$(CXXFLAGS)" == *GSL_DOUBLE* ]]; then \
+		cp libbct.a $(install_dir)/lib/libbct.a; \
+	elif [[ "$(CXXFLAGS)" == *GSL_LONG_DOUBLE* ]]; then \
+		cp libbct.a $(install_dir)/lib/libbct_long_double.a; \
+	fi
 
 uninstall:
 	-rm -rf $(install_dir)/include/bct
+	-rm $(install_dir)/lib/libbct_float.a
 	-rm $(install_dir)/lib/libbct.a
+	-rm $(install_dir)/lib/libbct_long_double.a
 
 clean:
 	-rm _bct_gsl.so _bct_py.so bct_gsl.py bct_py.py bct_gsl.pyc bct_py.pyc bct_gsl_wrap.cpp bct_py_wrap.cpp bct_gsl_wrap.o bct_py_wrap.o libbct.a $(objects)

assortativity.cpp

@@ -5,7 +5,7 @@ FP_T assortativity(const VECTOR_T*, const MATRIX_T*);
 /*
  * Computes assortativity for a directed graph.  Connection weights are ignored.
  */
-FP_T bct::assortativity_dir(const MATRIX_T* CIJ) {
+FP_T BCT_NAMESPACE::assortativity_dir(const MATRIX_T* CIJ) {
 	if (safe_mode) check_status(CIJ, SQUARE | DIRECTED, "assortativity_dir");
 
 	// [id,od,deg] = degrees_dir(CIJ);
@@ -26,7 +26,7 @@ FP_T bct::assortativity_dir(const MATRIX_T* CIJ) {
  * Computes assortativity for an undirected graph.  Connection weights are
  * ignored.
  */
-FP_T bct::assortativity_und(const MATRIX_T* CIJ) {
+FP_T BCT_NAMESPACE::assortativity_und(const MATRIX_T* CIJ) {
 	if (safe_mode) check_status(CIJ, SQUARE | UNDIRECTED, "assortativity_und");
 
 	// [deg] = degrees_und(m);
@@ -46,7 +46,7 @@ FP_T bct::assortativity_und(const MATRIX_T* CIJ) {
 }
 
 FP_T assortativity(const VECTOR_T* deg, const MATRIX_T* ij) {
-	using namespace bct;
+	using namespace BCT_NAMESPACE;
 
 	VECTOR_ID(const_view) i = MATRIX_ID(const_column)(ij, 0);
 	VECTOR_ID(const_view) j = MATRIX_ID(const_column)(ij, 1);

bct.h

@@ -1,15 +1,40 @@
-#ifndef BCT_H
+#include "precision.h"
+
+#undef SKIP
+
+#ifdef GSL_FLOAT
+#ifdef BCT_FLOAT_H
+#define SKIP
+#else
+#define BCT_FLOAT_H
+#endif
+#endif
+
+#ifdef GSL_DOUBLE
+#ifdef BCT_H
+#define SKIP
+#else
 #define BCT_H
+#endif
+#endif
 
-#include "precision.h"
+#ifdef GSL_LONG_DOUBLE
+#ifdef BCT_LONG_DOUBLE_H
+#define SKIP
+#else
+#define BCT_LONG_DOUBLE_H
+#endif
+#endif
+
+#ifndef SKIP
 
 #include <stdexcept>
 #include <vector>
 
 #include "matlab/matlab.h"
 
-namespace bct {
-	using namespace matlab;
+namespace BCT_NAMESPACE {
+	using namespace MATLAB_NAMESPACE;
 
 	class bct_exception : public std::runtime_error {
 	public:

betweenness_bin.cpp

@@ -3,7 +3,7 @@
 /*
  * Computes node betweenness for a binary graph.
  */
-VECTOR_T* bct::betweenness_bin(const MATRIX_T* G) {
+VECTOR_T* BCT_NAMESPACE::betweenness_bin(const MATRIX_T* G) {
 	VECTOR_T* BC;
 	MATRIX_T* EBC = edge_betweenness_bin(G, &BC);
 	MATRIX_ID(free)(EBC);
@@ -13,7 +13,7 @@ VECTOR_T* bct::betweenness_bin(const MATRIX_T* G) {
 /*
  * Computes node and edge betweenness for a binary graph.
  */
-MATRIX_T* bct::edge_betweenness_bin(const MATRIX_T* G, VECTOR_T** BC) {
+MATRIX_T* BCT_NAMESPACE::edge_betweenness_bin(const MATRIX_T* G, VECTOR_T** BC) {
 	if (safe_mode) check_status(G, SQUARE | BINARY, "edge_betweenness_bin");
 
 	// n=length(G);

betweenness_wei.cpp

@@ -5,7 +5,7 @@
 /*
  * Computes node betweenness for a weighted graph.
  */
-VECTOR_T* bct::betweenness_wei(const MATRIX_T* G) {
+VECTOR_T* BCT_NAMESPACE::betweenness_wei(const MATRIX_T* G) {
 	VECTOR_T* BC;
 	MATRIX_T* EBC = edge_betweenness_wei(G, &BC);
 	MATRIX_ID(free)(EBC);
@@ -15,7 +15,7 @@ VECTOR_T* bct::betweenness_wei(const MATRIX_T* G) {
 /*
  * Computes node and edge betweenness for a weighted graph.
  */
-MATRIX_T* bct::edge_betweenness_wei(const MATRIX_T* G, VECTOR_T** BC) {
+MATRIX_T* BCT_NAMESPACE::edge_betweenness_wei(const MATRIX_T* G, VECTOR_T** BC) {
 	if (safe_mode) check_status(G, SQUARE | WEIGHTED, "edge_betweenness_wei");
 
 	// n=length(G);

breadth.cpp

@@ -8,7 +8,7 @@
  * 0 precedes node i or that node i is unreachable.  Check distance(i) for
  * GSL_POSINF to differentiate between these two cases.
  */
-VECTOR_T* bct::breadth(const MATRIX_T* CIJ, int source, VECTOR_T** branch) {
+VECTOR_T* BCT_NAMESPACE::breadth(const MATRIX_T* CIJ, int source, VECTOR_T** branch) {
 	if (safe_mode) check_status(CIJ, SQUARE, "breadth");
 
 	// N = size(CIJ,1);

breadthdist.cpp

@@ -5,7 +5,7 @@
 /*
  * Computes reachability and distance matrices using breadth-first search.
  */
-MATRIX_T* bct::breadthdist(const MATRIX_T* CIJ, MATRIX_T** D) {
+MATRIX_T* BCT_NAMESPACE::breadthdist(const MATRIX_T* CIJ, MATRIX_T** D) {
 	if (safe_mode) check_status(CIJ, SQUARE, "breadthdist");
 
 	// N = size(CIJ,1);

cat.cpp

@@ -153,14 +153,14 @@ const FP_T cat_ctx[52 * 52] = {
 	0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 2, 1, 0, 1, 1, 0, 1, 2, 3, 2, 3, 2, 0
 };
 
-MATRIX_T* bct::get_cat_all() {
+MATRIX_T* BCT_NAMESPACE::get_cat_all() {
 	MATRIX_ID(const_view) mv = MATRIX_ID(const_view_array)(cat_all, 95, 95);
 	MATRIX_T* m = MATRIX_ID(alloc)(95, 95);
 	MATRIX_ID(memcpy)(m, &mv.matrix);
 	return m;
 }
 
-MATRIX_T* bct::get_cat_ctx() {
+MATRIX_T* BCT_NAMESPACE::get_cat_ctx() {
 	MATRIX_ID(const_view) mv = MATRIX_ID(const_view_array)(cat_ctx, 52, 52);
 	MATRIX_T* m = MATRIX_ID(alloc)(52, 52);
 	MATRIX_ID(memcpy)(m, &mv.matrix);

charpath.cpp

@@ -3,15 +3,15 @@
 #include "bct.h"
 
 /*
- * WARNING: bct::charpath_lambda takes a distance matrix, but
- * bct::capped_charpath_lambda takes a connection matrix.  Both should be
+ * WARNING: BCT_NAMESPACE::charpath_lambda takes a distance matrix, but
+ * BCT_NAMESPACE::capped_charpath_lambda takes a connection matrix.  Both should be
  * lengths, not weights (called distances in CalcMetric).
  */
 
 /*
  * Given a distance matrix, computes characteristic path length.
  */
-FP_T bct::charpath_lambda(const MATRIX_T* D) {
+FP_T BCT_NAMESPACE::charpath_lambda(const MATRIX_T* D) {
 	if (safe_mode) check_status(D, SQUARE, "charpath_lambda");
 
 	// lambda = sum(sum(D(D~=Inf)))/length(nonzeros(D~=Inf));
@@ -27,7 +27,7 @@ FP_T bct::charpath_lambda(const MATRIX_T* D) {
 /*
  * Given a connection matrix, computes capped characteristic path length.
  */
-FP_T bct::capped_charpath_lambda(const MATRIX_T* L) {
+FP_T BCT_NAMESPACE::capped_charpath_lambda(const MATRIX_T* L) {
 	if (safe_mode) check_status(L, SQUARE, "capped_charpath_lambda");
 	int N = L->size1;
 	int nonzeros = 0;
@@ -65,7 +65,7 @@ FP_T bct::capped_charpath_lambda(const MATRIX_T* L) {
 /*
  * Given a distance matrix, computes eccentricity, radius, and diameter.
  */
-VECTOR_T* bct::charpath_ecc(const MATRIX_T* D, FP_T* radius, FP_T* diameter) {
+VECTOR_T* BCT_NAMESPACE::charpath_ecc(const MATRIX_T* D, FP_T* radius, FP_T* diameter) {
 	if (safe_mode) check_status(D, SQUARE, "charpath_ecc");
 
 	// ecc = max(D.*(D~=Inf),[],2);

clustering_coef_bd.cpp

@@ -5,7 +5,7 @@
 /*
  * Computes clustering coefficient for a binary directed graph.
  */
-VECTOR_T* bct::clustering_coef_bd(const MATRIX_T* A) {
+VECTOR_T* BCT_NAMESPACE::clustering_coef_bd(const MATRIX_T* A) {
 	if (safe_mode) check_status(A, SQUARE | BINARY | DIRECTED, "clustering_coef_bd");
 
 	// S=A+A.';

clustering_coef_bu.cpp

@@ -3,7 +3,7 @@
 /*
  * Computes the clustering coefficient for a binary undirected graph.
  */
-VECTOR_T* bct::clustering_coef_bu(const MATRIX_T* G) {
+VECTOR_T* BCT_NAMESPACE::clustering_coef_bu(const MATRIX_T* G) {
 	if (safe_mode) check_status(G, SQUARE | BINARY | UNDIRECTED, "clustering_coef_bu");
 
 	// n=length(G);

clustering_coef_wd.cpp

@@ -5,7 +5,7 @@
 /*
  * Computes the clustering coefficient for a weighted directed graph.
  */
-VECTOR_T* bct::clustering_coef_wd(const MATRIX_T* W) {
+VECTOR_T* BCT_NAMESPACE::clustering_coef_wd(const MATRIX_T* W) {
 	if (safe_mode) check_status(W, SQUARE | WEIGHTED | DIRECTED, "clustering_coef_wd");
 
 	// A=W~=0;

clustering_coef_wu.cpp

@@ -5,7 +5,7 @@
 /*
  * Computes the clustering coefficient for a weighted undirected graph.
  */
-VECTOR_T* bct::clustering_coef_wu(const MATRIX_T* W) {
+VECTOR_T* BCT_NAMESPACE::clustering_coef_wu(const MATRIX_T* W) {
 	if (safe_mode) check_status(W, SQUARE | WEIGHTED | UNDIRECTED, "clustering_coef_wu");
 
 	// K=sum(W~=0,2);

connectivity_length.cpp

@@ -8,7 +8,7 @@
  * Marchiori and Latora (2000). Harmony in the small-world. Physica A 285:
  * 539-546.
  */
-FP_T bct::connectivity_length(const MATRIX_T* D) {
+FP_T BCT_NAMESPACE::connectivity_length(const MATRIX_T* D) {
 	if (safe_mode) check_status(D, SQUARE, "connectivity_length");
 	int N = D->size1;
 	FP_T sum = 0.0;

convert.cpp

@@ -5,7 +5,7 @@
 /*
  * Returns a copy of the given matrix with each nonzero element inverted.
  */
-MATRIX_T* bct::invert_elements(const MATRIX_T* m) {
+MATRIX_T* BCT_NAMESPACE::invert_elements(const MATRIX_T* m) {
 	MATRIX_T* inv_m = MATRIX_ID(alloc)(m->size1, m->size2);
 	for (int i = 0; i < (int)m->size1; i++) {
 		for (int j = 0; j < (int)m->size2; j++) {
@@ -23,7 +23,7 @@ MATRIX_T* bct::invert_elements(const MATRIX_T* m) {
 /*
  * Returns a copy of the given matrix with no loops.
  */
-MATRIX_T* bct::remove_loops(const MATRIX_T* m) {
+MATRIX_T* BCT_NAMESPACE::remove_loops(const MATRIX_T* m) {
 	MATRIX_T* nl_m = copy(m);
 	VECTOR_ID(view) diag_nl_m = MATRIX_ID(diagonal)(nl_m);
 	VECTOR_ID(set_zero)(&diag_nl_m.vector);
@@ -33,14 +33,14 @@ MATRIX_T* bct::remove_loops(const MATRIX_T* m) {
 /*
  * Returns a binary copy of the given matrix.
  */
-MATRIX_T* bct::to_binary(const MATRIX_T* m) {
+MATRIX_T* BCT_NAMESPACE::to_binary(const MATRIX_T* m) {
 	return compare_elements(m, fp_not_equal, 0.0);
 }
 
 /*
  * Returns a positive copy of the given matrix.
  */
-MATRIX_T* bct::to_positive(const MATRIX_T* m) {
+MATRIX_T* BCT_NAMESPACE::to_positive(const MATRIX_T* m) {
 	MATRIX_T* pos_m = MATRIX_ID(alloc)(m->size1, m->size2);
 	for (int i = 0; i < (int)m->size1; i++) {
 		for (int j = 0; j < (int)m->size2; j++) {
@@ -55,7 +55,7 @@ MATRIX_T* bct::to_positive(const MATRIX_T* m) {
  * nodes, if either m(i, j) or m(j, i) is nonzero, then both m(i, j) and m(j, i)
  * are set to one.  Otherwise, both are set to zero.
  */
-MATRIX_T* bct::to_undirected_bin(const MATRIX_T* m) {
+MATRIX_T* BCT_NAMESPACE::to_undirected_bin(const MATRIX_T* m) {
 	MATRIX_T* und_m = MATRIX_ID(calloc)(m->size1, m->size2);
 	for (int i = 0; i < (int)m->size1; i++) {
 		for (int j = i; j < (int)m->size2; j++) {
@@ -74,7 +74,7 @@ MATRIX_T* bct::to_undirected_bin(const MATRIX_T* m) {
  * Returns an undirected copy of the given weighted matrix.  For every pair of
  * nodes, m(i, j) and m(j, i) are both set to the average of their two values.
  */
-MATRIX_T* bct::to_undirected_wei(const MATRIX_T* m) {
+MATRIX_T* BCT_NAMESPACE::to_undirected_wei(const MATRIX_T* m) {
 	MATRIX_T* und_m = MATRIX_ID(calloc)(m->size1, m->size2);
 	for (int i = 0; i < (int)m->size1; i++) {
 		for (int j = i; j < (int)m->size2; j++) {

cycprob.cpp

@@ -3,7 +3,7 @@
 /*
  * Computes the fraction of all paths that are cycles.
  */
-VECTOR_T* bct::cycprob_fcyc(const std::vector<MATRIX_T*>& Pq) {
+VECTOR_T* BCT_NAMESPACE::cycprob_fcyc(const std::vector<MATRIX_T*>& Pq) {
 
 	// fcyc = zeros(1,size(Pq,3));
 	VECTOR_T* fcyc = zeros_vector(Pq.size());
@@ -31,7 +31,7 @@ VECTOR_T* bct::cycprob_fcyc(const std::vector<MATRIX_T*>& Pq) {
  * Computes the probability that a non-cyclic path of length (q - 1) can be
  * extended to form a cycle of length q.
  */
-VECTOR_T* bct::cycprob_pcyc(const std::vector<MATRIX_T*>& Pq) {
+VECTOR_T* BCT_NAMESPACE::cycprob_pcyc(const std::vector<MATRIX_T*>& Pq) {
 
 	// pcyc = zeros(1,size(Pq,3));
 	VECTOR_T* pcyc = zeros_vector(Pq.size());

debug.cpp

@@ -6,7 +6,7 @@
  * Prints a vector using the given format for each element.  This is only
  * provided for debugging purposes.  In other cases, use gsl_vector_fprintf.
  */
-void bct::printf(const VECTOR_T* v, const std::string& format) {
+void BCT_NAMESPACE::printf(const VECTOR_T* v, const std::string& format) {
 	for (int i = 0; i < (int)v->size; i++) {
 		std::printf(format.c_str(), VECTOR_ID(get)(v, i));
 		std::printf(" ");
@@ -18,7 +18,7 @@ void bct::printf(const VECTOR_T* v, const std::string& format) {
  * Prints a matrix using the given format for each element.  This is only
  * provided for debugging purposes.  In other cases, use gsl_matrix_fprintf.
  */
-void bct::printf(const MATRIX_T* m, const std::string& format) {
+void BCT_NAMESPACE::printf(const MATRIX_T* m, const std::string& format) {
 	for (int i = 0; i < (int)m->size1; i++) {
 		for (int j = 0; j < (int)m->size2; j++) {
 			std::printf(format.c_str(), MATRIX_ID(get)(m, i, j));
@@ -33,7 +33,7 @@ void bct::printf(const MATRIX_T* m, const std::string& format) {
  * provided for debugging purposes.  In other cases, use
  * gsl_permutation_fprintf.
  */
-void bct::printf(const gsl_permutation* p, const std::string& format) {
+void BCT_NAMESPACE::printf(const gsl_permutation* p, const std::string& format) {
 	for (int i = 0; i < (int)p->size; i++) {
 		std::printf(format.c_str(), gsl_permutation_get(p, i));
 		std::printf(" ");

degrees_dir.cpp

@@ -4,7 +4,7 @@
  * Computes degree, in-degree, and out-degree for a directed graph.  Connection
  * weights are ignored.
  */
-VECTOR_T* bct::degrees_dir(const MATRIX_T* CIJ, VECTOR_T** id, VECTOR_T** od) {
+VECTOR_T* BCT_NAMESPACE::degrees_dir(const MATRIX_T* CIJ, VECTOR_T** id, VECTOR_T** od) {
 	if (safe_mode) check_status(CIJ, SQUARE | DIRECTED, "degrees_dir");
 
 	// CIJ = FP_T(CIJ~=0);

degrees_und.cpp

@@ -3,7 +3,7 @@
 /*
  * Computes degree for an undirected graph.  Connection weights are ignored.
  */
-VECTOR_T* bct::degrees_und(const MATRIX_T* CIJ) {
+VECTOR_T* BCT_NAMESPACE::degrees_und(const MATRIX_T* CIJ) {
 	if (safe_mode) check_status(CIJ, SQUARE | UNDIRECTED, "degrees_und");
 
 	// CIJ = FP_T(CIJ~=0);

density_dir.cpp

@@ -3,7 +3,7 @@
 /*
  * Computes density for a directed graph.  Connection weights are ignored.
  */
-FP_T bct::density_dir(const MATRIX_T* CIJ) {
+FP_T BCT_NAMESPACE::density_dir(const MATRIX_T* CIJ) {
 	if (safe_mode) check_status(CIJ, SQUARE | DIRECTED, "density_dir");
 
 	// N = size(CIJ,1);

density_und.cpp

@@ -3,7 +3,7 @@
 /*
  * Computes density for an undirected graph.  Connection weights are ignored.
  */
-FP_T bct::density_und(const MATRIX_T* CIJ) {
+FP_T BCT_NAMESPACE::density_und(const MATRIX_T* CIJ) {
 	if (safe_mode) check_status(CIJ, SQUARE | UNDIRECTED, "density_und");
 
 	// N = size(CIJ,1);