added a temporary peice of junk code to save the state of individual domain_newton interations

2 files changed, 20 insertions, 9 deletions

diff --git a/src/domain_minimize.cpp b/src/domain_minimize.cpp
index f316e4b..57f3923 100644
--- a/src/domain_minimize.cpp
+++ b/src/domain_minimize.cpp
@@ -63,7 +63,7 @@ int domain_minimize_naked(gsl_vector *z, unsigned n, double c, double eps, unsig
 	nakedgethess hess(c, n);
 	nakedgetenergy energy(c, n);
 
-	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, eta0, 0.1, 100, verb);
+	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, eta0, 0.1, 100, verb, false);
 }
 
 struct fixedgetgrad {
@@ -102,7 +102,7 @@ int domain_minimize_fixed(gsl_vector *z, unsigned n, double c, double eps, unsig
 	fixedgethess hess(c, n);
 	fixedgetenergy energy(c, n);
 
-	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, 0, 0, 0.1, 10, true);
+	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, 0, 0, 0.1, 10, true, false);
 }
 
 struct randgetgrad {
@@ -153,7 +153,7 @@ int domain_minimize_rand(gsl_vector *z, unsigned n, double c, unsigned ord, cons
 	randgethess hess(c, n, ord, k, a, phi);
 	randgetenergy energy(c, n, ord, k, a, phi);
 
-	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, bound, 0.1, 2, verb);
+	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, bound, 0.1, 2, verb, true);
 }
 
 struct nakedwellgetgrad {
@@ -198,7 +198,7 @@ int domain_minimize_nakedWell(gsl_vector *z, unsigned n, double c, double w, dou
 	nakedwellgethess hess(c, n, w, ss);
 	nakedwellgetenergy energy(c, n, w, ss);
 
-	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, eta0, 0.1, 100, verb);
+	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, eta0, 0.1, 100, verb, false);
 }
 
 
@@ -256,7 +256,7 @@ int domain_minimize_randWell(gsl_vector *z, unsigned n, double c, unsigned ord,
 	randwellgethess hess(c, n, ord, k, a, phi, w, ss);
 	randwellgetenergy energy(c, n, ord, k, a, phi, w, ss);
 
-	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, eta0, 0.1, 100, verb);
+	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, eta0, 0.1, 100, verb, false);
 }
 
 
@@ -296,5 +296,5 @@ int domain_minimize_fixedmin(gsl_vector *z, unsigned n, double c, double eps, un
 	fixedmingethess hess(c, n);
 	fixedmingetenergy energy(c, n);
 
-	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, bound, 0.1, 10, verb);
+	return domain_newton(z, size, params, energy, grad, hess, eps, N, beta, s, sigma, gamma, bound, 0.1, 10, verb, false);
 }
diff --git a/src/domain_newton.h b/src/domain_newton.h
index abb7ea3..3c806c0 100644
--- a/src/domain_newton.h
+++ b/src/domain_newton.h
@@ -2,6 +2,8 @@
 #define DOMAIN_NEWTON_H
 
 #include <math.h>
+#include <iostream>
+#include <string>
 
 // GSL includes.
 #include <gsl/gsl_math.h>
@@ -24,7 +26,8 @@ template <class energy_func, class grad_func, class hess_func>
 int domain_newton(gsl_vector *state, unsigned size, unsigned params,
 		energy_func get_energy, grad_func get_grad, hess_func get_hess, double
 		epsilon, unsigned max_iterations, double beta, double s, double sigma,
-		double gamma, double eta_0, double delta, double bound, bool verbose) {
+		double gamma, double eta_0, double delta, double bound, bool verbose, bool
+		save_states) {
 /* The function domain_newton carries out a modified version of Newton's
  * method.  On success, 0 is returned.  On failure, 1 is returned.
  * 
@@ -192,8 +195,8 @@ int domain_newton(gsl_vector *state, unsigned size, unsigned params,
 		if (fabs(norm - old_norm) < gamma * eta) eta *= delta;
 
 		// Prints several useful statistics for debugging purposes.
-		if (verbose) printf("m was %i, grad_norm was %e, eta was %e, energy was %e\n",
-				m, norm, eta, energy);
+		if (verbose) printf("NEWTON STEP %06d: m %i, grad_norm %e, eta %e, energy %e\n",
+				iterations, m, norm, eta, energy);
 
 		// Determines if the process has converged to acceptable precision.
 		if (norm < epsilon) {
@@ -207,6 +210,14 @@ int domain_newton(gsl_vector *state, unsigned size, unsigned params,
 		// Reset the norm for the next iteration.
 		old_norm = norm;
 
+		if (save_states) {
+			char str[40];
+			sprintf(str, "states/state-%06d.dat", iterations);
+			FILE *fout = fopen(str, "w");
+			gsl_vector_fprintf(fout, state, "%.15e");
+			fclose(fout);
+		}
+
 		// Increment the counter.
 		iterations++;
 	}