summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--src/domain_minimize.cpp12
-rw-r--r--src/domain_newton.h17
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++;
}