22#include <gsl/gsl_errno.h>
23#include <gsl/gsl_math.h>
24#include <gsl/gsl_min.h>
25#include <gsl/gsl_multimin.h>
37typedef void (*
func_d)(
const int *,
const double *,
double *,
const int *,
44static void my_fdf(
const gsl_vector *v,
void *params,
double *
f,
46 double *x, *gradient, ff2;
53 x = (
double *)malloc(dim *
sizeof(
double));
54 gradient = (
double *)malloc(dim *
sizeof(
double));
56 for (
i = 0;
i < dim;
i++)
57 x[
i] = gsl_vector_get(v,
i);
58 getgrad = (df == NULL) ? 0 : 1;
60 p->
func(&dim, x, &ff2, &getgrad, gradient);
66 for (
i = 0;
i < dim;
i++)
67 gsl_vector_set(df,
i, gradient[
i]);
74static double my_f(
const gsl_vector *v,
void *params) {
77 my_fdf(v, params, &val, NULL);
82static void my_df(
const gsl_vector *v,
void *params, gsl_vector *df) {
86typedef void (*
print_f_ptr)(
const int *,
const int *,
const double *,
87 const double *,
const double *,
const double *);
90 OCT_MINIMIZE)(
const int *method,
const int *dim,
double *point,
91 const double *step,
const double *line_tol,
92 const double *tolgrad,
const double *toldr,
97 double maxgrad, maxdr;
102 const gsl_multimin_fdfminimizer_type *T = NULL;
103 gsl_multimin_fdfminimizer *
s;
105 gsl_vector *absgrad, *absdr;
106 gsl_multimin_function_fdf my_func;
112 oldpoint = (
double *)malloc(*dim *
sizeof(
double));
113 grad = (
double *)malloc(*dim *
sizeof(
double));
119 my_func.params = (
void *)&p;
122 x = gsl_vector_alloc(*dim);
123 for (
i = 0;
i < *dim;
i++)
124 gsl_vector_set(x,
i, point[
i]);
127 absgrad = gsl_vector_alloc(*dim);
128 absdr = gsl_vector_alloc(*dim);
133 T = gsl_multimin_fdfminimizer_steepest_descent;
136 T = gsl_multimin_fdfminimizer_conjugate_fr;
139 T = gsl_multimin_fdfminimizer_conjugate_pr;
142 T = gsl_multimin_fdfminimizer_vector_bfgs;
145 T = gsl_multimin_fdfminimizer_vector_bfgs2;
149 s = gsl_multimin_fdfminimizer_alloc(T, *dim);
151 gsl_multimin_fdfminimizer_set(
s, &my_func, x, *step, *line_tol);
154 for (
i = 0;
i < *dim;
i++)
155 oldpoint[
i] = point[
i];
158 status = gsl_multimin_fdfminimizer_iterate(
s);
161 *
minimum = gsl_multimin_fdfminimizer_minimum(
s);
162 for (
i = 0;
i < *dim;
i++)
163 point[
i] = gsl_vector_get(gsl_multimin_fdfminimizer_x(
s),
i);
164 for (
i = 0;
i < *dim;
i++)
165 grad[
i] = gsl_vector_get(gsl_multimin_fdfminimizer_gradient(
s),
i);
168 for (
i = 0;
i < *dim;
i++)
169 gsl_vector_set(absdr,
i,
fabs(point[
i] - oldpoint[
i]));
170 maxdr = gsl_vector_max(absdr);
171 for (
i = 0;
i < *dim;
i++)
172 gsl_vector_set(absgrad,
i,
fabs(grad[
i]));
173 maxgrad = gsl_vector_max(absgrad);
176 write_info(&iter, dim,
minimum, &maxdr, &maxgrad, point);
179 for (
i = 0;
i < *dim;
i++)
180 oldpoint[
i] = point[
i];
185 if ((maxgrad <= *tolgrad) || (maxdr <= *toldr))
186 status = GSL_SUCCESS;
188 status = GSL_CONTINUE;
189 }
while (status == GSL_CONTINUE && iter <= *
maxiter);
191 if (status == GSL_CONTINUE)
194 gsl_multimin_fdfminimizer_free(
s);
196 gsl_vector_free(absgrad);
197 gsl_vector_free(absdr);
212typedef void (*
func1)(
const double *,
double *);
217double fn1(
double x,
void *params) {
224void FC_FUNC_(oct_1dminimize, OCT_1DMINIMIZE)(
double *a,
double *b,
double *m,
228 const gsl_min_fminimizer_type *T;
229 gsl_min_fminimizer *
s;
236 F.params = (
void *)&p;
238 T = gsl_min_fminimizer_brent;
239 s = gsl_min_fminimizer_alloc(T);
241 *status = gsl_min_fminimizer_set(
s, &F, *m, *a, *b);
243 gsl_set_error_handler_off();
247 *status = gsl_min_fminimizer_iterate(
s);
249 *m = gsl_min_fminimizer_x_minimum(
s);
250 *a = gsl_min_fminimizer_x_lower(
s);
251 *b = gsl_min_fminimizer_x_upper(
s);
253 *status = gsl_min_test_interval(*a, *b, 0.00001, 0.0);
257 }
while (*status == GSL_CONTINUE && iter < max_iter);
258 gsl_min_fminimizer_free(
s);
268typedef void (*
funcn)(
int *,
double *,
double *);
273double fn(
const gsl_vector *v,
void *params) {
281 x = (
double *)malloc(dim *
sizeof(
double));
283 for (
i = 0;
i < dim;
i++)
284 x[
i] = gsl_vector_get(v,
i);
285 p->
func(&dim, x, &val);
291typedef void (*
print_f_fn_ptr)(
const int *,
const int *,
const double *,
292 const double *,
const double *);
295 OCT_MINIMIZE_DIRECT)(
const int *method,
const int *dim,
296 double *point,
const double *step,
297 const double *toldr,
const int *
maxiter,
300 int iter = 0, status,
i;
303 const gsl_multimin_fminimizer_type *T = NULL;
304 gsl_multimin_fminimizer *
s = NULL;
306 gsl_multimin_function my_func;
313 my_func.params = (
void *)&p;
316 ss = gsl_vector_alloc(*dim);
317 gsl_vector_set_all(ss, *step);
320 x = gsl_vector_alloc(*dim);
321 for (
i = 0;
i < *dim;
i++)
322 gsl_vector_set(x,
i, point[
i]);
326 T = gsl_multimin_fminimizer_nmsimplex;
330 s = gsl_multimin_fminimizer_alloc(T, *dim);
331 gsl_multimin_fminimizer_set(
s, &my_func, x, ss);
335 status = gsl_multimin_fminimizer_iterate(
s);
341 for (
i = 0;
i < *dim;
i++)
342 point[
i] = gsl_vector_get(gsl_multimin_fminimizer_x(
s),
i);
344 size = gsl_multimin_fminimizer_size(
s);
345 status = gsl_multimin_test_size(size, *toldr);
347 write_info(&iter, dim,
minimum, &size, point);
349 }
while (status == GSL_CONTINUE && iter < *
maxiter);
351 if (status == GSL_CONTINUE)
356 gsl_multimin_fminimizer_free(
s);
real(real64) function func(r1, rn, n, a)
double fn(const gsl_vector *v, void *params)
static void my_fdf(const gsl_vector *v, void *params, double *f, gsl_vector *df)
double fn1(double x, void *params)
double fabs(double __x) __attribute__((__nothrow__
void(* print_f_fn_ptr)(const int *, const int *, const double *, const double *, const double *)
void(* func_d)(const int *, const double *, double *, const int *, double *)
void(* func1)(const double *, double *)
int FC_FUNC_(oct_minimize, OCT_MINIMIZE) const
static double my_f(const gsl_vector *v, void *params)
void(* funcn)(int *, double *, double *)
static void my_df(const gsl_vector *v, void *params, gsl_vector *df)
void(* print_f_ptr)(const int *, const int *, const double *, const double *, const double *, const double *)
integer, parameter, public minimum
real(real64) function s()
static double f(double w, void *p)