ext/gsl_native/interp.c in gsl-1.16.0.6 vs ext/gsl_native/interp.c in gsl-2.1.0

@@ -220,9 +220,23 @@
       for (i = 0; (int) i < na->total; i++)
         ptr[i] = (*eval)(rgi->p, ptrx, ptry, ptrz[i], rgi->a);
       return ary;
     }
 #endif
+#ifdef HAVE_NMATRIX_H
+    if (NM_IsNMatrix(xx)) {
+      NM_DENSE_STORAGE *nm;
+      double *ptrz = NULL, *ptr = NULL;
+      size_t n = NM_DENSE_COUNT(xx);
+      nm = NM_STORAGE_DENSE(xx);
+      ptrz = (double*) nm->elements;
+      ary = rb_nmatrix_dense_create(FLOAT64, nm->shape, nm->dim, nm->elements, n);
+      ptr = (double*)NM_DENSE_ELEMENTS(ary);
+      for (i = 0; (int) i < n; i++)
+        ptr[i] = (*eval)(rgi->p, ptrx, ptry, ptrz[i], rgi->a);
+      return ary;
+    }
+#endif
     if (VECTOR_P(xx)) {
       Data_Get_Struct(xx, gsl_vector, v);
       vnew = gsl_vector_alloc(v->size);
       for (i = 0; i < v->size; i++) {
         val = (*eval)(rgi->p, ptrx, ptry, gsl_vector_get(v, i), rgi->a);