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 
46 #ifndef DOXYGEN_NO_O2NS
47 namespace o2scl {
48 #endif
49 
50  /** \brief Multi-dimensional interpolation by kriging
51 
52  \note This class assumes that the function specified in the
53  call to set_data() is the same as that passed to the
54  eval() functions. If this is not the case, the
55  behavior of this class is undefined.
56 
57  \note Experimental.
58  */
59  template<class vec_t=boost::numeric::ublas::vector<double>,
60  class mat_col_t=boost::numeric::ublas::matrix_column<
61  boost::numeric::ublas::matrix<double> >,
62  class covar_func_t=std::function<double(const vec_t &,const vec_t &)> >
63  class interpm_krige {
64 
65  public:
66 
70 
71  protected:
72 
73  /** \brief Inverse covariance matrix times function vector
74  */
75  std::vector<ubvector> Kinvf;
76 
77  public:
78 
79  interpm_krige() {
80  data_set=false;
81  verbose=0;
82  }
83 
84  /** \brief Verbosity parameter (default 0)
85  */
86  int verbose;
87 
88  /** \brief Initialize the data for the interpolation
89  */
90  template<class vec_vec_t, class vec_vec2_t>
91  void set_data(size_t n_in, size_t n_out, size_t n_points,
92  vec_vec_t &x, vec_vec2_t &y,
93  std::vector<covar_func_t> &fcovar) {
94 
95  if (n_points<3) {
96  O2SCL_ERR2("Must provide at least three points in ",
97  "interpm_krige::set_data()",exc_efailed);
98  }
99  if (n_in<1) {
100  O2SCL_ERR2("Must provide at least one input column in ",
101  "interpm_krige::set_data()",exc_efailed);
102  }
103  if (n_out<1) {
104  O2SCL_ERR2("Must provide at least one output column in ",
105  "interpm_krige::set_data()",exc_efailed);
106  }
107  np=n_points;
108  nd_in=n_in;
109  nd_out=n_out;
110  ptrs_x.resize(n_points);
111  for(size_t i=0;i<n_points;i++) {
112  if (x[i].size()!=n_in) {
113  O2SCL_ERR2("Size of x not correct in ",
114  "interpm_krige::set_data().",o2scl::exc_efailed);
115  }
116  std::swap(ptrs_x[i],x[i]);
117  }
118  data_set=true;
119 
120  if (verbose>0) {
121  std::cout << "interpm_krige::set_data() : Using " << n_points
122  << " points with " << nd_in << " input variables and\n\t"
123  << nd_out << " output variables." << std::endl;
124  }
125 
126  Kinvf.resize(n_out);
127 
128  // Loop over all output functions
129  for(size_t iout=0;iout<n_out;iout++) {
130 
131  // Construct the KXX matrix
132  ubmatrix KXX(n_points,n_points);
133  for(size_t irow=0;irow<n_points;irow++) {
134  for(size_t icol=0;icol<n_points;icol++) {
135  if (irow>icol) {
136  KXX(irow,icol)=KXX(icol,irow);
137  } else {
138  KXX(irow,icol)=fcovar[iout](ptrs_x[irow],ptrs_x[icol]);
139  }
140  }
141  }
142 
143  // Construct the inverse of KXX
144  ubmatrix inv_KXX(n_points,n_points);
145  o2scl::permutation p(n_points);
146  int signum;
147  if (verbose>0) {
148  std::cout << "interpm_krige::set_data() : "
149  << "LU decompose and invert " << iout+1 << " of " << n_out
150  << std::endl;
151  }
152  o2scl_linalg::LU_decomp(n_points,KXX,p,signum);
153  if (o2scl_linalg::diagonal_has_zero(n_points,KXX)) {
154  O2SCL_ERR2("KXX matrix is singular in ",
155  "interpm_krige::set_data().",
157  }
158  o2scl_linalg::LU_invert<ubmatrix,ubmatrix,mat_col_t>
159  (n_points,KXX,p,inv_KXX);
160 
161  // Inverse covariance matrix times function vector
162  Kinvf[iout].resize(n_points);
163  boost::numeric::ublas::axpy_prod(inv_KXX,y[iout],Kinvf[iout],true);
164 
165  }
166 
167  return;
168  }
169 
170  /** \brief Perform the interpolation
171  */
172  template<class vec2_t, class vec3_t>
173  void eval(const vec2_t &x, vec3_t &y,
174  std::vector<covar_func_t> &fcovar) const {
175 
176  if (data_set==false) {
177  O2SCL_ERR("Data not set in interpm_krige::eval().",
178  exc_einval);
179  }
180 
181  y.resize(nd_out);
182  for(size_t iout=0;iout<nd_out;iout++) {
183  y[iout]=0.0;
184  for(size_t ipoints=0;ipoints<np;ipoints++) {
185  y[iout]+=fcovar[iout](x,ptrs_x[ipoints])*Kinvf[iout][ipoints];
186  }
187  }
188 
189  return;
190 
191  }
192 
193 #ifndef DOXYGEN_INTERNAL
194 
195  protected:
196 
197  /// The number of points
198  size_t np;
199  /// The number of dimensions of the inputs
200  size_t nd_in;
201  /// The number of dimensions of the outputs
202  size_t nd_out;
203  /// A vector of pointers holding the data
204  std::vector<vec_t> ptrs_x;
205  /// True if the data has been specified
206  bool data_set;
207  /// Number of points to include in each interpolation (default 3)
208  size_t order;
209 
210 #endif
211 
212  };
213 
214  /** \brief Multi-dimensional interpolation by kriging with
215  nearest-neighbor
216 
217  \note This class assumes that the function specified in the
218  call to set_data() is the same as that passed to the
219  eval() functions. If this is not the case, the
220  behavior of this class is undefined.
221 
222  \note Experimental.
223  */
224  template<class vec_t=boost::numeric::ublas::vector<double>,
225  class mat_col_t=boost::numeric::ublas::matrix_column<
226  boost::numeric::ublas::matrix<double> >,
227  class covar_func_t=std::function<double(const vec_t &,const vec_t &)> >
229 
230  public:
231 
235 
236  interpm_krige_nn() {
237  data_set=false;
238  verbose=0;
239  }
240 
241  /** \brief Verbosity parameter (default 0)
242  */
243  int verbose;
244 
245  /** \brief Initialize the data for the interpolation
246  */
247  template<class vec_vec_t, class vec_vec2_t>
248  void set_data(size_t n_in, size_t n_out, size_t n_points,
249  vec_vec_t &x, vec_vec2_t &y,
250  std::vector<covar_func_t> &fcovar, size_t order) {
251 
252  if (n_points<3) {
253  O2SCL_ERR2("Must provide at least three points in ",
254  "interpm_krige::set_data()",exc_efailed);
255  }
256  if (n_in<1) {
257  O2SCL_ERR2("Must provide at least one input column in ",
258  "interpm_krige::set_data()",exc_efailed);
259  }
260  if (n_out<1) {
261  O2SCL_ERR2("Must provide at least one output column in ",
262  "interpm_krige::set_data()",exc_efailed);
263  }
264  np=n_points;
265  nd_in=n_in;
266  nd_out=n_out;
267  ptrs_x.resize(n_points);
268  norder=order;
269  for(size_t i=0;i<n_points;i++) {
270  if (x[i].size()!=n_in) {
271  O2SCL_ERR2("Size of x not correct in ",
272  "interpm_krige_nn::set_data().",o2scl::exc_efailed);
273  }
274  std::swap(ptrs_x[i],x[i]);
275  }
276  ptrs_y.resize(n_out);
277  for(size_t i=0;i<n_out;i++) {
278  if (y[i].size()!=n_points) {
279  O2SCL_ERR2("Size of y not correct in ",
280  "interpm_krige_nn::set_data().",o2scl::exc_efailed);
281  }
282  std::swap(ptrs_y[i],y[i]);
283  }
284  data_set=true;
285 
286  if (verbose>0) {
287  std::cout << "interpm_krige_nn::set_data() : Using " << n_points
288  << " points with " << nd_in << " input variables and\n\t"
289  << nd_out << " output variables and order "
290  << norder << " ." << std::endl;
291  }
292 
293  return;
294  }
295 
296  /** \brief Perform the interpolation
297  */
298  template<class vec2_t, class vec3_t>
299  void eval(const vec2_t &x, vec3_t &y,
300  std::vector<covar_func_t> &fcovar) const {
301 
302  if (data_set==false) {
303  O2SCL_ERR("Data not set in interpm_krige::eval().",
304  exc_einval);
305  }
306 
307  y.resize(nd_out);
308 
309  // Loop over all output functions
310  for(size_t iout=0;iout<nd_out;iout++) {
311 
312 
313  // Find points closest to requested point, as defined
314  // by the negative covariance for this output function
315  ubvector dists(np);
316  for(size_t ip=0;ip<np;ip++) {
317  dists[ip]=-fcovar[iout](x,ptrs_x[ip]);
318  }
319 
320  // Empty index vector (resized by the vector_smallest_index
321  // function)
322  ubvector_size_t index;
323  o2scl::vector_smallest_index<ubvector,double,ubvector_size_t>
324  (np,dists,norder,index);
325 
326  // Construct subset of function values for nearest neighbors
327  ubvector func(norder);
328  for(size_t io=0;io<norder;io++) {
329  func[io]=ptrs_y[iout][index[io]];
330  }
331 
332  // Construct the nearest neighbor KXX matrix
333  ubmatrix KXX(norder,norder);
334  for(size_t irow=0;irow<norder;irow++) {
335  for(size_t icol=0;icol<norder;icol++) {
336  if (irow>icol) {
337  KXX(irow,icol)=KXX(icol,irow);
338  } else {
339  KXX(irow,icol)=fcovar[iout](ptrs_x[index[irow]],
340  ptrs_x[index[icol]]);
341  }
342  }
343  }
344 
345  // Construct the inverse of KXX
346  ubmatrix inv_KXX(norder,norder);
347  o2scl::permutation p(norder);
348  int signum;
349  o2scl_linalg::LU_decomp(norder,KXX,p,signum);
350  if (o2scl_linalg::diagonal_has_zero(norder,KXX)) {
351  O2SCL_ERR2("KXX matrix is singular in ",
352  "interpm_krige_nn::eval().",
354  }
355  o2scl_linalg::LU_invert<ubmatrix,ubmatrix,mat_col_t>
356  (norder,KXX,p,inv_KXX);
357 
358  // Inverse covariance matrix times function vector
359  ubvector Kinvf(norder);
360  boost::numeric::ublas::axpy_prod(inv_KXX,func,Kinvf,true);
361 
362  // Comput the final result
363  y[iout]=0.0;
364  for(size_t ipoints=0;ipoints<norder;ipoints++) {
365  y[iout]+=-dists[index[ipoints]]*Kinvf[ipoints];
366  }
367 
368  }
369 
370  return;
371 
372  }
373 
374  /** \brief Find a set of linearly independent points
375  */
376  template<class vec2_t>
377  void find_lin_indep(const vec2_t &x, size_t iout,
378  std::vector<covar_func_t> &fcovar,
379  ubvector_size_t &index,
380  ubvector_size_t &indep) const {
381 
382  if (x.size()<nd_in || index.size()<np || indep.size()<norder
383  || iout>=nd_out) {
384  std::cout << x.size() << " " << nd_in << std::endl;
385  std::cout << index.size() << " " << np << std::endl;
386  std::cout << indep.size() << " " << norder << std::endl;
387  std::cout << iout << " " << nd_out << std::endl;
388  O2SCL_ERR("Vectors not of correct size in find_lin_indep().",
390  }
391 
392  bool done=false;
393  while (done==false) {
394 
395  // Construct subset of function values for nearest neighbors
396  ubvector func(norder);
397  for(size_t io=0;io<norder;io++) {
398  func[io]=ptrs_y[iout][index[indep[io]]];
399  }
400 
401  // Construct the nearest neighbor KXX matrix
402  ubmatrix KXX(norder,norder);
403  for(size_t irow=0;irow<norder;irow++) {
404  for(size_t icol=0;icol<norder;icol++) {
405  if (irow>icol) {
406  KXX(irow,icol)=KXX(icol,irow);
407  } else {
408  KXX(irow,icol)=fcovar[iout](ptrs_x[index[indep[irow]]],
409  ptrs_x[index[indep[icol]]]);
410  }
411  }
412  }
413 
414  // Construct the inverse of KXX
415  o2scl::permutation p(norder);
416  int signum;
417  o2scl_linalg::LU_decomp(norder,KXX,p,signum);
418  if (!o2scl_linalg::diagonal_has_zero(norder,KXX)) {
419  done=true;
420  } else {
421  if (verbose>1) {
422  std::cout << "Finding new independent rows." << std::endl;
423  for(size_t j=0;j<norder;j++) {
424  std::cout << indep[j] << " " << KXX(j,j) << std::endl;
425  }
426  }
427  size_t max=o2scl::vector_max_value<ubvector_size_t,
428  double>(indep);
429  if (verbose>1) {
430  std::cout << "Max is: " << max << std::endl;
431  }
432  for(size_t j=0;j<norder;j++) {
433  if (KXX(j,j)==0.0) {
434  if (max+1<np) {
435  if (verbose>1) {
436  std::cout << "Entry " << j << " is zero so replacing "
437  << "entry with " << max+1 << std::endl;
438  }
439  indep[j]=max+1;
440  max++;
441  } else {
442  O2SCL_ERR3("Failed to find set of independent points ",
443  "in interpm_krige_nn::find_lin_indep",
444  "(const vec2_t &, size_t).",
446  }
447  }
448  }
449  }
450 
451  }
452  return;
453  }
454 
455  /** \brief Perform the interpolation
456  */
457  template<class vec2_t>
458  double eval(const vec2_t &x, size_t iout,
459  std::vector<covar_func_t> &fcovar) const {
460 
461  if (data_set==false) {
462  O2SCL_ERR("Data not set in interpm_krige_nn::eval().",
463  exc_einval);
464  }
465 
466  double ret;
467 
468  // Find points closest to requested point, as defined
469  // by the negative covariance for this output function
470  ubvector dists(np);
471  for(size_t ip=0;ip<np;ip++) {
472  dists[ip]=-fcovar[iout](x,ptrs_x[ip]);
473  }
474 
475  ubvector_size_t index(np);
476  o2scl::vector_sort_index<ubvector,ubvector_size_t>(np,dists,index);
477 
478  // Vector for storing the indexes in the index array which
479  // will store the closest norder points which are
480  // linearly independent
481  ubvector_size_t indep(norder);
482  for(size_t io=0;io<norder;io++) {
483  indep[io]=io;
484  }
485 
486  find_lin_indep(x,iout,fcovar,index,indep);
487 
488  // Construct subset of function values for nearest neighbors
489  ubvector func(norder);
490  for(size_t io=0;io<norder;io++) {
491  func[io]=ptrs_y[iout][index[indep[io]]];
492  }
493 
494  // Construct the nearest neighbor KXX matrix
495  ubmatrix KXX(norder,norder);
496  for(size_t irow=0;irow<norder;irow++) {
497  for(size_t icol=0;icol<norder;icol++) {
498  if (irow>icol) {
499  KXX(irow,icol)=KXX(icol,irow);
500  } else {
501  KXX(irow,icol)=fcovar[iout](ptrs_x[index[indep[irow]]],
502  ptrs_x[index[indep[icol]]]);
503  }
504  }
505  }
506 
507  // Construct the inverse of KXX
508  ubmatrix inv_KXX(norder,norder);
509  o2scl::permutation p(norder);
510  int signum;
511  o2scl_linalg::LU_decomp(norder,KXX,p,signum);
512  if (o2scl_linalg::diagonal_has_zero(norder,KXX)) {
513  O2SCL_ERR2("KXX matrix is singular in ",
514  "interpm_krige_nn::eval().",
516  }
517  o2scl_linalg::LU_invert<ubmatrix,ubmatrix,mat_col_t>
518  (norder,KXX,p,inv_KXX);
519 
520  // Inverse covariance matrix times function vector
521  ubvector Kinvf(norder);
522  boost::numeric::ublas::axpy_prod(inv_KXX,func,Kinvf,true);
523 
524  // Comput the final result
525  ret=0.0;
526  for(size_t ipoints=0;ipoints<norder;ipoints++) {
527  ret+=-dists[index[indep[ipoints]]]*Kinvf[ipoints];
528  }
529 
530  return ret;
531 
532  }
533 
534  /** \brief Use jackknife to match interpolated to real function
535  values
536  */
537  template<class vec2_t>
538  double eval_jackknife(const vec2_t &x, size_t iout,
539  std::vector<covar_func_t> &fcovar) const {
540 
541  if (data_set==false) {
542  O2SCL_ERR("Data not set in interpm_krige::eval().",
543  exc_einval);
544  }
545 
546  double qual=0.0;
547 
548  // Interpolated function value inside jackknife loop
549  double ytmp;
550 
551  // For a distance measurement, just use the the negative
552  // covariance for this output function
553  ubvector dists(np);
554  for(size_t ip=0;ip<np;ip++) {
555  dists[ip]=-fcovar[iout](x,ptrs_x[ip]);
556  }
557 
558  // Create an index array which sorts by distance
559  ubvector_size_t index(np);
560  o2scl::vector_sort_index<ubvector,ubvector_size_t>(np,dists,index);
561 
562  // Vector for storing the indexes in the index array which
563  // will store the closest norder points which are
564  // linearly independent
565  ubvector_size_t indep(norder);
566  for(size_t io=0;io<norder;io++) {
567  indep[io]=io;
568  }
569 
570  // -------------------------------------------------------------
571  // Before the jackknife loop, we want to create a full
572  // set of norder linearly independent points
573 
574  find_lin_indep(x,iout,fcovar,index,indep);
575 
576  // -------------------------------------------------------------
577  // Now, the jackknife loop, removing one point at a time
578 
579  for(size_t jk=0;jk<norder;jk++) {
580 
581  if (verbose>0) {
582  std::cout << "Jackknife: " << jk << " matching function value "
583  << ptrs_y[iout][index[jk]] << std::endl;
584  }
585 
586  ubvector_size_t indep_jk;
587  vector_copy_jackknife(indep,jk,indep_jk);
588 
589  // Construct subset of function values for nearest neighbors
590  ubvector func(norder-1);
591  for(size_t io=0;io<norder-1;io++) {
592  func[io]=ptrs_y[iout][index[indep_jk[io]]];
593  }
594 
595  // Construct the nearest neighbor KXX matrix
596  ubmatrix KXX(norder-1,norder-1);
597  for(size_t irow=0;irow<norder-1;irow++) {
598  for(size_t icol=0;icol<norder-1;icol++) {
599  if (irow>icol) {
600  KXX(irow,icol)=KXX(icol,irow);
601  } else {
602  KXX(irow,icol)=fcovar[iout](ptrs_x[index[indep_jk[irow]]],
603  ptrs_x[index[indep_jk[icol]]]);
604  }
605  }
606  }
607 
608  // Construct the inverse of KXX
609  o2scl::permutation p(norder-1);
610  ubmatrix inv_KXX(norder-1,norder-1);
611  int signum;
612  o2scl_linalg::LU_decomp(norder-1,KXX,p,signum);
613  int cnt=0;
614  while (o2scl_linalg::diagonal_has_zero(norder-1,KXX)) {
615  std::cout << "Second fli run starting." << std::endl;
616  find_lin_indep(x,iout,fcovar,index,indep);
617  std::cout << "Second fli run done." << std::endl;
618  cnt++;
619  if (cnt==10) {
620  O2SCL_ERR3("Failed to find set of independent points ",
621  "in interpm_krige_nn::eval_jackknife",
622  "(const vec2_t &, size_t).",
624  }
625  }
626  o2scl_linalg::LU_invert<ubmatrix,ubmatrix,mat_col_t>
627  (norder-1,KXX,p,inv_KXX);
628 
629  // Inverse covariance matrix times function vector
630  ubvector Kinvf(norder-1);
631  boost::numeric::ublas::axpy_prod(inv_KXX,func,Kinvf,true);
632 
633  // Comput the final result
634  ytmp=0.0;
635  for(size_t ipoints=0;ipoints<norder-1;ipoints++) {
636  ytmp+=-dists[index[indep_jk[ipoints]]]*Kinvf[ipoints];
637  }
638 
639  // Add the squared deviation to y[iout]
640  qual+=pow(ptrs_y[iout][index[jk]]-ytmp,2.0);
641 
642  if (verbose>0) {
643  std::cout << "Original value: "
644  << ptrs_y[iout][index[jk]] << " interpolated: "
645  << ytmp << std::endl;
646  }
647 
648  // End of jackknife loop
649  }
650 
651  return qual;
652  }
653 
654 
655 #ifndef DOXYGEN_INTERNAL
656 
657  protected:
658 
659  /// Desc
660  size_t norder;
661  /// The number of points
662  size_t np;
663  /// The number of dimensions of the inputs
664  size_t nd_in;
665  /// The number of dimensions of the outputs
666  size_t nd_out;
667  /// A vector of pointers holding the data
668  std::vector<vec_t> ptrs_x;
669  /// A vector of pointers holding the data
670  std::vector<vec_t> ptrs_y;
671  /// True if the data has been specified
672  bool data_set;
673  /// Number of points to include in each interpolation (default 3)
674  size_t order;
675 
676 #endif
677 
678  };
679 
680 #ifndef DOXYGEN_NO_O2NS
681 }
682 #endif
683 
684 #endif
685 
686 
687 
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
Perform the interpolation.
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:2299
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
Use jackknife to match interpolated to real function values.
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:63
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:91
invalid argument supplied by user
Definition: err_hnd.h:59
A class for representing permutations.
Definition: permutation.h:70
int LU_decomp(const size_t N, mat_t &A, o2scl::permutation &p, int &signum)
Compute the LU decomposition of the matrix A.
Definition: lu_base.h:86
int verbose
Verbosity parameter (default 0)
Definition: interpm_krige.h:86
size_t nd_in
The number of dimensions of the inputs.
std::vector< vec_t > ptrs_x
A vector of pointers holding the data.
double eval(const vec2_t &x, size_t iout, std::vector< covar_func_t > &fcovar) const
Perform the interpolation.
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
size_t order
Number of points to include in each interpolation (default 3)
void eval(const vec2_t &x, vec3_t &y, std::vector< covar_func_t > &fcovar) const
Perform the interpolation.
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:1078
size_t order
Number of points to include in each interpolation (default 3)
#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.
#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:75
size_t np
The number of points.
int diagonal_has_zero(const size_t N, mat_t &A)
Return 1 if at least one of the elements in the diagonal is zero.
Definition: lu_base.h:54
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).