Hi Brian,
Attached is the patch. As there are many parameters for miser and vegas,
I added a structure which holds all the parameters and added functions
which copy the values in and out of the structure. Otherwise we would
need over 10 getter/setter functions for each method.
Let me know if you like this approach, I'll update the manual then.
All the best,
Szymon
On Thu, 2009-07-09 at 21:06 +0100, Brian Gough wrote:
At Wed, 8 Jul 2009 06:54:17 -0400 (EDT),
Szymon Jaroszewicz wrote:
> Should the API be updated to allow setting at least some of the
> parameters without messing with state fields directly? At the
> minimum I think there should be a way to get the chisq field for
> Vegas which is used to test convergence.
This makes sense, would you like to send me a patch for the functions
you want.
diff --git a/monte/demo.c b/monte/demo.c
index ed41151..6b58e6b 100644
--- a/monte/demo.c
+++ b/monte/demo.c
@@ -81,9 +81,9 @@ main ()
{
gsl_monte_vegas_integrate (&G, xl, xu, 3, calls, r, s, &res, &err);
printf ("result = % .6f sigma = % .6f chisq/dof = %.1f\n",
- res, err, s->chisq);
+ res, err, gsl_monte_vegas_chisq (s));
}
- while (fabs (s->chisq - 1.0) > 0.5);
+ while (fabs (gsl_monte_vegas_chisq (s) - 1.0) > 0.5);
display_results ("vegas final", res, err);
diff --git a/monte/gsl_monte_miser.h b/monte/gsl_monte_miser.h
index b053063..76dcb3f 100644
--- a/monte/gsl_monte_miser.h
+++ b/monte/gsl_monte_miser.h
@@ -77,6 +77,17 @@ int gsl_monte_miser_init(gsl_monte_miser_state* state);
void gsl_monte_miser_free(gsl_monte_miser_state* state);
+typedef struct {
+ double estimate_frac;
+ size_t min_calls;
+ size_t min_calls_per_bisection;
+ double alpha;
+ double dither;
+} gsl_monte_miser_params;
+
+void gsl_monte_miser_get_params(const gsl_monte_miser_state* state,
gsl_monte_miser_params* params);
+
+void gsl_monte_miser_set_params(gsl_monte_miser_state* state, const
gsl_monte_miser_params* params);
__END_DECLS
diff --git a/monte/gsl_monte_vegas.h b/monte/gsl_monte_vegas.h
index 8ecef46..35b8090 100644
--- a/monte/gsl_monte_vegas.h
+++ b/monte/gsl_monte_vegas.h
@@ -22,6 +22,7 @@
#ifndef __GSL_MONTE_VEGAS_H__
#define __GSL_MONTE_VEGAS_H__
+#include <stdlib.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_monte.h>
@@ -99,6 +100,21 @@ int gsl_monte_vegas_init(gsl_monte_vegas_state* state);
void gsl_monte_vegas_free (gsl_monte_vegas_state* state);
+double gsl_monte_vegas_chisq (const gsl_monte_vegas_state* state);
+
+typedef struct {
+ double alpha;
+ size_t iterations;
+ int stage;
+ int mode;
+ int verbose;
+ FILE * ostream;
+} gsl_monte_vegas_params;
+
+void gsl_monte_vegas_get_params(const gsl_monte_vegas_state* state,
gsl_monte_vegas_params* params);
+
+void gsl_monte_vegas_set_params(gsl_monte_vegas_state* state, const
gsl_monte_vegas_params* params);
+
__END_DECLS
#endif /* __GSL_MONTE_VEGAS_H__ */
diff --git a/monte/miser.c b/monte/miser.c
index f7b29c3..c468a64 100644
--- a/monte/miser.c
+++ b/monte/miser.c
@@ -567,6 +567,27 @@ gsl_monte_miser_free (gsl_monte_miser_state * s)
free (s);
}
+void
+gsl_monte_miser_get_params(const gsl_monte_miser_state* s,
gsl_monte_miser_params* p)
+{
+ p->estimate_frac = s->estimate_frac;
+ p->min_calls = s->min_calls;
+ p->min_calls_per_bisection = s->min_calls_per_bisection;
+ p->alpha = s->alpha;
+ p->dither = s->dither;
+}
+
+void
+gsl_monte_miser_set_params(gsl_monte_miser_state* s, const
gsl_monte_miser_params* p)
+{
+ s->estimate_frac = p->estimate_frac;
+ s->min_calls = p->min_calls;
+ s->min_calls_per_bisection = p->min_calls_per_bisection;
+ s->alpha = p->alpha;
+ s->dither = p->dither;
+}
+
+
static int
estimate_corrmc (gsl_monte_function * f,
const double xl[], const double xu[],
diff --git a/monte/test.c b/monte/test.c
index c1dd0a7..72d2810 100644
--- a/monte/test.c
+++ b/monte/test.c
@@ -282,6 +282,26 @@ main (void)
#undef MONTE_SPEEDUP
#endif
+#ifdef MISER
+#define NAME "miser(params)"
+#define MONTE_STATE gsl_monte_miser_state
+#define MONTE_ALLOC gsl_monte_miser_alloc
+#define MONTE_PARAMS gsl_monte_miser_params
+#define MONTE_INTEGRATE(f,xl,xu,dim,calls,r,s,res,err) {
gsl_monte_miser_get_params(s, ¶ms) ; params.alpha = 1.5 ;
gsl_monte_miser_set_params(s, ¶ms) ;
gsl_monte_miser_integrate(f,xl,xu,dim,calls,r,s,res,err); }
+#define MONTE_FREE gsl_monte_miser_free
+#define MONTE_SPEEDUP 2
+#define MONTE_ERROR_TEST(err,expected) gsl_test(err > 5.0 * expected, NAME ",
%s, abserr[%d] (obs %g vs plain %g)", I->description, i, err, expected)
+#include "test_main.c"
+#undef NAME
+#undef MONTE_STATE
+#undef MONTE_ALLOC
+#undef MONTE_PARAMS
+#undef MONTE_INTEGRATE
+#undef MONTE_FREE
+#undef MONTE_ERROR_TEST
+#undef MONTE_SPEEDUP
+#endif
+
#ifdef VEGAS
#define NAME "vegas"
#define MONTE_STATE gsl_monte_vegas_state
@@ -289,7 +309,7 @@ main (void)
#define MONTE_INTEGRATE(f,xl,xu,dim,calls,r,s,res,err) {
gsl_monte_vegas_integrate(f,xl,xu,dim,calls,r,s,res,err) ; }
#define MONTE_FREE gsl_monte_vegas_free
#define MONTE_SPEEDUP 3
-#define MONTE_ERROR_TEST(err,expected) gsl_test(err > 3.0 * (expected == 0 ?
1.0/(I->calls/MONTE_SPEEDUP) : expected), NAME ", %s, abserr[%d] (obs %g vs exp
%g)", I->description, i, err, expected) ; gsl_test(s->chisq < 0, NAME " returns
valid chisq (%g)", s->chisq)
+#define MONTE_ERROR_TEST(err,expected) gsl_test(err > 3.0 * (expected == 0 ?
1.0/(I->calls/MONTE_SPEEDUP) : expected), NAME ", %s, abserr[%d] (obs %g vs exp
%g)", I->description, i, err, expected) ; gsl_test(gsl_monte_vegas_chisq(s) <
0, NAME " returns valid chisq (%g)", gsl_monte_vegas_chisq(s))
#include "test_main.c"
#undef NAME
#undef MONTE_STATE
@@ -308,11 +328,32 @@ main (void)
#define MONTE_INTEGRATE(f,xl,xu,dim,calls,r,s,res,err) {
gsl_monte_vegas_integrate(f,xl,xu,dim,calls,r,s,res,err) ;
gsl_monte_vegas_integrate(f,xl,xu,dim,calls,r,s,res,err); }
#define MONTE_FREE gsl_monte_vegas_free
#define MONTE_SPEEDUP 3
-#define MONTE_ERROR_TEST(err,expected) gsl_test(err > 3.0 * (expected == 0 ?
1.0/(I->calls/MONTE_SPEEDUP) : expected), NAME ", %s, abserr[%d] (obs %g vs exp
%g)", I->description, i, err, expected); gsl_test(s->chisq < 0, NAME " returns
valid chisq (%g)", s->chisq)
+#define MONTE_ERROR_TEST(err,expected) gsl_test(err > 3.0 * (expected == 0 ?
1.0/(I->calls/MONTE_SPEEDUP) : expected), NAME ", %s, abserr[%d] (obs %g vs exp
%g)", I->description, i, err, expected); gsl_test(gsl_monte_vegas_chisq(s) < 0,
NAME " returns valid chisq (%g)", gsl_monte_vegas_chisq(s))
+#include "test_main.c"
+#undef NAME
+#undef MONTE_STATE
+#undef MONTE_ALLOC
+#undef MONTE_INTEGRATE
+#undef MONTE_FREE
+#undef MONTE_ERROR_TEST
+#undef MONTE_SPEEDUP
+#endif
+
+
+#ifdef VEGAS
+#define NAME "vegas(params)"
+#define MONTE_STATE gsl_monte_vegas_state
+#define MONTE_ALLOC gsl_monte_vegas_alloc
+#define MONTE_PARAMS gsl_monte_vegas_params
+#define MONTE_INTEGRATE(f,xl,xu,dim,calls,r,s,res,err) {
gsl_monte_vegas_get_params(s, ¶ms) ; params.alpha = 2 ; params.iterations =
3 ; gsl_monte_vegas_set_params(s, ¶ms) ;
gsl_monte_vegas_integrate(f,xl,xu,dim,calls,r,s,res,err); }
+#define MONTE_FREE gsl_monte_vegas_free
+#define MONTE_SPEEDUP 3
+#define MONTE_ERROR_TEST(err,expected) gsl_test(err > 3.0 * (expected == 0 ?
1.0/(I->calls/MONTE_SPEEDUP) : expected), NAME ", %s, abserr[%d] (obs %g vs exp
%g)", I->description, i, err, expected); gsl_test(gsl_monte_vegas_chisq(s) < 0,
NAME " returns valid chisq (%g)", gsl_monte_vegas_chisq(s))
#include "test_main.c"
#undef NAME
#undef MONTE_STATE
#undef MONTE_ALLOC
+#undef MONTE_PARAMS
#undef MONTE_INTEGRATE
#undef MONTE_FREE
#undef MONTE_ERROR_TEST
diff --git a/monte/test_main.c b/monte/test_main.c
index a797e8a..7cd71d1 100644
--- a/monte/test_main.c
+++ b/monte/test_main.c
@@ -15,7 +15,10 @@ for (I = problems ; I->f != 0; I++)
for (i = 0; i < TRIALS ; i++)
{
MONTE_STATE *s = MONTE_ALLOC (I->dim);
-
+#ifdef MONTE_PARAMS
+ MONTE_PARAMS params;
+#endif
+
I->f->dim = I->dim;
MONTE_INTEGRATE (I->f, I->xl, I->xu,
diff --git a/monte/vegas.c b/monte/vegas.c
index e8c4d38..87d91fa 100644
--- a/monte/vegas.c
+++ b/monte/vegas.c
@@ -520,6 +520,32 @@ gsl_monte_vegas_free (gsl_monte_vegas_state * s)
free (s);
}
+double
+gsl_monte_vegas_chisq (const gsl_monte_vegas_state* s)
+{
+ return s->chisq;
+}
+
+void gsl_monte_vegas_get_params(const gsl_monte_vegas_state* s,
gsl_monte_vegas_params* p)
+{
+ p->alpha = s->alpha;
+ p->iterations = s->iterations;
+ p->stage = s->stage;
+ p->mode = s->mode;
+ p->verbose = s->verbose;
+ p->ostream = s->ostream;
+}
+
+void gsl_monte_vegas_set_params(gsl_monte_vegas_state* s, const
gsl_monte_vegas_params* p)
+{
+ s->alpha = p->alpha;
+ s->iterations = p->iterations;
+ s->stage = p->stage;
+ s->mode = p->mode;
+ s->verbose = p->verbose;
+ s->ostream = p->ostream;
+}
+
static void
init_box_coord (gsl_monte_vegas_state * s, coord box[])
{
