Markus Mottl avatar Markus Mottl committed 5c3e770

Copied latest C-file changes from regular ocamlgsl SVN tree

Comments (0)

Files changed (5)

lib/mlgsl_error.c

   return copy_string(gsl_strerror(gsl_errno));
 }
 
-static inline int conv_err_code(int gsl_errno)
-{
-  if(gsl_errno < 0)
-    return gsl_errno + 2 ;
-  else
-    return gsl_errno + 1 ;
-}
-
-static value       *ml_gsl_exn;
-
-static void ml_gsl_raise_exn(const char *msg, int gsl_errno)
-{
-  CAMLlocal2(exn_msg, exn_arg);
-  exn_msg = copy_string(msg);
-  exn_arg = alloc_small(2, 0);
-  Field(exn_arg, 0) = Val_int(conv_err_code(gsl_errno));
-  Field(exn_arg, 1) = exn_msg;
-  if(ml_gsl_exn != NULL)
-    raise_with_arg(*ml_gsl_exn, exn_arg);
-  else
-    failwith("GSL error");
-}
+static value *ml_gsl_err_handler;
 
 static void ml_gsl_error_handler(const char *reason, const char *file,
 				 int line, int gsl_errno)
 {
-  ml_gsl_raise_exn(reason, gsl_errno);
+  CAMLparam0();
+  CAMLlocal2(exn_msg, exn_arg);
+  int ml_errno;
+
+  if (0 < gsl_errno && gsl_errno <= GSL_EOF)
+    ml_errno = gsl_errno + 1;
+  else if (GSL_CONTINUE <= gsl_errno && gsl_errno <= GSL_FAILURE)
+    ml_errno = gsl_errno + 2;
+  else
+    failwith("invalid GSL error code");
+
+  exn_msg = copy_string(reason);
+  exn_arg = alloc_small(2, 0);
+  Field(exn_arg, 0) = Val_int(ml_errno);
+  Field(exn_arg, 1) = exn_msg;
+  callback(Field(*ml_gsl_err_handler, 0), exn_arg);
+
+  CAMLreturn0;
 }
 
 CAMLprim value ml_gsl_error_init(value init)
 {
   static gsl_error_handler_t *old;
-  if(ml_gsl_exn == NULL) 
-    ml_gsl_exn = caml_named_value("mlgsl_exn");
-  if(Bool_val(init))
-    old = gsl_set_error_handler(&ml_gsl_error_handler);
+  if(ml_gsl_err_handler == NULL)
+    ml_gsl_err_handler = caml_named_value("mlgsl_err_handler");
+
+  if (Bool_val(init)) {
+    gsl_error_handler_t *prev;
+    prev = gsl_set_error_handler(&ml_gsl_error_handler);
+    if (prev != ml_gsl_error_handler)
+      old = prev;
+  }
   else
     gsl_set_error_handler(old);
+
   return Val_unit;
 }
 
 CAMLprim value ml_gsl_fft_complex_rad2_forward(value dif, value stride, value data)
 {
-  size_t N = Double_array_length(data);
+  size_t N = Double_array_length(data) / 2;
   size_t c_stride = Opt_arg(stride, Int_val, 1);
   int c_dif = Opt_arg(dif, Bool_val, 0);
 
 CAMLprim value ml_gsl_fft_complex_rad2_transform(value dif, value stride, 
 					value data, value sign)
 {
-  size_t N = Double_array_length(data);
+  size_t N = Double_array_length(data) / 2;
   size_t c_stride = Opt_arg(stride, Int_val, 1);
   int c_dif = Opt_arg(dif, Bool_val, 0);
   gsl_fft_direction c_sign =
 
 CAMLprim value ml_gsl_fft_complex_rad2_backward(value dif, value stride, value data)
 {
-  size_t N = Double_array_length(data);
+  size_t N = Double_array_length(data) / 2;
   size_t c_stride = Opt_arg(stride, Int_val, 1);
   int c_dif = Opt_arg(dif, Bool_val, 0);
 
 
 CAMLprim value ml_gsl_fft_complex_rad2_inverse(value dif, value stride, value data)
 {
-  size_t N = Double_array_length(data);
+  size_t N = Double_array_length(data) / 2;
   size_t c_stride = Opt_arg(stride, Int_val, 1);
   int c_dif = Opt_arg(dif, Bool_val, 0);
 
   struct callback_params *p=params;
   value res;
   value v_x = copy_double(x);
-  res=callback_exn(p->closure, v_x);
-  if(Is_exception_result(res))
-    return GSL_NAN;
+  res=callback(p->closure, v_x);
   return Double_val(res);
 }
 
+/* FDF CALLBACKS */
 double gslfun_callback_indir(double x, void *params)
 {
   value res;
   value v_x = copy_double(x);
   value *closure = params;
-  res=callback_exn(*closure, v_x);
-  if(Is_exception_result(res))
-    return GSL_NAN;
+  res=callback(*closure, v_x);
   return Double_val(res);
 }
  
-/* FDF CALLBACKS */
 double gslfun_callback_f(double x, void *params)
 {
   struct callback_params *p=params;
   value res;
   value v_x=copy_double(x);
-  res=callback_exn(Field(p->closure, 0), v_x);
-  if(Is_exception_result(res))
-    return GSL_NAN;
+  res=callback(Field(p->closure, 0), v_x);
   return Double_val(res);
 }
 
   struct callback_params *p=params;
   value res;
   value v_x=copy_double(x);
-  res=callback_exn(Field(p->closure, 1), v_x);
-  if(Is_exception_result(res))
-    return GSL_NAN;
+  res=callback(Field(p->closure, 1), v_x);
   return Double_val(res);
 }
 
   struct callback_params *p=params;
   value res;
   value v_x=copy_double(x);
-  res=callback_exn(Field(p->closure, 2), v_x);
-  if(Is_exception_result(res)){
-    *f=GSL_NAN; *df=GSL_NAN;
-    return;
-  }
+  res=callback(Field(p->closure, 2), v_x);
   *f =Double_val(Field(res, 0));
   *df=Double_val(Field(res, 1));
 }
   value res;
 
   memcpy(Double_array_val(p->dbl), x_arr, dim*sizeof(double));
-  res=callback_exn(p->closure, p->dbl);
-  if(Is_exception_result(res))
-    return GSL_NAN;
+  res=callback(p->closure, p->dbl);
   return Double_val(res);
 }
 
   struct callback_params *p=params;
   value res;
 
-  res=callback_exn(p->closure, (value)x_arr);
-  if(Is_exception_result(res))
-    return GSL_NAN;
+  res=callback(p->closure, (value)x_arr);
   return Double_val(res);
 }
 
 /* MULTIROOT CALLBACKS */
 int gsl_multiroot_callback(const gsl_vector *x, void *params, gsl_vector *F)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr, f_barr;
   int len = x->size;
-  LOCALARRAY(double, x_arr, len); 
-  LOCALARRAY(double, f_arr, len); 
   gsl_vector_view x_v, f_v;
-  value res;
 
-  x_v = gsl_vector_view_array(x_arr, len);
-  f_v = gsl_vector_view_array(f_arr, len);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  f_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+  f_v = gsl_vector_view_array(Data_bigarray_val(f_barr), len);
+
   gsl_vector_memcpy(&x_v.vector, x);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, len);
-  f_barr = alloc_bigarray_dims(barr_flags, 1, f_arr, len);
-  res=callback2_exn(p->closure, x_barr, f_barr);
-  if(Is_exception_result(res))
-    return GSL_FAILURE;
+  callback2(p->closure, x_barr, f_barr);
   gsl_vector_memcpy(F, &f_v.vector);
   return GSL_SUCCESS;
 }
 
 int gsl_multiroot_callback_f(const gsl_vector *x, void *params, gsl_vector *F)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr, f_barr;
   int len = x->size;
-  LOCALARRAY(double, x_arr, len); 
-  LOCALARRAY(double, f_arr, len); 
   gsl_vector_view x_v, f_v;
-  value res;
 
-  x_v = gsl_vector_view_array(x_arr, len);
-  f_v = gsl_vector_view_array(f_arr, len);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  f_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+  f_v = gsl_vector_view_array(Data_bigarray_val(f_barr), len);
+
   gsl_vector_memcpy(&x_v.vector, x);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, len);
-  f_barr = alloc_bigarray_dims(barr_flags, 1, f_arr, len);
-  res=callback2_exn(Field(p->closure, 0), x_barr, f_barr);
-  if(Is_exception_result(res))
-    return GSL_FAILURE;
+  callback2(Field(p->closure, 0), x_barr, f_barr);
   gsl_vector_memcpy(F, &f_v.vector);
   return GSL_SUCCESS;
 }
 
 int gsl_multiroot_callback_df(const gsl_vector *x, void *params, gsl_matrix *J)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr, j_barr;
   int len = x->size;
-  LOCALARRAY(double, x_arr, len); 
-  LOCALARRAY(double, j_arr, len*len); 
   gsl_vector_view x_v;
   gsl_matrix_view j_v;
-  value res;
 
-  x_v = gsl_vector_view_array(x_arr, len);
-  j_v = gsl_matrix_view_array(j_arr, len, len);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  j_barr = alloc_bigarray_dims(barr_flags, 2, NULL, len, len);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+  j_v = gsl_matrix_view_array(Data_bigarray_val(j_barr), len, len);
+
   gsl_vector_memcpy(&x_v.vector, x);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, len);
-  j_barr = alloc_bigarray_dims(barr_flags, 2, j_arr, len, len);
-  res=callback2_exn(Field(p->closure, 1), x_barr, j_barr);
-  if(Is_exception_result(res))
-    return GSL_FAILURE;
+  callback2(Field(p->closure, 1), x_barr, j_barr);
   gsl_matrix_memcpy(J, &j_v.matrix);
   return GSL_SUCCESS;
 }
 int gsl_multiroot_callback_fdf(const gsl_vector *x, void *params, 
 			   gsl_vector *F, gsl_matrix *J)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr, f_barr, j_barr;
   int len = x->size;
-  LOCALARRAY(double, x_arr, len); 
-  LOCALARRAY(double, f_arr, len); 
-  LOCALARRAY(double, j_arr, len*len); 
   gsl_vector_view x_v, f_v;
   gsl_matrix_view j_v;
-  value res;
   
-  x_v = gsl_vector_view_array(x_arr, len);
-  f_v = gsl_vector_view_array(f_arr, len);
-  j_v = gsl_matrix_view_array(j_arr, len, len);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  f_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  j_barr = alloc_bigarray_dims(barr_flags, 2, NULL, len, len);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+  f_v = gsl_vector_view_array(Data_bigarray_val(f_barr), len);
+  j_v = gsl_matrix_view_array(Data_bigarray_val(j_barr), len, len);
+
   gsl_vector_memcpy(&x_v.vector, x);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, len);
-  f_barr = alloc_bigarray_dims(barr_flags, 1, f_arr, len);
-  j_barr = alloc_bigarray_dims(barr_flags, 2, j_arr, len, len);
-  res=callback3_exn(Field(p->closure, 2), x_barr, f_barr, j_barr);
-  if(Is_exception_result(res))
-    return GSL_FAILURE;
+  callback3(Field(p->closure, 2), x_barr, f_barr, j_barr);
   gsl_vector_memcpy(F, &f_v.vector);
   gsl_matrix_memcpy(J, &j_v.matrix);
   return GSL_SUCCESS;
   int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr;
-  intnat len = x->size;
+  int len = x->size;
   gsl_vector_view x_v;
   value res;
 
   x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
   x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+
   gsl_vector_memcpy(&x_v.vector, x);
-  res=callback_exn(p->closure, x_barr);
-
-  if(Is_exception_result(res)) {
-    return GSL_NAN;
-  }
-
+  res=callback(p->closure, x_barr);
   return Double_val(res);
 }
 
 double gsl_multimin_callback_f(const gsl_vector *x, void *params)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr;
   int len = x->size;
-
-/* CR mmottl: Stack allocation may segfault with large lengths, hence
-   malloc.  Note that the OCaml callbacks may put the bigarrays into
-   references.  This is evil, and these bindings should really get fixed
-   to avoid problems of that sort. */
-  double *x_arr = malloc(sizeof(double) * len);
-#if 0
-  LOCALARRAY(double, x_arr, len);
-#endif
-
   gsl_vector_view x_v;
   value res;
 
-  x_v = gsl_vector_view_array(x_arr, len);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+
   gsl_vector_memcpy(&x_v.vector, x);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, len);
-  res=callback_exn(Field(p->closure, 0), x_barr);
-
-  /* CR mmottl: need to free malloced memory now */
-  free(x_arr);
-
-  if(Is_exception_result(res))
-    return GSL_NAN;
-
+  res=callback(Field(p->closure, 0), x_barr);
   return Double_val(res);
 }
 
 void gsl_multimin_callback_df(const gsl_vector *x, void *params, gsl_vector *G)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr, g_barr;
   int len = x->size;
+  gsl_vector_view x_v, g_v;
 
-/* CR mmottl: Stack allocation may segfault with large lengths, hence
-   malloc.  Note that the OCaml callbacks may put the bigarrays into
-   references.  This is evil, and these bindings should really get fixed
-   to avoid problems of that sort. */
-  double *x_arr = malloc(sizeof(double) * len);
-  double *g_arr = malloc(sizeof(double) * len);
-#if 0
-  LOCALARRAY(double, x_arr, len);
-  LOCALARRAY(double, g_arr, len); 
-#endif
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  g_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+  g_v = gsl_vector_view_array(Data_bigarray_val(g_barr), len);
 
-  gsl_vector_view x_v, g_v;
-  value res;
-
-  x_v = gsl_vector_view_array(x_arr, len);
-  g_v = gsl_vector_view_array(g_arr, len);
   gsl_vector_memcpy(&x_v.vector, x);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, len);
-  g_barr = alloc_bigarray_dims(barr_flags, 1, g_arr, len);
-  res=callback2_exn(Field(p->closure, 1), x_barr, g_barr);
-
-  if(Is_exception_result(res)){
-    /* CR mmottl: need to free malloced memory now */
-    free(x_arr);
-    free(g_arr);
-
-    /* the caml functions raised an exception but there's no way we can
-       indicate this to GSL since the return type is void.
-       So we set the out param G to NaN. */
-    gsl_vector_set_all(G, GSL_NAN);
-    return;
-  }
-
+  callback2(Field(p->closure, 1), x_barr, g_barr);
   gsl_vector_memcpy(G, &g_v.vector);
-
-  /* CR mmottl: need to free malloced memory now */
-  free(x_arr);
-  free(g_arr);
-
 }
 
 void gsl_multimin_callback_fdf(const gsl_vector *x, void *params, 
 			       double *f, gsl_vector *G)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *p=params;
   value x_barr, g_barr;
   int len = x->size;
-
-/* CR mmottl: Stack allocation may segfault with large lengths, hence
-   malloc.  Note that the OCaml callbacks may put the bigarrays into
-   references.  This is evil, and these bindings should really get fixed
-   to avoid problems of that sort. */
-  double *x_arr = malloc(sizeof(double) * len);
-  double *g_arr = malloc(sizeof(double) * len);
-#if 0
-  LOCALARRAY(double, x_arr, len);
-  LOCALARRAY(double, g_arr, len); 
-#endif
-
   gsl_vector_view x_v, g_v;
   value res;
   
-  x_v = gsl_vector_view_array(x_arr, len);
-  g_v = gsl_vector_view_array(g_arr, len);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  g_barr = alloc_bigarray_dims(barr_flags, 1, NULL, len);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), len);
+  g_v = gsl_vector_view_array(Data_bigarray_val(g_barr), len);
+
   gsl_vector_memcpy(&x_v.vector, x);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, len);
-  g_barr = alloc_bigarray_dims(barr_flags, 1, g_arr, len);
-  res=callback2_exn(Field(p->closure, 2), x_barr, g_barr);
-
-  if(Is_exception_result(res)){
-    /* CR mmottl: need to free malloced memory now */
-    free(x_arr);
-    free(g_arr);
-
-    *f=GSL_NAN;
-    return;
-  }
+  res=callback2(Field(p->closure, 2), x_barr, g_barr);
   gsl_vector_memcpy(G, &g_v.vector);
-
-  /* CR mmottl: need to free malloced memory now */
-  free(x_arr);
-  free(g_arr);
-
   *f=Double_val(res);
 }
 
 /* MULTIFIT CALLBACKS */
 int gsl_multifit_callback_f(const gsl_vector *X, void *params, gsl_vector *F)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *parms=params;
   value x_barr, f_barr;
   size_t p = X->size;
   size_t n = F->size;
-  LOCALARRAY(double, x_arr, p); 
-  LOCALARRAY(double, f_arr, n); 
   gsl_vector_view x_v, f_v;
-  value res;
 
-  x_v = gsl_vector_view_array(x_arr, p);
-  f_v = gsl_vector_view_array(f_arr, n);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, p);
+  f_barr = alloc_bigarray_dims(barr_flags, 1, NULL, n);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), p);
+  f_v = gsl_vector_view_array(Data_bigarray_val(f_barr), n);
+
   gsl_vector_memcpy(&x_v.vector, X);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, p);
-  f_barr = alloc_bigarray_dims(barr_flags, 1, f_arr, n);
-  res=callback2_exn(Field(parms->closure, 0), x_barr, f_barr);
-  if(Is_exception_result(res))
-    return GSL_FAILURE;
+  callback2(Field(parms->closure, 0), x_barr, f_barr);
   gsl_vector_memcpy(F, &f_v.vector);
   return GSL_SUCCESS;
 }
 
 int gsl_multifit_callback_df(const gsl_vector *X, void *params, gsl_matrix *J)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *parms=params;
   value x_barr, j_barr;
   size_t p = X->size;
   size_t n = J->size1;
-  LOCALARRAY(double, x_arr, p); 
-  LOCALARRAY(double, j_arr, n*p); 
   gsl_vector_view x_v;
   gsl_matrix_view j_v;
   value res;
 
-  x_v = gsl_vector_view_array(x_arr, p);
-  j_v = gsl_matrix_view_array(j_arr, n, p);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, p);
+  j_barr = alloc_bigarray_dims(barr_flags, 2, NULL, n, p);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), p);
+  j_v = gsl_matrix_view_array(Data_bigarray_val(j_barr), n, p);
+
   gsl_vector_memcpy(&x_v.vector, X);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, p);
-  j_barr = alloc_bigarray_dims(barr_flags, 2, j_arr, n, p);
-  res=callback2_exn(Field(parms->closure, 1), x_barr, j_barr);
+  res=callback2(Field(parms->closure, 1), x_barr, j_barr);
   if(Is_exception_result(res))
     return GSL_FAILURE;
   gsl_matrix_memcpy(J, &j_v.matrix);
 int gsl_multifit_callback_fdf(const gsl_vector *X, void *params, 
 			      gsl_vector *F, gsl_matrix *J)
 {
-  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT | BIGARRAY_EXTERNAL;
+  int barr_flags = BIGARRAY_FLOAT64 | BIGARRAY_C_LAYOUT;
   struct callback_params *parms=params;
   value x_barr, f_barr, j_barr;
   size_t p = X->size;
   size_t n = F->size;
-  LOCALARRAY(double, x_arr, p); 
-  LOCALARRAY(double, f_arr, n); 
-  LOCALARRAY(double, j_arr, n*p); 
   gsl_vector_view x_v, f_v;
   gsl_matrix_view j_v;
-  value res;
   
-  x_v = gsl_vector_view_array(x_arr, p);
-  f_v = gsl_vector_view_array(f_arr, n);
-  j_v = gsl_matrix_view_array(j_arr, n, p);
+  x_barr = alloc_bigarray_dims(barr_flags, 1, NULL, p);
+  f_barr = alloc_bigarray_dims(barr_flags, 1, NULL, n);
+  j_barr = alloc_bigarray_dims(barr_flags, 2, NULL, n, p);
+  x_v = gsl_vector_view_array(Data_bigarray_val(x_barr), p);
+  f_v = gsl_vector_view_array(Data_bigarray_val(f_barr), n);
+  j_v = gsl_matrix_view_array(Data_bigarray_val(j_barr), n, p);
+
   gsl_vector_memcpy(&x_v.vector, X);
-  x_barr = alloc_bigarray_dims(barr_flags, 1, x_arr, p);
-  f_barr = alloc_bigarray_dims(barr_flags, 1, f_arr, n);
-  j_barr = alloc_bigarray_dims(barr_flags, 2, j_arr, n, p);
-  res=callback3_exn(Field(parms->closure, 2), x_barr, f_barr, j_barr);
-  if(Is_exception_result(res))
-    return GSL_FAILURE;
+  callback3(Field(parms->closure, 2), x_barr, f_barr, j_barr);
   gsl_vector_memcpy(F, &f_v.vector);
   gsl_matrix_memcpy(J, &j_v.matrix);
   return GSL_SUCCESS;

lib/mlgsl_permut.c

 			      Double_array_length(arr)/2);
   return Val_unit;
 }
+
+CAMLprim value ml_gsl_permute_mul(value p, value pa, value pb)
+{
+  GSL_PERMUT_OF_BIGARRAY(p);
+  GSL_PERMUT_OF_BIGARRAY(pa);
+  GSL_PERMUT_OF_BIGARRAY(pb);
+  gsl_permutation_mul(&perm_p, &perm_pa, &perm_pb);
+  return Val_unit;
+}
+
+CAMLprim value ml_gsl_permute_linear_to_canonical(value q, value p)
+{
+  GSL_PERMUT_OF_BIGARRAY(q);
+  GSL_PERMUT_OF_BIGARRAY(p);
+  gsl_permutation_linear_to_canonical (&perm_q, &perm_p);
+  return Val_unit;
+}
+
+CAMLprim value ml_gsl_permute_canonical_to_linear(value p, value q)
+{
+  GSL_PERMUT_OF_BIGARRAY(p);
+  GSL_PERMUT_OF_BIGARRAY(q);
+  gsl_permutation_canonical_to_linear (&perm_p, &perm_q);
+  return Val_unit;
+}
+
+CAMLprim value ml_gsl_permute_inversions(value p) 
+{
+  size_t inv;
+  GSL_PERMUT_OF_BIGARRAY(p);
+  inv = gsl_permutation_inversions (&perm_p);
+  return Val_long(inv);
+}
+
+CAMLprim value ml_gsl_permute_canonical_cycles(value q) 
+{
+  size_t c;
+  GSL_PERMUT_OF_BIGARRAY(q);
+  c = gsl_permutation_canonical_cycles (&perm_q);
+  return Val_long(c);
+}
+
+CAMLprim value ml_gsl_permute_linear_cycles(value p) 
+{
+  size_t c;
+  GSL_PERMUT_OF_BIGARRAY(p);
+  c = gsl_permutation_linear_cycles (&perm_p);
+  return Val_long(c);
+}

lib/mlgsl_stats.c

 CAMLprim value ml_gsl_stats_correlation(value data1, value data2)
 {
   size_t len = Double_array_length(data1);
+  double r;
   check_array_size(data1, data2);
-  double r = gsl_stats_correlation(Double_array_val(data1), 1,
-                                   Double_array_val(data2), 1, len);
+  r = gsl_stats_correlation(Double_array_val(data1), 1,
+                            Double_array_val(data2), 1, len);
   return copy_double(r);
 }
Tip: Filter by directory path e.g. /media app.js to search for public/media/app.js.
Tip: Use camelCasing e.g. ProjME to search for ProjectModifiedEvent.java.
Tip: Filter by extension type e.g. /repo .js to search for all .js files in the /repo directory.
Tip: Separate your search with spaces e.g. /ssh pom.xml to search for src/ssh/pom.xml.
Tip: Use ↑ and ↓ arrow keys to navigate and return to view the file.
Tip: You can also navigate files with Ctrl+j (next) and Ctrl+k (previous) and view the file with Ctrl+o.
Tip: You can also navigate files with Alt+j (next) and Alt+k (previous) and view the file with Alt+o.