Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Adapt to GSL 2.0 #31

Open
wants to merge 5 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 8 additions & 0 deletions ext/gsl_native/extconf.rb
Original file line number Diff line number Diff line change
Expand Up @@ -115,4 +115,12 @@ def gsl_gem_config(target, dir = 'ext')
gsl_have_header('tamuanova', 'tamu_anova/tamu_anova.h')
end

have_struct_member('gsl_multifit_fdfsolver', 'J', 'gsl/gsl_multifit_nlin.h')
have_func('gsl_sf_mathieu_a_e', 'gsl/gsl_sf_mathieu.h');
have_func('gsl_sf_mathieu_b_e', 'gsl/gsl_sf_mathieu.h');
have_func('gsl_sf_mathieu_ce_e', 'gsl/gsl_sf_mathieu.h');
have_func('gsl_sf_mathieu_se_e', 'gsl/gsl_sf_mathieu.h');
have_func('gsl_sf_mathieu_Mc_e', 'gsl/gsl_sf_mathieu.h');
have_func('gsl_sf_mathieu_Ms_e', 'gsl/gsl_sf_mathieu.h');

create_makefile('gsl_native')
42 changes: 42 additions & 0 deletions ext/gsl_native/histogram.c
Original file line number Diff line number Diff line change
Expand Up @@ -1169,6 +1169,9 @@ static VALUE rb_gsl_histogram_fit_gaussian(int argc, VALUE *argv, VALUE obj)
size_t n, dof; /* # of data points */
size_t p = 3; /* # of fitting parameters */
gsl_multifit_function_fdf f;
#ifndef GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
gsl_matrix *covar = NULL;
gsl_vector *x = NULL;
double sigma, mean, height, errs, errm, errh, chi2;
Expand Down Expand Up @@ -1197,6 +1200,9 @@ static VALUE rb_gsl_histogram_fit_gaussian(int argc, VALUE *argv, VALUE obj)
hh.binend = binend;
n = binend - binstart + 1;

#ifndef GSL_MULTIFIT_FDFSOLVER_J
J = gsl_matrix_alloc(n, p);
#endif
covar = gsl_matrix_alloc(p, p);

f.f = Gaussian_f;
Expand All @@ -1219,7 +1225,12 @@ static VALUE rb_gsl_histogram_fit_gaussian(int argc, VALUE *argv, VALUE obj)
sigma = sqrt(gsl_vector_get(s->x, 0));
mean = gsl_vector_get(s->x, 1);
height = gsl_vector_get(s->x, 2)*sigma*sqrt(2*M_PI);
#ifdef GSL_MULTIFIT_FDFSOLVER_J
gsl_multifit_covar(s->J, 0.0, covar);
#else
gsl_multifit_fdfsolver_jac(s, J);
gsl_multifit_covar(J, 0.0, covar);
#endif
chi2 = gsl_pow_2(gsl_blas_dnrm2(s->f)); /* not reduced chi-square */
dof = n - p;
errs = sqrt(chi2/dof*gsl_matrix_get(covar, 0, 0))/sigma/2;
Expand All @@ -1228,6 +1239,9 @@ static VALUE rb_gsl_histogram_fit_gaussian(int argc, VALUE *argv, VALUE obj)

gsl_multifit_fdfsolver_free(s);
gsl_vector_free(x);
#ifndef GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix_free(J);
#endif
gsl_matrix_free(covar);
return rb_ary_new3(8, rb_float_new(sigma), rb_float_new(mean),
rb_float_new(height), rb_float_new(errs),
Expand Down Expand Up @@ -1305,6 +1319,9 @@ static VALUE rb_gsl_histogram_fit_rayleigh(int argc, VALUE *argv, VALUE obj)
size_t n, dof; /* # of data points */
size_t p = 2; /* # of fitting parameters */
gsl_multifit_function_fdf f;
#ifndef GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
gsl_matrix *covar = NULL;
gsl_vector *x = NULL;
double sigma, height, errs, errh, chi2;
Expand Down Expand Up @@ -1332,6 +1349,9 @@ static VALUE rb_gsl_histogram_fit_rayleigh(int argc, VALUE *argv, VALUE obj)
hh.binend = binend;
n = binend - binstart + 1;

#ifndef GSL_MULTIFIT_FDFSOLVER_J
J = gsl_matrix_alloc(n, p);
#endif
covar = gsl_matrix_alloc(p, p);

f.f = Rayleigh_f;
Expand All @@ -1353,14 +1373,22 @@ static VALUE rb_gsl_histogram_fit_rayleigh(int argc, VALUE *argv, VALUE obj)
} while (status == GSL_CONTINUE);
sigma = sqrt(gsl_vector_get(s->x, 0));
height = gsl_vector_get(s->x, 1)*sigma*sigma;
#ifdef GSL_MULTIFIT_FDFSOLVER_J
gsl_multifit_covar(s->J, 0.0, covar);
#else
gsl_multifit_fdfsolver_jac(s, J);
gsl_multifit_covar(J, 0.0, covar);
#endif
chi2 = gsl_pow_2(gsl_blas_dnrm2(s->f)); /* not reduced chi-square */
dof = n - p;
errs = sqrt(chi2/dof*gsl_matrix_get(covar, 0, 0))/sigma/2;
errh = sqrt(chi2/dof*gsl_matrix_get(covar, 1, 1));

gsl_multifit_fdfsolver_free(s);
gsl_vector_free(x);
#ifndef GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix_free(J);
#endif
gsl_matrix_free(covar);
return rb_ary_new3(6, rb_float_new(sigma),
rb_float_new(height), rb_float_new(errs),
Expand Down Expand Up @@ -1441,6 +1469,9 @@ static VALUE rb_gsl_histogram_fit_xexponential(int argc, VALUE *argv, VALUE obj)
size_t n, dof; /* # of data points */
size_t p = 2; /* # of fitting parameters */
gsl_multifit_function_fdf f;
#ifndef GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
gsl_matrix *covar = NULL;
gsl_vector *x = NULL;
double b, height, errs, errh, chi2;
Expand Down Expand Up @@ -1468,6 +1499,9 @@ static VALUE rb_gsl_histogram_fit_xexponential(int argc, VALUE *argv, VALUE obj)
hh.binend = binend;
n = binend - binstart + 1;

#ifndef GSL_MULTIFIT_FDFSOLVER_J
J = gsl_matrix_alloc(n, p);
#endif
covar = gsl_matrix_alloc(p, p);

f.f = xExponential_f;
Expand All @@ -1489,14 +1523,22 @@ static VALUE rb_gsl_histogram_fit_xexponential(int argc, VALUE *argv, VALUE obj)
} while (status == GSL_CONTINUE);
b = gsl_vector_get(s->x, 0);
height = gsl_vector_get(s->x, 1);
#ifdef GSL_MULTIFIT_FDFSOLVER_J
gsl_multifit_covar(s->J, 0.0, covar);
#else
gsl_multifit_fdfsolver_jac(s, J);
gsl_multifit_covar(J, 0.0, covar);
#endif
chi2 = gsl_pow_2(gsl_blas_dnrm2(s->f)); /* not reduced chi-square */
dof = n - p;
errs = sqrt(chi2/dof*gsl_matrix_get(covar, 0, 0));
errh = sqrt(chi2/dof*gsl_matrix_get(covar, 1, 1));

gsl_multifit_fdfsolver_free(s);
gsl_vector_free(x);
#ifndef GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix_free(J);
#endif
gsl_matrix_free(covar);
return rb_ary_new3(6, rb_float_new(b),
rb_float_new(height), rb_float_new(errs),
Expand Down
71 changes: 71 additions & 0 deletions ext/gsl_native/multifit.c
Original file line number Diff line number Diff line change
Expand Up @@ -324,14 +324,24 @@ static VALUE rb_gsl_multifit_fdfsolver_test_gradient(int argc, VALUE *argv, VALU
{
gsl_multifit_fdfsolver *solver = NULL;
gsl_vector *g = NULL;
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
int status;
double epsabs;
Data_Get_Struct(obj, gsl_multifit_fdfsolver, solver);
switch (argc) {
case 1:
Need_Float(argv[0]);
g = gsl_vector_alloc(solver->x->size);
#ifdef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_multifit_gradient(solver->J, solver->f, g);
#else
J = gsl_matrix_alloc(solver->fdf->n, solver->fdf->p);
gsl_multifit_fdfsolver_jac(solver, J);
gsl_multifit_gradient(J, solver->f, g);
gsl_matrix_free(J);
#endif
epsabs = NUM2DBL(argv[0]);
status = gsl_multifit_test_gradient(g, epsabs);
gsl_vector_free(g);
Expand All @@ -352,23 +362,47 @@ static VALUE rb_gsl_multifit_fdfsolver_test_gradient(int argc, VALUE *argv, VALU
static VALUE rb_gsl_multifit_fdfsolver_gradient(int argc, VALUE *argv, VALUE obj)
{
gsl_multifit_fdfsolver *solver = NULL;
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
gsl_vector *g = NULL;
// local variable "status" declared and set, but never used
//int status;
Data_Get_Struct(obj, gsl_multifit_fdfsolver, solver);
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix_alloc(solver->fdf->n, solver->fdf->p);
gsl_multifit_fdfsolver_jac(solver, J);
#endif
if (argc == 1) {
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
int retval = 0;
#endif
Data_Get_Vector(argv[0], g);
#ifdef HAVE_GSL_MULTIFIT_FDFSOLVER_J
return INT2FIX(gsl_multifit_gradient(solver->J, solver->f, g));
#else
retval = gsl_multifit_gradient(J, solver->f, g);
gsl_matrix_free(J);
return INT2FIX(retval);
#endif
} else {
g = gsl_vector_alloc(solver->x->size);
#ifdef HAVE_GSL_MULTIFIT_FDFSOLVER_J
/*status =*/ gsl_multifit_gradient(solver->J, solver->f, g);
#else
/*status =*/ gsl_multifit_gradient(J, solver->f, g);
gsl_matrix_free(J);
#endif
return Data_Wrap_Struct(cgsl_vector, 0, gsl_vector_free, g);
}
}

static VALUE rb_gsl_multifit_fdfsolver_covar(int argc, VALUE *argv, VALUE obj)
{
gsl_multifit_fdfsolver *solver = NULL;
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
gsl_matrix *covar = NULL;
double epsrel;
// local variable "status" declared and set, but never used
Expand All @@ -378,14 +412,32 @@ static VALUE rb_gsl_multifit_fdfsolver_covar(int argc, VALUE *argv, VALUE obj)
Data_Get_Struct(obj, gsl_multifit_fdfsolver, solver);
epsrel = NUM2DBL(argv[0]);
switch (argc) {
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
int retval = 0;
#endif
case 1:
covar = gsl_matrix_alloc(solver->x->size, solver->x->size);
#ifdef HAVE_GSL_MULTIFIT_FDFSOLVER_J
/*status =*/ gsl_multifit_covar(solver->J, epsrel, covar);
#else
J = gsl_matrix_alloc(solver->fdf->n, solver->fdf->p);
gsl_multifit_fdfsolver_jac(solver, J);
/*status =*/ gsl_multifit_covar(J, epsrel, covar);
gsl_matrix_free(J);
#endif
return Data_Wrap_Struct(cgsl_matrix, 0, gsl_matrix_free, covar);
break;
case 2:
Data_Get_Matrix(argv[1], covar);
#ifdef HAVE_GSL_MULTIFIT_FDFSOLVER_J
return INT2FIX(gsl_multifit_covar(solver->J, epsrel, covar));
#else
J = gsl_matrix_alloc(solver->fdf->n, solver->fdf->p);
gsl_multifit_fdfsolver_jac(solver, J);
retval = gsl_multifit_covar(J, epsrel, covar);
gsl_matrix_free(J);
return INT2FIX(retval);
#endif
break;
default:
rb_raise(rb_eArgError, "wrong number of arguments");
Expand Down Expand Up @@ -417,8 +469,17 @@ static VALUE rb_gsl_multifit_fdfsolver_f(VALUE obj)
static VALUE rb_gsl_multifit_fdfsolver_J(VALUE obj)
{
gsl_multifit_fdfsolver *solver = NULL;
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
Data_Get_Struct(obj, gsl_multifit_fdfsolver, solver);
#ifdef HAVE_GSL_MULTIFIT_FDFSOLVER_J
return Data_Wrap_Struct(cgsl_matrix_view_ro, 0, NULL, solver->J);
#else
J = gsl_matrix_alloc(solver->fdf->n, solver->fdf->p);
gsl_multifit_fdfsolver_jac(solver, J);
return Data_Wrap_Struct(cgsl_matrix_view_ro, 0, gsl_matrix_free, J);
#endif
}

/* singleton */
Expand Down Expand Up @@ -1626,6 +1687,9 @@ static VALUE rb_gsl_multifit_fit(int argc, VALUE *argv, VALUE klass)
size_t n, dof; /* # of data points */
size_t p; /* # of fitting parameters */
gsl_multifit_function_fdf f;
#ifndef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_matrix *J = NULL;
#endif
gsl_matrix *covar = NULL;
gsl_vector *v = NULL;
gsl_vector *x, *y, *w = NULL;
Expand Down Expand Up @@ -1699,7 +1763,14 @@ static VALUE rb_gsl_multifit_fit(int argc, VALUE *argv, VALUE klass)
covar = gsl_matrix_alloc(p, p);
chi2 = gsl_pow_2(gsl_blas_dnrm2(solver->f)); /* not reduced chi-square */
dof = n - p;
#ifdef HAVE_GSL_MULTIFIT_FDFSOLVER_J
gsl_multifit_covar(solver->J, 0.0, covar);
#else
J = gsl_matrix_alloc(n, p);
gsl_multifit_fdfsolver_jac(solver, J);
gsl_multifit_covar(J, 0.0, covar);
gsl_matrix_free(J);
#endif
for (i = 0; i < p; i++)
gsl_vector_set(verr, i, sqrt(chi2/dof*gsl_matrix_get(covar, i, i)));
gsl_matrix_free(covar);
Expand Down
12 changes: 10 additions & 2 deletions ext/gsl_native/sf_ellint.c
Original file line number Diff line number Diff line change
Expand Up @@ -82,17 +82,25 @@ static VALUE rb_gsl_sf_ellint_P_e(VALUE obj, VALUE phi, VALUE k,
static VALUE rb_gsl_sf_ellint_D(int argc, VALUE *argv, VALUE obj)
{
if (argc == 3)
#if GSL_MAJOR_VERSION > 1
return rb_gsl_sf_eval_double2_m(gsl_sf_ellint_D, argv[0], argv[1],
#else
return rb_gsl_sf_eval_double3_m(gsl_sf_ellint_D, argv[0], argv[1], argv[2],
#endif
INT2FIX(GSL_PREC_DOUBLE));
else
#if GSL_MAJOR_VERSION > 1
return rb_gsl_sf_eval_double2_m(gsl_sf_ellint_D, argv[0], argv[1],
#else
return rb_gsl_sf_eval_double3_m(gsl_sf_ellint_D, argv[0], argv[1], argv[2],
#endif
argv[3]);
}

static VALUE rb_gsl_sf_ellint_D_e(VALUE obj, VALUE phi, VALUE k,
VALUE n, VALUE m)
VALUE m)
{
return rb_gsl_sf_eval_e_double3_m(gsl_sf_ellint_D_e, phi, k, n, m);
return rb_gsl_sf_eval_e_double2_m(gsl_sf_ellint_D_e, phi, k, m);
}

static VALUE rb_gsl_sf_ellint_RC(int argc, VALUE *argv, VALUE obj)
Expand Down
Loading