interpm_krige.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2017-2018, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_INTERPM_KRIGE_H
24 #define O2SCL_INTERPM_KRIGE_H
25 
26 /** \file interpm_krige.h
27  \brief File defining \ref o2scl::interpm_krige and
28  \ref o2scl::interpm_krige_nn
29 */
30 
31 #include <iostream>
32 #include <string>
33 #include <cmath>
34 
35 #include <boost/numeric/ublas/matrix.hpp>
36 #include <boost/numeric/ublas/operation.hpp>
37 
38 #include <gsl/gsl_combination.h>
39 
40 #include <o2scl/err_hnd.h>
41 #include <o2scl/vector.h>
42 #include <o2scl/vec_stats.h>
43 #include <o2scl/linear_solver.h>
44 #include <o2scl/columnify.h>
45 #include <o2scl/cholesky.h>
46 
47 #ifndef DOXYGEN_NO_O2NS
48 namespace o2scl {
49 #endif
50 
51  /** \brief Multi-dimensional interpolation by kriging
52 
53  \note This class assumes that the function specified in the
54  call to set_data() is the same as that passed to the
55  eval() functions. If this is not the case, the
56  behavior of this class is undefined.
57 
58  \note Experimental.
59  */
60  template<class vec_t=boost::numeric::ublas::vector<double>,
61  class mat_col_t=boost::numeric::ublas::matrix_column<
62  boost::numeric::ublas::matrix<double> >,
63  class covar_func_t=std::function<double(const vec_t &,const vec_t &)> >
64  class interpm_krige {
65 
66  public:
67 
71 
72  protected:
73 
74  /** \brief Inverse covariance matrix times function vector
75  */
76  std::vector<ubvector> Kinvf;
77 
78  public:
79 
80  interpm_krige() {
81  data_set=false;
82  verbose=0;
83  }
84 
85  /** \brief Verbosity parameter (default 0)
86  */
87  int verbose;
88 
89  /** \brief Initialize the data for the interpolation
90  */
91  template<class vec_vec_t, class vec_vec2_t>
92  void set_data(size_t n_in, size_t n_out, size_t n_points,
93  vec_vec_t &x, vec_vec2_t &y,
94  std::vector<covar_func_t> &fcovar) {
95 
96  if (n_points<3) {
97  O2SCL_ERR2("Must provide at least three points in ",
98  "interpm_krige::set_data()",exc_efailed);
99  }
100  if (n_in<1) {
101  O2SCL_ERR2("Must provide at least one input column in ",
102  "interpm_krige::set_data()",exc_efailed);
103  }
104  if (n_out<1) {
105  O2SCL_ERR2("Must provide at least one output column in ",
106  "interpm_krige::set_data()",exc_efailed);
107  }
108  np=n_points;
109  nd_in=n_in;
110  nd_out=n_out;
111  ptrs_x.resize(n_points);
112  for(size_t i=0;i<n_points;i++) {
113  if (x[i].size()!=n_in) {
114  O2SCL_ERR2("Size of x not correct in ",
115  "interpm_krige::set_data().",o2scl::exc_efailed);
116  }
117  std::swap(ptrs_x[i],x[i]);
118  }
119  data_set=true;
120 
121  if (verbose>0) {
122  std::cout << "interpm_krige::set_data() : Using " << n_points
123  << " points with " << nd_in << " input variables and\n\t"
124  << nd_out << " output variables." << std::endl;
125  }
126 
127  Kinvf.resize(n_out);
128 
129  // Loop over all output functions
130  for(size_t iout=0;iout<n_out;iout++) {
131 
132  // Construct the KXX matrix
133  ubmatrix KXX(n_points,n_points);
134  for(size_t irow=0;irow<n_points;irow++) {
135  for(size_t icol=0;icol<n_points;icol++) {
136  if (irow>icol) {
137  KXX(irow,icol)=KXX(icol,irow);
138  } else {
139  KXX(irow,icol)=fcovar[iout](ptrs_x[irow],ptrs_x[icol]);
140  }
141  }
142  }
143 
144  o2scl_linalg::cholesky_decomp(n_points,KXX);
145  ubmatrix &inv_KXX=KXX;
146  o2scl_linalg::cholesky_invert<ubmatrix>(n_points,inv_KXX);
147 
148  // Inverse covariance matrix times function vector
149  Kinvf[iout].resize(n_points);
150  boost::numeric::ublas::axpy_prod(inv_KXX,y[iout],Kinvf[iout],true);
151 
152  }
153 
154  return;
155  }
156 
157  /** \brief Given covariance function \c fcovar and input vector \c x
158  store the result of the interpolation in \c y
159  */
160  template<class vec2_t, class vec3_t>
161  void eval(const vec2_t &x, vec3_t &y,
162  std::vector<covar_func_t> &fcovar) const {
163 
164  if (data_set==false) {
165  O2SCL_ERR("Data not set in interpm_krige::eval().",
166  exc_einval);
167  }
168 
169  y.resize(nd_out);
170  for(size_t iout=0;iout<nd_out;iout++) {
171  y[iout]=0.0;
172  for(size_t ipoints=0;ipoints<np;ipoints++) {
173  y[iout]+=fcovar[iout](x,ptrs_x[ipoints])*Kinvf[iout][ipoints];
174  }
175  }
176 
177  return;
178 
179  }
180 
181 #ifndef DOXYGEN_INTERNAL
182 
183  protected:
184 
185  /// The number of points
186  size_t np;
187  /// The number of dimensions of the inputs
188  size_t nd_in;
189  /// The number of dimensions of the outputs
190  size_t nd_out;
191  /// A vector of pointers holding the data
192  std::vector<vec_t> ptrs_x;
193  /// True if the data has been specified
194  bool data_set;
195 
196 #endif
197 
198  };
199 
200  /** \brief Multi-dimensional interpolation by kriging with
201  nearest-neighbor
202 
203  \note This class assumes that the function specified in the
204  call to set_data() is the same as that passed to the
205  eval() functions. If this is not the case, the
206  behavior of this class is undefined.
207 
208  \note Experimental.
209  */
210  template<class vec_t=boost::numeric::ublas::vector<double>,
211  class mat_col_t=boost::numeric::ublas::matrix_column<
212  boost::numeric::ublas::matrix<double> >,
213  class covar_func_t=std::function<double(const vec_t &,const vec_t &)> >
215 
216  public:
217 
221 
222  interpm_krige_nn() {
223  data_set=false;
224  verbose=0;
225  }
226 
227  /** \brief Verbosity parameter (default 0)
228  */
229  int verbose;
230 
231  /** \brief Initialize the data for the interpolation
232  */
233  template<class vec_vec_t, class vec_vec2_t>
234  void set_data(size_t n_in, size_t n_out, size_t n_points,
235  vec_vec_t &x, vec_vec2_t &y,
236  std::vector<covar_func_t> &fcovar, size_t order) {
237 
238  if (n_points<3) {
239  O2SCL_ERR2("Must provide at least three points in ",
240  "interpm_krige_nn::set_data()",exc_efailed);
241  }
242  if (n_in<1) {
243  O2SCL_ERR2("Must provide at least one input column in ",
244  "interpm_krige_nn::set_data()",exc_efailed);
245  }
246  if (n_out<1) {
247  O2SCL_ERR2("Must provide at least one output column in ",
248  "interpm_krige_nn::set_data()",exc_efailed);
249  }
250  np=n_points;
251  nd_in=n_in;
252  nd_out=n_out;
253  ptrs_x.resize(n_points);
254  n_order=order;
255  for(size_t i=0;i<n_points;i++) {
256  if (x[i].size()!=n_in) {
257  O2SCL_ERR2("Size of x not correct in ",
258  "interpm_krige_nn::set_data().",o2scl::exc_efailed);
259  }
260  std::swap(ptrs_x[i],x[i]);
261  }
262  ptrs_y.resize(n_out);
263  for(size_t i=0;i<n_out;i++) {
264  if (y[i].size()!=n_points) {
265  O2SCL_ERR2("Size of y not correct in ",
266  "interpm_krige_nn::set_data().",o2scl::exc_efailed);
267  }
268  std::swap(ptrs_y[i],y[i]);
269  }
270  data_set=true;
271 
272  if (verbose>0) {
273  std::cout << "interpm_krige_nn::set_data() : Using " << n_points
274  << " points with " << nd_in << " input variables and\n\t"
275  << nd_out << " output variables and order "
276  << n_order << " ." << std::endl;
277  }
278 
279  return;
280  }
281 
282  /** \brief Given covariance function \c fcovar and input vector \c x
283  store the result of the interpolation in \c y
284  */
285  template<class vec2_t, class vec3_t>
286  void eval(const vec2_t &x, vec3_t &y,
287  std::vector<covar_func_t> &fcovar) const {
288 
289  if (data_set==false) {
290  O2SCL_ERR("Data not set in interpm_krige_nn::eval().",
291  exc_einval);
292  }
293 
294  y.resize(nd_out);
295 
296  // Loop over all output functions
297  for(size_t iout=0;iout<nd_out;iout++) {
298 
299  // Find points closest to requested point, as defined
300  // by the negative covariance for this output function
301  ubvector dists(np);
302  for(size_t ip=0;ip<np;ip++) {
303  dists[ip]=-fcovar[iout](x,ptrs_x[ip]);
304  }
305 
306  // Empty index vector (resized by the vector_smallest_index
307  // function)
308  ubvector_size_t index;
309  o2scl::vector_smallest_index<ubvector,double,ubvector_size_t>
310  (np,dists,n_order,index);
311 
312  // Construct subset of function values for nearest neighbors
313  ubvector func(n_order);
314  for(size_t io=0;io<n_order;io++) {
315  func[io]=ptrs_y[iout][index[io]];
316  }
317 
318  // Construct the nearest neighbor KXX matrix
319  ubmatrix KXX(n_order,n_order);
320  for(size_t irow=0;irow<n_order;irow++) {
321  for(size_t icol=0;icol<n_order;icol++) {
322  if (irow>icol) {
323  KXX(irow,icol)=KXX(icol,irow);
324  } else {
325  KXX(irow,icol)=fcovar[iout](ptrs_x[index[irow]],
326  ptrs_x[index[icol]]);
327  }
328  }
329  }
330 
331  // Construct the inverse of KXX
332  o2scl_linalg::cholesky_decomp(n_order,KXX);
333  ubmatrix &inv_KXX=KXX;
334  o2scl_linalg::cholesky_invert<ubmatrix>(n_order,inv_KXX);
335 
336  // Inverse covariance matrix times function vector
337  ubvector Kinvf(n_order);
338  boost::numeric::ublas::axpy_prod(inv_KXX,func,Kinvf,true);
339 
340  // Comput the final result
341  y[iout]=0.0;
342  for(size_t ipoints=0;ipoints<n_order;ipoints++) {
343  y[iout]+=-dists[index[ipoints]]*Kinvf[ipoints];
344  }
345 
346  }
347 
348  return;
349 
350  }
351 
352  /** \brief Find a set of linearly independent points
353 
354  Given a point \c x, a covariance function
355  \c fcovar, the index of the output function
356  \c iout, and an array specifying the closest points
357  \c index, this function produces a list of
358  */
359  template<class vec2_t>
360  void find_lin_indep(const vec2_t &x, size_t iout,
361  std::vector<covar_func_t> &fcovar,
362  ubvector_size_t &index,
363  ubvector_size_t &indep) const {
364 
365  if (x.size()<nd_in || index.size()<np || indep.size()<n_order
366  || iout>=nd_out) {
367  std::cout << x.size() << " " << nd_in << std::endl;
368  std::cout << index.size() << " " << np << std::endl;
369  std::cout << indep.size() << " " << n_order << std::endl;
370  std::cout << iout << " " << nd_out << std::endl;
371  O2SCL_ERR("Vectors not of correct size in find_lin_indep().",
373  }
374 
375  bool done=false;
376  while (done==false) {
377 
378  // Construct subset of function values for nearest neighbors
379  ubvector func(n_order);
380  for(size_t io=0;io<n_order;io++) {
381  func[io]=ptrs_y[iout][index[indep[io]]];
382  }
383 
384  // Construct the nearest neighbor KXX matrix
385  ubmatrix KXX(n_order,n_order);
386  for(size_t irow=0;irow<n_order;irow++) {
387  for(size_t icol=0;icol<n_order;icol++) {
388  if (irow>icol) {
389  KXX(irow,icol)=KXX(icol,irow);
390  } else {
391  KXX(irow,icol)=fcovar[iout](ptrs_x[index[indep[irow]]],
392  ptrs_x[index[indep[icol]]]);
393  }
394  }
395  }
396 
397  // Construct the inverse of KXX
398  int cret=o2scl_linalg::cholesky_decomp(n_order,KXX);
399  if (cret==0) {
400  done=true;
401  } else {
402  if (verbose>1) {
403  std::cout << "Finding new independent rows." << std::endl;
404  for(size_t j=0;j<n_order;j++) {
405  std::cout << indep[j] << " " << KXX(j,j) << std::endl;
406  }
407  }
408  size_t max=o2scl::vector_max_value<ubvector_size_t,
409  double>(indep);
410  if (verbose>1) {
411  std::cout << "Max is: " << max << std::endl;
412  }
413  for(size_t j=0;j<n_order;j++) {
414  if (KXX(j,j)==0.0) {
415  if (max+1<np) {
416  if (verbose>1) {
417  std::cout << "Entry " << j << " is zero so replacing "
418  << "entry with " << max+1 << std::endl;
419  }
420  indep[j]=max+1;
421  max++;
422  } else {
423  O2SCL_ERR3("Failed to find set of independent points ",
424  "in interpm_krige_nn::find_lin_indep",
425  "(const vec2_t &, size_t).",
427  }
428  }
429  }
430  }
431 
432  }
433  return;
434  }
435 
436  /** \brief Given covariance function \c fcovar and input vector \c x
437  return the result of the interpolation for function with
438  index \c iout
439  */
440  template<class vec2_t>
441  double eval(const vec2_t &x, size_t iout,
442  std::vector<covar_func_t> &fcovar) const {
443 
444  if (data_set==false) {
445  O2SCL_ERR("Data not set in interpm_krige_nn::eval().",
446  exc_einval);
447  }
448 
449  double ret;
450 
451  // Find points closest to requested point, as defined
452  // by the negative covariance for this output function
453  ubvector dists(np);
454  for(size_t ip=0;ip<np;ip++) {
455  dists[ip]=-fcovar[iout](x,ptrs_x[ip]);
456  }
457 
458  ubvector_size_t index(np);
459  o2scl::vector_sort_index<ubvector,ubvector_size_t>(np,dists,index);
460 
461  // Vector for storing the indexes in the index array which
462  // will store the closest n_order points which are
463  // linearly independent
464  ubvector_size_t indep(n_order);
465  for(size_t io=0;io<n_order;io++) {
466  indep[io]=io;
467  }
468 
469  find_lin_indep(x,iout,fcovar,index,indep);
470 
471  // Construct subset of function values for nearest neighbors
472  ubvector func(n_order);
473  for(size_t io=0;io<n_order;io++) {
474  func[io]=ptrs_y[iout][index[indep[io]]];
475  }
476 
477  // Construct the nearest neighbor KXX matrix
478  ubmatrix KXX(n_order,n_order);
479  for(size_t irow=0;irow<n_order;irow++) {
480  for(size_t icol=0;icol<n_order;icol++) {
481  if (irow>icol) {
482  KXX(irow,icol)=KXX(icol,irow);
483  } else {
484  KXX(irow,icol)=fcovar[iout](ptrs_x[index[indep[irow]]],
485  ptrs_x[index[indep[icol]]]);
486  }
487  }
488  }
489 
490  // Construct the inverse of KXX
491  o2scl_linalg::cholesky_decomp(n_order,KXX);
492  ubmatrix &inv_KXX=KXX;
493  o2scl_linalg::cholesky_invert<ubmatrix>(n_order,inv_KXX);
494 
495  // Inverse covariance matrix times function vector
496  ubvector Kinvf(n_order);
497  boost::numeric::ublas::axpy_prod(inv_KXX,func,Kinvf,true);
498 
499  // Comput the final result
500  ret=0.0;
501  for(size_t ipoints=0;ipoints<n_order;ipoints++) {
502  ret+=-dists[index[indep[ipoints]]]*Kinvf[ipoints];
503  }
504 
505  return ret;
506 
507  }
508 
509  /** \brief Compute a quality factor for interpolation
510  using jackknife resampling
511  */
512  template<class vec2_t>
513  double eval_jackknife(const vec2_t &x, size_t iout,
514  std::vector<covar_func_t> &fcovar) const {
515 
516  if (data_set==false) {
517  O2SCL_ERR("Data not set in interpm_krige_nn::eval_jackknife().",
518  exc_einval);
519  }
520 
521  double qual=0.0;
522 
523  // Interpolated function value inside jackknife loop
524  double ytmp;
525 
526  // For a distance measurement, just use the the negative
527  // covariance for this output function
528  ubvector dists(np);
529  for(size_t ip=0;ip<np;ip++) {
530  dists[ip]=-fcovar[iout](x,ptrs_x[ip]);
531  }
532 
533  // Create an index array which sorts by distance
534  ubvector_size_t index(np);
535  o2scl::vector_sort_index<ubvector,ubvector_size_t>(np,dists,index);
536 
537  // Vector for storing the indexes in the index array which
538  // will store the closest n_order points which are
539  // linearly independent
540  ubvector_size_t indep(n_order);
541  for(size_t io=0;io<n_order;io++) {
542  indep[io]=io;
543  }
544 
545  // -------------------------------------------------------------
546  // Before the jackknife loop, we want to create a full
547  // set of n_order linearly independent points
548 
549  find_lin_indep(x,iout,fcovar,index,indep);
550 
551  // -------------------------------------------------------------
552  // Now, the jackknife loop, removing one point at a time
553 
554  for(size_t jk=0;jk<n_order;jk++) {
555 
556  if (verbose>0) {
557  std::cout << "Jackknife: " << jk << " matching function value "
558  << ptrs_y[iout][index[jk]] << std::endl;
559  }
560 
561  ubvector_size_t indep_jk;
562  vector_copy_jackknife(indep,jk,indep_jk);
563 
564  // Construct subset of function values for nearest neighbors
565  ubvector func(n_order-1);
566  for(size_t io=0;io<n_order-1;io++) {
567  func[io]=ptrs_y[iout][index[indep_jk[io]]];
568  }
569 
570  // Construct the nearest neighbor KXX matrix
571  ubmatrix KXX(n_order-1,n_order-1);
572  for(size_t irow=0;irow<n_order-1;irow++) {
573  for(size_t icol=0;icol<n_order-1;icol++) {
574  if (irow>icol) {
575  KXX(irow,icol)=KXX(icol,irow);
576  } else {
577  KXX(irow,icol)=fcovar[iout](ptrs_x[index[indep_jk[irow]]],
578  ptrs_x[index[indep_jk[icol]]]);
579  }
580  }
581  }
582 
583  // Construct the inverse of KXX
584  o2scl_linalg::cholesky_decomp(n_order-1,KXX);
585  ubmatrix &inv_KXX=KXX;
586  o2scl_linalg::cholesky_invert<ubmatrix>(n_order-1,inv_KXX);
587 
588  // Inverse covariance matrix times function vector
589  ubvector Kinvf(n_order-1);
590  boost::numeric::ublas::axpy_prod(inv_KXX,func,Kinvf,true);
591 
592  // Comput the final result
593  ytmp=0.0;
594  for(size_t ipoints=0;ipoints<n_order-1;ipoints++) {
595  ytmp+=-dists[index[indep_jk[ipoints]]]*Kinvf[ipoints];
596  }
597 
598  // Add the squared deviation to y[iout]
599  qual+=pow(ptrs_y[iout][index[jk]]-ytmp,2.0);
600 
601  if (verbose>0) {
602  std::cout << "Original value: "
603  << ptrs_y[iout][index[jk]] << " interpolated: "
604  << ytmp << std::endl;
605  }
606 
607  // End of jackknife loop
608  }
609 
610  return qual;
611  }
612 
613 
614 #ifndef DOXYGEN_INTERNAL
615 
616  protected:
617 
618  /// The order of the interpolation (specified by \ref set_data() )
619  size_t n_order;
620  /// The number of points
621  size_t np;
622  /// The number of dimensions of the inputs
623  size_t nd_in;
624  /// The number of dimensions of the outputs
625  size_t nd_out;
626  /// A vector of pointers holding the data
627  std::vector<vec_t> ptrs_x;
628  /// A vector of pointers holding the data
629  std::vector<vec_t> ptrs_y;
630  /// True if the data has been specified
631  bool data_set;
632 
633 #endif
634 
635  };
636 
637 #ifndef DOXYGEN_NO_O2NS
638 }
639 #endif
640 
641 #endif
642 
643 
644 
bool data_set
True if the data has been specified.
size_t nd_out
The number of dimensions of the outputs.
void eval(const vec2_t &x, vec3_t &y, std::vector< covar_func_t > &fcovar) const
Given covariance function fcovar and input vector x store the result of the interpolation in y...
void vector_copy_jackknife(const vec_t &src, size_t iout, vec2_t &dest)
From a given vector, create a new vector by removing a specified element.
Definition: vector.h:2367
Multi-dimensional interpolation by kriging with nearest-neighbor.
double eval_jackknife(const vec2_t &x, size_t iout, std::vector< covar_func_t > &fcovar) const
Compute a quality factor for interpolation using jackknife resampling.
The main O<span style=&#39;position: relative; top: 0.3em; font-size: 0.8em&#39;>2</span>scl O$_2$scl names...
Definition: anneal.h:42
int verbose
Verbosity parameter (default 0)
void set_data(size_t n_in, size_t n_out, size_t n_points, vec_vec_t &x, vec_vec2_t &y, std::vector< covar_func_t > &fcovar, size_t order)
Initialize the data for the interpolation.
void find_lin_indep(const vec2_t &x, size_t iout, std::vector< covar_func_t > &fcovar, ubvector_size_t &index, ubvector_size_t &indep) const
Find a set of linearly independent points.
size_t nd_out
The number of dimensions of the outputs.
Multi-dimensional interpolation by kriging.
Definition: interpm_krige.h:64
void set_data(size_t n_in, size_t n_out, size_t n_points, vec_vec_t &x, vec_vec2_t &y, std::vector< covar_func_t > &fcovar)
Initialize the data for the interpolation.
Definition: interpm_krige.h:92
invalid argument supplied by user
Definition: err_hnd.h:59
int verbose
Verbosity parameter (default 0)
Definition: interpm_krige.h:87
size_t nd_in
The number of dimensions of the inputs.
std::vector< vec_t > ptrs_x
A vector of pointers holding the data.
size_t n_order
The order of the interpolation (specified by set_data() )
double eval(const vec2_t &x, size_t iout, std::vector< covar_func_t > &fcovar) const
Given covariance function fcovar and input vector x return the result of the interpolation for functi...
generic failure
Definition: err_hnd.h:61
#define O2SCL_ERR3(d, d2, d3, n)
Set an error, three-string version.
Definition: err_hnd.h:286
void eval(const vec2_t &x, vec3_t &y, std::vector< covar_func_t > &fcovar) const
Given covariance function fcovar and input vector x store the result of the interpolation in y...
data_t vector_max_value(size_t n, const vec_t &data)
Compute the maximum of the first n elements of a vector.
Definition: vector.h:1124
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
Definition: err_hnd.h:281
bool data_set
True if the data has been specified.
int cholesky_decomp(const size_t M, mat_t &A, bool err_on_fail=true)
Compute the in-place Cholesky decomposition of a symmetric positive-definite square matrix...
Definition: cholesky_base.h:62
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
Definition: err_hnd.h:273
std::vector< vec_t > ptrs_y
A vector of pointers holding the data.
std::vector< ubvector > Kinvf
Inverse covariance matrix times function vector.
Definition: interpm_krige.h:76
size_t np
The number of points.
size_t np
The number of points.
std::vector< vec_t > ptrs_x
A vector of pointers holding the data.
size_t nd_in
The number of dimensions of the inputs.

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).