aboutsummaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
Diffstat (limited to 'src')
-rw-r--r--src/config.h.in6
-rw-r--r--src/configure.ac5
-rw-r--r--src/mtx_bessel.c25
3 files changed, 29 insertions, 7 deletions
diff --git a/src/config.h.in b/src/config.h.in
index ea56ead..c86dafa 100644
--- a/src/config.h.in
+++ b/src/config.h.in
@@ -9,9 +9,15 @@
/* can gsl compute complex eigenvalues? */
#undef HAVE_GSL_EIGEN_NONSYMM
+/* can gsl compute Bessel functions? */
+#undef HAVE_GSL_BESSEL
+
/* Define to 1 if you have the <gsl/gsl_linalg.h> header file. */
#undef HAVE_GSL_GSL_LINALG_H
+/* Define to 1 if you have the <gsl/gsl_sf_bessel.h> header file. */
+#undef HAVE_GSL_GSL_LINALG_H
+
/* Define to 1 if you have the <inttypes.h> header file. */
#undef HAVE_INTTYPES_H
diff --git a/src/configure.ac b/src/configure.ac
index b6285f9..2a0a7be 100644
--- a/src/configure.ac
+++ b/src/configure.ac
@@ -69,9 +69,10 @@ dnl AC_CHECK_LIB([m], [cos])
AC_CHECK_LIB([gslcblas], [cblas_dgemm])
AC_CHECK_LIB([gsl], [gsl_blas_dgemm])
AC_CHECK_LIB([gsl], [gsl_eigen_nonsymm],AC_DEFINE(HAVE_GSL_EIGEN_NONSYMM,1,[can gsl compute complex eigenvalues?]))
+AC_CHECK_LIB([gsl], [gsl_sf_bessel_Jn],AC_DEFINE(HAVE_GSL_BESSEL,1,[can gsl compute Bessel functions?]))
dnl for math.h Bessel/Neumann functions
-AC_CHECK_LIB([m], [jnf],AC_DEFINE(HAVE_MATH_BESSEL,1,[can math compute Bessel functions?]))
+AC_CHECK_LIB([m], [jn],AC_DEFINE(HAVE_MATH_BESSEL,1,[can math compute Bessel functions?]))
if test "x$with_pd" != "x"; then
@@ -100,7 +101,7 @@ AC_CHECK_LIB(pd, nullfn)
dnl Checks for header files.
AC_HEADER_STDC
-AC_CHECK_HEADERS(stdlib.h stdio.h string.h math.h time.h sys/time.h fftw3.h sndfile.h gsl/gsl_linalg.h)
+AC_CHECK_HEADERS(stdlib.h stdio.h string.h math.h time.h sys/time.h fftw3.h sndfile.h gsl/gsl_linalg.h gsl/gsl_sf_bessel.h)
dnl Checks for typedefs, structures, and compiler characteristics.
AC_HEADER_TIME
diff --git a/src/mtx_bessel.c b/src/mtx_bessel.c
index f17342f..458100c 100644
--- a/src/mtx_bessel.c
+++ b/src/mtx_bessel.c
@@ -16,6 +16,9 @@
#include "iemmatrix.h"
#include "math.h"
#include <stdlib.h>
+#ifdef HAVE_GSL_BESSEL
+#include <gsl/gsl_sf_bessel.h>
+#endif
static t_class *mtx_bessel_class;
@@ -121,7 +124,7 @@ static void mTXBesselMatrix (MTXBessel *x, t_symbol *s,
int n,m,ofs;
-#ifdef HAVE_MATH_BESSEL
+#if defined HAVE_MATH_BESSEL || defined HAVE_GSL_BESSEL
/* size check */
if (!size)
@@ -139,16 +142,28 @@ static void mTXBesselMatrix (MTXBessel *x, t_symbol *s,
for (n=0;n<x->l;n++) {
x->kr[n]=(double) atom_getfloat(argv+n);
}
-
+
+#ifdef HAVE_GSL_BESSEL
if (x->h_re!=0)
for (m=0;m<x->l;m++)
for (n=0;n<x->nmax+1;n++)
- x->h_re[n+m*(x->nmax+1)]=jnf(n,x->kr[m]);
+ x->h_re[n+m*(x->nmax+1)]=gsl_sf_bessel_Jn(n,x->kr[m]);
if (x->h_im!=0)
for (m=0;m<x->l;m++)
for (n=0;n<x->nmax+1;n++)
- x->h_im[n+m*(x->nmax+1)]=ynf(n,x->kr[m]);
+ x->h_im[n+m*(x->nmax+1)]=gsl_sf_bessel_Yn(n,x->kr[m]);
+#else
+ if (x->h_re!=0)
+ for (m=0;m<x->l;m++)
+ for (n=0;n<x->nmax+1;n++)
+ x->h_re[n+m*(x->nmax+1)]=jn(n,x->kr[m]);
+
+ if (x->h_im!=0)
+ for (m=0;m<x->l;m++)
+ for (n=0;n<x->nmax+1;n++)
+ x->h_im[n+m*(x->nmax+1)]=yn(n,x->kr[m]);
+#endif
if (x->h_re!=0) {
@@ -169,7 +184,7 @@ static void mTXBesselMatrix (MTXBessel *x, t_symbol *s,
}
#else
- post("mtx_bessel: implementation requires math.h that implements jnf and ynf Bessel functions");
+ post("mtx_bessel: implementation requires math.h that implements jn and yn Bessel functions or gsl_sf_bessel.h");
#endif
}