This is the mail archive of the
gsl-discuss@sources.redhat.com
mailing list for the GSL project.
Re: Additional multiroots functions.
Brian Gough <bjg@network-theory.co.uk> writes:
> It probably makes sense for alloc and set to become orthogonal, so
> that alloc only allocates memory, which would also solve your problem
> I think.
Please find below a patch which does this. I've also updated test.c
to reflect the change in semantics, and adjusted doc/multiroots.texi
to correctly describe the functions.
Thanks,
Toby
diff -ur gsl.orig/doc/multiroots.texi gsl/doc/multiroots.texi
--- gsl.orig/doc/multiroots.texi Sun May 14 21:10:13 2000
+++ gsl/doc/multiroots.texi Wed Feb 14 12:15:16 2001
@@ -94,15 +94,15 @@
@node Initializing the Multidimensional Solver
@section Initializing the Solver
-@deftypefun {gsl_multiroot_fsolver *} gsl_multiroot_fsolver_alloc (const gsl_multiroot_fsolver_type * @var{T}, gsl_multiroot_function * @var{f}, gsl_vector * @var{x})
+@deftypefun {gsl_multiroot_fsolver *} gsl_multiroot_fsolver_alloc (const gsl_multiroot_fsolver_type * @var{T}, size_t @var{n})
This function returns a pointer to a a newly allocated instance of a
-solver of type @var{T} for the function @var{f}, with an initial bracket
-on the root of @var{x}. For example, the following code creates an
-instance of a Brent solver,
+solver of type @var{T} for a system of @var{n} dimensions.
+For example, the following code creates an instance of a Brent solver,
+to solve a 3-dimensional system of equations.
@example
gsl_multiroot_fsolver * s =
- gsl_multiroot_fsolver_alloc (gsl_multiroot_fsolver_brent, f, x);
+ gsl_multiroot_fsolver_alloc (gsl_multiroot_fsolver_brent, 3);
@end example
If there is insufficient memory to create the solver then the function
@@ -110,15 +110,15 @@
code of @code{GSL_ENOMEM}.
@end deftypefun
-@deftypefun {gsl_multiroot_fdfsolver *} gsl_multiroot_fdfsolver_alloc (const gsl_multiroot_fdfsolver_type * @var{T}, gsl_multiroot_function_fdf * @var{fdf}, gsl_vector * @var{x})
+@deftypefun {gsl_multiroot_fdfsolver *} gsl_multiroot_fdfsolver_alloc (const gsl_multiroot_fdfsolver_type * @var{T}, size_t @var{n})
This function returns a pointer to a a newly allocated instance of a
-derivative solver of type @var{T} for the function @var{fdf}, with an
-initial guess for the root of @var{x}. For example, the following code
-creates an instance of a Newton-Raphson solver,
+derivative solver of type @var{T} for a system of size_t @var{n} dimensions.
+For example, the following code creates an instance of a Newton-Raphson solver,
+for a 2-dimensional system of equations.
@example
gsl_multiroot_fdfsolver * s =
- gsl_multiroot_fdfsolver_alloc (gsl_multiroot_fdfsolver_newton, fdf, x);
+ gsl_multiroot_fdfsolver_alloc (gsl_multiroot_fdfsolver_newton, 2);
@end example
If there is insufficient memory to create the solver then the function
@@ -127,12 +127,12 @@
@end deftypefun
@deftypefun int gsl_multiroot_fsolver_set (gsl_multiroot_fsolver * @var{s}, gsl_multiroot_function * @var{f}, gsl_vector * @var{x})
-This function reinitializes an existing solver @var{s} to use the
+This function sets, or resets, an existing solver @var{s} to use the
function @var{f} and the initial guess @var{x}.
@end deftypefun
@deftypefun int gsl_multiroot_fdfsolver_set (gsl_multiroot_fdfsolver * @var{s}, gsl_function_fdf * @var{fdf}, gsl_vector * @var{x})
-This function reinitializes an existing solver @var{s} to use the
+This function sets, or resets, an existing solver @var{s} to use the
function and derivative @var{fdf} and the initial guess @var{x}.
@end deftypefun
@@ -778,7 +778,8 @@
gsl_vector_set (x, 1, x_init[1]);
T = gsl_multiroot_fsolver_hybrids;
- s = gsl_multiroot_fsolver_alloc (T, &f, x);
+ s = gsl_multiroot_fsolver_alloc (T, 2);
+ gsl_multiroot_fsolver_set (s, &f, x);
print_state (iter, s);
@@ -919,7 +920,8 @@
gsl_vector_set (x, 1, x_init[1]);
T = gsl_multiroot_fdfsolver_gnewton;
- s = gsl_multiroot_fdfsolver_alloc (T, &f, x);
+ s = gsl_multiroot_fdfsolver_alloc (T, n);
+ gsl_multiroot_fdfsolver_alloc (s, &f, x);
print_state (iter, s);
diff -ur gsl.orig/multiroots/fdfsolver.c gsl/multiroots/fdfsolver.c
--- gsl.orig/multiroots/fdfsolver.c Mon Jun 5 20:27:37 2000
+++ gsl/multiroots/fdfsolver.c Wed Feb 14 12:07:12 2001
@@ -25,21 +25,12 @@
gsl_multiroot_fdfsolver *
gsl_multiroot_fdfsolver_alloc (const gsl_multiroot_fdfsolver_type * T,
- gsl_multiroot_function_fdf * f,
- gsl_vector * x)
+ size_t n)
{
int status;
gsl_multiroot_fdfsolver * s;
- const size_t n = f->n ;
-
- if (x->size != n)
- {
- GSL_ERROR_VAL ("vector length not compatible with function",
- GSL_EBADLEN, 0);
- }
-
s = (gsl_multiroot_fdfsolver *) malloc (sizeof (gsl_multiroot_fdfsolver));
if (s == 0)
@@ -116,21 +107,8 @@
GSL_ERROR_VAL ("failed to set solver", status, 0);
}
- status = gsl_multiroot_fdfsolver_set (s, f, x); /* seed the generator */
+ s->fdf = NULL;
- if (status != GSL_SUCCESS)
- {
- (s->type->free)(s->state);
- free (s->state);
- gsl_vector_free (s->dx);
- gsl_vector_free (s->x);
- gsl_vector_free (s->f);
- gsl_matrix_free (s->J);
- free (s); /* exception in constructor, avoid memory leak */
-
- GSL_ERROR_VAL ("failed to set solver", status, 0);
- }
-
return s;
}
@@ -139,6 +117,17 @@
gsl_multiroot_function_fdf * f,
gsl_vector * x)
{
+ if (s->x->size != f->n)
+ {
+ GSL_ERROR("function incompatible with solver size", GSL_EBADLEN);
+ }
+
+ if (x->size != f->n)
+ {
+ GSL_ERROR_VAL ("vector length not compatible with function",
+ GSL_EBADLEN, 0);
+ }
+
s->fdf = f;
gsl_vector_memcpy(s->x,x);
@@ -148,6 +137,11 @@
int
gsl_multiroot_fdfsolver_iterate (gsl_multiroot_fdfsolver * s)
{
+ if (!s->fdf)
+ {
+ GSL_ERROR ("fdfsolver not set", GSL_EINVAL);
+ }
+
return (s->type->iterate) (s->state, s->fdf, s->x, s->f, s->J, s->dx);
}
diff -ur gsl.orig/multiroots/fsolver.c gsl/multiroots/fsolver.c
--- gsl.orig/multiroots/fsolver.c Mon Jun 5 20:27:37 2000
+++ gsl/multiroots/fsolver.c Wed Feb 14 12:06:49 2001
@@ -23,22 +23,14 @@
#include <gsl/gsl_errno.h>
#include <gsl/gsl_multiroots.h>
-gsl_multiroot_fsolver *
-gsl_multiroot_fsolver_alloc (const gsl_multiroot_fsolver_type * T,
- gsl_multiroot_function * f, gsl_vector * x)
+gsl_multiroot_fsolver *
+gsl_multiroot_fsolver_alloc (const gsl_multiroot_fsolver_type * T,
+ size_t n)
{
int status;
gsl_multiroot_fsolver * s;
- const size_t n = f->n ;
-
- if (x->size != n)
- {
- GSL_ERROR_VAL ("vector length not compatible with function",
- GSL_EBADLEN, 0);
- }
-
s = (gsl_multiroot_fsolver *) malloc (sizeof (gsl_multiroot_fsolver));
if (s == 0)
@@ -102,20 +94,9 @@
GSL_ERROR_VAL ("failed to set solver", status, 0);
}
-
- status = gsl_multiroot_fsolver_set (s, f, x); /* seed the generator */
-
- if (status != GSL_SUCCESS)
- {
- free (s->state);
- gsl_vector_free (s->dx);
- gsl_vector_free (s->x);
- gsl_vector_free (s->f);
- free (s); /* exception in constructor, avoid memory leak */
-
- GSL_ERROR_VAL ("failed to set solver", status, 0);
- }
+ s->function = NULL;
+
return s;
}
@@ -124,6 +105,17 @@
gsl_multiroot_function * f,
gsl_vector * x)
{
+ if (s->x->size != f->n)
+ {
+ GSL_ERROR("function incompatible with solver size", GSL_EBADLEN);
+ }
+
+ if (x->size != f->n)
+ {
+ GSL_ERROR_VAL ("vector length not compatible with function",
+ GSL_EBADLEN, 0);
+ }
+
s->function = f;
gsl_vector_memcpy(s->x,x);
@@ -133,6 +125,11 @@
int
gsl_multiroot_fsolver_iterate (gsl_multiroot_fsolver * s)
{
+ if (!s->function)
+ {
+ GSL_ERROR ("fsolver not set", GSL_EINVAL);
+ }
+
return (s->type->iterate) (s->state, s->function, s->x, s->f, s->dx);
}
diff -ur gsl.orig/multiroots/gsl_multiroots.h gsl/multiroots/gsl_multiroots.h
--- gsl.orig/multiroots/gsl_multiroots.h Thu May 4 12:25:05 2000
+++ gsl/multiroots/gsl_multiroots.h Wed Feb 14 12:05:53 2001
@@ -77,9 +77,10 @@
}
gsl_multiroot_fsolver;
-gsl_multiroot_fsolver *
+gsl_multiroot_fsolver *
gsl_multiroot_fsolver_alloc (const gsl_multiroot_fsolver_type * T,
- gsl_multiroot_function * f, gsl_vector * x);
+ size_t n);
+
void gsl_multiroot_fsolver_free (gsl_multiroot_fsolver * s);
int gsl_multiroot_fsolver_set (gsl_multiroot_fsolver * s,
@@ -131,11 +132,9 @@
}
gsl_multiroot_fdfsolver;
-
gsl_multiroot_fdfsolver *
-gsl_multiroot_fdfsolver_alloc (const gsl_multiroot_fdfsolver_type * T,
- gsl_multiroot_function_fdf * fdf,
- gsl_vector * x);
+gsl_multiroot_fdfsolver_alloc (const gsl_multiroot_fdfsolver_type * T,
+ size_t n);
int
gsl_multiroot_fdfsolver_set (gsl_multiroot_fdfsolver * s,
diff -ur gsl.orig/multiroots/test.c gsl/multiroots/test.c
--- gsl.orig/multiroots/test.c Thu May 4 12:25:05 2000
+++ gsl/multiroots/test.c Wed Feb 14 12:07:58 2001
@@ -134,7 +134,8 @@
if (factor != 1.0) scale(x, factor);
- s = gsl_multiroot_fdfsolver_alloc (T, function, x);
+ s = gsl_multiroot_fdfsolver_alloc (T, n);
+ gsl_multiroot_fdfsolver_set (s, function, x);
do
{
@@ -220,7 +221,8 @@
if (factor != 1.0) scale(x, factor);
- s = gsl_multiroot_fsolver_alloc (T, &function, x);
+ s = gsl_multiroot_fsolver_alloc (T, n);
+ gsl_multiroot_fsolver_set (s, &function, x);
/* printf("x "); gsl_vector_fprintf (stdout, s->x, "%g"); printf("\n"); */
/* printf("f "); gsl_vector_fprintf (stdout, s->f, "%g"); printf("\n"); */
--
Toby White, University Chemical Lab., Lensfield Road, Cambridge. CB2 1EW. U.K.
Email: <tow@theor.ch.cam.ac.uk> GPG Key ID: 1DE9DE75
Web: <URL:http://ket.ch.cam.ac.uk/people/tow/index.html>
Tel: +44 1223 336423
Fax: +44 1223 336362