Jpp  18.2.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JHermiteSpline.hh
Go to the documentation of this file.
1 #ifndef __JTOOLS__JHERMITESPLINE__
2 #define __JTOOLS__JHERMITESPLINE__
3 
4 #include <utility>
5 #include <sstream>
6 #include <cmath>
7 
8 #include "JMath/JZero.hh"
9 #include "JLang/JException.hh"
10 #include "JLang/JClass.hh"
11 #include "JTools/JFunctional.hh"
12 #include "JTools/JDistance.hh"
13 #include "JTools/JResult.hh"
14 #include "JTools/JMapCollection.hh"
15 
16 
17 /**
18  * \author mdejong
19  */
20 
21 namespace JTOOLS {}
22 namespace JPP { using namespace JTOOLS; }
23 
24 namespace JTOOLS {
25 
26  using JLANG::JNoValue;
30 
31 
32  /**
33  * Template base class spline interpolations.
34  *
35  * This class implements the JFunctional interface.
36  *
37  * Note that the data structure of the elements in the collection should have the additional methods:
38  * <pre>
39  * ordinate_type getU() const;
40  * void setU(ordinate_type u);
41  * </pre>
42  * to get and set the derivatives, respectively.
43  *
44  * Note that -by default- the compilation is for a monotonous interpolation.
45  */
46  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
48  public JCollection_t<JElement_t, JDistance_t>,
49  public virtual JFunctional<>
50  {
51  public:
52 
53  typedef JCollection_t<JElement_t, JDistance_t> collection_type;
54 
55  typedef typename collection_type::abscissa_type abscissa_type;
56  typedef typename collection_type::ordinate_type ordinate_type;
57  typedef typename collection_type::value_type value_type;
58  typedef typename collection_type::distance_type distance_type;
59 
60  typedef typename collection_type::const_iterator const_iterator;
61  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
62  typedef typename collection_type::iterator iterator;
63  typedef typename collection_type::reverse_iterator reverse_iterator;
64 
66 
67 
68  /**
69  * Determination of derivatives.
70  *
71  * \param monotone monotone
72  */
73  void compile(const bool monotone)
74  {
75  using namespace std;
76 
77  if (this->size() >= 2u) {
78 
79  {
80  iterator j = this->begin(), i = j++;
81 
82  i->setU((j->getY() - i->getY()) / this->getDistance(i->getX(), j->getX()));
83  }
84 
85  {
86  reverse_iterator j = this->rbegin(), i = j++;
87 
88  i->setU((j->getY() - i->getY()) / this->getDistance(i->getX(), j->getX()));
89  }
90 
91  for (iterator k = this->begin(), i = k++, j = k++; k != this->end(); ++i, ++j, ++k) {
92  j->setU(0.5 * ((j->getY() - i->getY()) / this->getDistance(i->getX(), j->getX()) +
93  (k->getY() - j->getY()) / this->getDistance(j->getX(), k->getX())));
94  }
95 
96  if (monotone) {
97 
98  for (iterator j = this->begin(), i = j++; j != this->end(); ++i, ++j) {
99  if (i->getY() == j->getY()) {
100  j->setU(JMATH::zero);
101  }
102  }
103 
104  for (iterator j = this->begin(), i = j++; j != this->end(); ++i, ++j) {
105 
106  const ordinate_type u = (j->getY() - i->getY()) / this->getDistance(i->getX(), j->getX());
107  const ordinate_type w = (i->getU()*i->getU() + j->getU()*j->getU());
108 
109  if (w > 9.0*u*u) {
110 
111  const ordinate_type v = 3.0*u/sqrt(w);
112 
113  i->setU(v*i->getU());
114  j->setU(v*j->getU());
115  }
116  }
117  }
118  }
119  }
120 
121 
122  protected:
123 
124  static abscissa_type h00 (abscissa_type t) { return (1.0 + 2*t) * (1.0 - t) * (1.0 - t); }
125  static abscissa_type h10 (abscissa_type t) { return t * (1.0 - t) * (1.0 - t); }
126  static abscissa_type h01 (abscissa_type t) { return t * t * (3.0 - 2*t); }
127  static abscissa_type h11 (abscissa_type t) { return t * t * (t - 1.0); }
128 
129  static abscissa_type h00p(abscissa_type t) { return 6 * t * (t - 1.0); }
130  static abscissa_type h10p(abscissa_type t) { return t * (3*t - 4.0) + 1.0; }
131  static abscissa_type h01p(abscissa_type t) { return 6 * t * (1.0 -t); }
132  static abscissa_type h11p(abscissa_type t) { return t * (3*t - 2.0); }
133 
134  static abscissa_type H00 (abscissa_type t) { return t * (t * t * (0.5*t - 1.0) + 1.0); }
135  static abscissa_type H10 (abscissa_type t) { return t * t * (t * (0.25*t - 2.0/3.0) + 0.5); }
136  static abscissa_type H01 (abscissa_type t) { return t * t * t * (1.0 - 0.5*t); }
137  static abscissa_type H11 (abscissa_type t) { return t * t * t * (0.25*t - 1.0/3.0); }
138 
139 
140  /**
141  * Default constructor.
142  */
144  {}
145 
146 
147  /**
148  * Determination of derivatives.
149  */
150  virtual void do_compile() override
151  {
152  compile(true);
153  }
154  };
155 
156 
157  /**
158  * Template definition for functional collection with spline interpolation.
159  */
160  template<class JElement_t,
161  template<class, class> class JCollection_t,
162  class JResult_t,
163  class JDistance_t>
165 
166 
167  /**
168  * Template specialisation for functional collection with spline interpolation.
169  */
170  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
171  class JHermiteSplineFunction<JElement_t,
172  JCollection_t,
173  typename JResultType<typename JElement_t::ordinate_type>::result_type,
174  JDistance_t> :
175  public JHermiteSplineCollection<JElement_t, JCollection_t, JDistance_t>,
176  public virtual JFunction<typename JElement_t::abscissa_type,
177  typename JResultType<typename JElement_t::ordinate_type>::result_type>
178  {
179  public:
180 
182 
183  typedef typename collection_type::abscissa_type abscissa_type;
184  typedef typename collection_type::ordinate_type ordinate_type;
185  typedef typename collection_type::value_type value_type;
186  typedef typename collection_type::distance_type distance_type;
187 
188  typedef typename collection_type::const_iterator const_iterator;
189  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
190  typedef typename collection_type::iterator iterator;
191  typedef typename collection_type::reverse_iterator reverse_iterator;
192 
195 
199 
200 
201  /**
202  * Default constructor.
203  */
205  {}
206 
207 
208  /**
209  * Recursive interpolation method implementation.
210  *
211  * \param pX pointer to abscissa values
212  * \return function value
213  */
214  virtual result_type evaluate(const argument_type* pX) const override
215  {
216  const argument_type x = *pX;
217 
218  if (this->size() <= 1u) {
219 
220  std::ostringstream os;
221 
222  os << __FILE__ << ':' << __LINE__ << " not enough data " << STREAM("?") << x;
223 
224  return this->getExceptionHandler().action(JFunctionalException(os.str()));
225  }
226 
227  const_iterator p = this->lower_bound(x);
228 
229  if ((p == this->begin() && this->getDistance(x, (p++)->getX()) > distance_type::precision) ||
230  (p == this->end() && this->getDistance((--p)->getX(), x) > distance_type::precision)) {
231 
232  std::ostringstream os;
233 
234  os << __FILE__ << ':' << __LINE__ << " abscissa out of range "
235  << STREAM("?") << x << " <> "
236  << STREAM("?") << this->begin() ->getX() << ' '
237  << STREAM("?") << this->rbegin()->getX();
238 
239  return this->getExceptionHandler().action(JValueOutOfRange(os.str()));
240  }
241 
242  const_iterator q = p--;
243 
244  const double dx = this->getDistance(p->getX(), q->getX());
245  const double t = this->getDistance(p->getX(), x) / dx;
246 
247  return h00(t)*p->getY() + h10(t)*p->getU()*dx + h01(t)*q->getY() + h11(t)*q->getU()*dx;
248  }
249 
250  protected:
251 
252  using collection_type::h00;
253  using collection_type::h10;
254  using collection_type::h01;
255  using collection_type::h11;
256  };
257 
258 
259  /**
260  * Template specialisation for spline interpolation method with returning JResultDerivative data structure.
261  */
262  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
263  class JHermiteSplineFunction<JElement_t,
264  JCollection_t,
265  JResultDerivative<typename JResultType<typename JElement_t::ordinate_type>::result_type>,
266  JDistance_t> :
267  public JHermiteSplineCollection<JElement_t, JCollection_t, JDistance_t>,
268  public virtual JFunction<typename JElement_t::abscissa_type,
269  JResultDerivative<typename JResultType<typename JElement_t::ordinate_type>::result_type> >
270  {
271  public:
272 
274 
275  typedef typename collection_type::abscissa_type abscissa_type;
276  typedef typename collection_type::ordinate_type ordinate_type;
277  typedef typename collection_type::value_type value_type;
278  typedef typename collection_type::distance_type distance_type;
279 
280  typedef typename collection_type::const_iterator const_iterator;
281  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
282  typedef typename collection_type::iterator iterator;
283  typedef typename collection_type::reverse_iterator reverse_iterator;
284 
287 
291 
292 
293  /**
294  * Default constructor.
295  */
297  {}
298 
299 
300  /**
301  * Recursive interpolation method implementation.
302  *
303  * \param pX pointer to abscissa values
304  * \return function value
305  */
306  virtual result_type evaluate(const argument_type* pX) const override
307  {
308  const argument_type x = *pX;
309 
310  if (this->size() <= 1u) {
311 
312  std::ostringstream os;
313 
314  os << __FILE__ << ':' << __LINE__ << " not enough data " << STREAM("?") << x;
315 
316  return this->getExceptionHandler().action(JFunctionalException(os.str()));
317  }
318 
319  const_iterator p = this->lower_bound(x);
320 
321 
322  if ((p == this->begin() && this->getDistance(x, (p++)->getX()) > distance_type::precision) ||
323  (p == this->end() && this->getDistance((--p)->getX(), x) > distance_type::precision)) {
324 
325  std::ostringstream os;
326 
327  os << __FILE__ << ':' << __LINE__ << " abscissa out of range "
328  << STREAM("?") << x << " <> "
329  << STREAM("?") << this->begin() ->getX() << ' '
330  << STREAM("?") << this->rbegin()->getX();
331 
332  return this->getExceptionHandler().action(JValueOutOfRange(os.str()));
333  }
334 
335  const_iterator q = p--;
336 
337  const double dx = this->getDistance(p->getX(), q->getX());
338  const double t = this->getDistance(p->getX(), x) / dx;
339 
340  result.f = h00 (t)*p->getY() + h10 (t)*p->getU()*dx + h01 (t)*q->getY() + h11 (t)*q->getU()*dx;
341  result.fp = h00p(t)*p->getY()/dx + h10p(t)*p->getU() + h01p(t)*q->getY()/dx + h11p(t)*q->getU();
342 
343  return result;
344  }
345 
346 
347  protected:
348 
349  using collection_type::h00;
350  using collection_type::h10;
351  using collection_type::h01;
352  using collection_type::h11;
353 
354  using collection_type::h00p;
355  using collection_type::h10p;
356  using collection_type::h01p;
357  using collection_type::h11p;
358 
359  private:
361  };
362 
363 
364  /**
365  * Template specialisation for spline interpolation method with returning JResultPDF data structure.
366  *
367  * Note that the data structure of the elements in the collection should have the additional methods:
368  * <pre>
369  * ordinate_type getIntegral() const;
370  * void setIntegral(ordinate_type v);
371  * </pre>
372  * to get and set the integral values, respectively.
373  */
374  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
375  class JHermiteSplineFunction<JElement_t,
376  JCollection_t,
377  JResultPDF<typename JResultType<typename JElement_t::ordinate_type>::result_type>,
378  JDistance_t> :
379  public JHermiteSplineCollection<JElement_t, JCollection_t, JDistance_t>,
380  public virtual JFunction<typename JElement_t::abscissa_type,
381  JResultPDF<typename JResultType<typename JElement_t::ordinate_type>::result_type> >
382  {
383  public:
384 
386 
387  typedef typename collection_type::abscissa_type abscissa_type;
388  typedef typename collection_type::ordinate_type ordinate_type;
389  typedef typename collection_type::value_type value_type;
390  typedef typename collection_type::distance_type distance_type;
391 
392  typedef typename collection_type::const_iterator const_iterator;
393  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
394  typedef typename collection_type::iterator iterator;
395  typedef typename collection_type::reverse_iterator reverse_iterator;
396 
399 
403 
404 
405  /**
406  * Default constructor.
407  */
409  {}
410 
411 
412  /**
413  * Recursive interpolation method implementation.
414  *
415  * \param pX pointer to abscissa values
416  * \return function value
417  */
418  virtual result_type evaluate(const argument_type* pX) const override
419  {
420  const argument_type x = *pX;
421 
422  if (this->size() <= 1u) {
423 
424  std::ostringstream os;
425 
426  os << __FILE__ << ':' << __LINE__ << " not enough data " << STREAM("?") << x;
427 
428  return this->getExceptionHandler().action(JFunctionalException(os.str()));
429  }
430 
431  const_iterator p = this->lower_bound(x);
432 
433  if (p == this->begin() && this->getDistance(x, (p++)->getX()) > distance_type::precision) {
434 
435  try {
436 
437  std::ostringstream os;
438 
439  os << __FILE__ << ':' << __LINE__ << " abscissa out of range " << STREAM("?") << x << " < " << STREAM("?") << this->begin() ->getX();
440 
441  result = this->getExceptionHandler().action(JValueOutOfRange(os.str()));
442 
443  // overwrite integral values
444 
445  result.v = 0;
446  result.V = this->rbegin()->getIntegral();
447 
448  } catch(const JValueOutOfRange& exception) {
449  throw exception;
450  }
451 
452  return result;
453 
454  } else if (p == this->end() && this->getDistance((--p)->getX(), x) > distance_type::precision) {
455 
456  try {
457 
458  std::ostringstream os;
459 
460  os << __FILE__ << ':' << __LINE__ << " abscissa out of range " << STREAM("?") << x << " > " << STREAM("?") << this->rbegin()->getX();
461 
462  result = this->getExceptionHandler().action(JValueOutOfRange(os.str()));
463 
464  // overwrite integral values
465 
466  result.v = this->rbegin()->getIntegral();
467  result.V = this->rbegin()->getIntegral();
468 
469  } catch(const JValueOutOfRange& exception) {
470  throw exception;
471  }
472 
473  return result;
474  }
475 
476  const_iterator q = p--;
477 
478  const double dx = this->getDistance(p->getX(), q->getX());
479  const double t = this->getDistance(p->getX(), x) / dx;
480 
481  result.f = h00 (t)*p->getY() + h10 (t)*p->getU()*dx + h01 (t)*q->getY() + h11 (t)*q->getU()*dx;
482  result.fp = h00p(t)*p->getY()/dx + h10p(t)*p->getU() + h01p(t)*q->getY()/dx + h11p(t)*q->getU();
483  result.v = (p->getIntegral() +
484  (H00 (t)*p->getY() + H10 (t)*p->getU()*dx + H01 (t)*q->getY() + H11 (t)*q->getU()*dx)*dx);
485  result.V = this->rbegin()->getIntegral();
486 
487  return result;
488  }
489 
490 
491  protected:
492 
493  using collection_type::h00;
494  using collection_type::h10;
495  using collection_type::h01;
496  using collection_type::h11;
497 
498  using collection_type::h00p;
499  using collection_type::h10p;
500  using collection_type::h01p;
501  using collection_type::h11p;
502 
503  using collection_type::H00;
504  using collection_type::H10;
505  using collection_type::H01;
506  using collection_type::H11;
507 
508  /**
509  * Determination of derivatives.
510  */
511  virtual void do_compile() override
512  {
513  if (!this->empty()) {
514 
516 
517  this->begin()->setIntegral(JMATH::zero);
518 
519  for (iterator j = this->begin(), i = j++; j != this->end(); ++i, ++j) {
520 
521  const double dx = this->getDistance(i->getX(), j->getX());
522  const ordinate_type y = i->getY() + j->getY();
523  const ordinate_type z = i->getU() - j->getU();
524 
525  const ordinate_type v = dx * 0.50 * y;
526  const ordinate_type w = dx * 1.00 * z*dx/12.0;
527 
528  j->setIntegral(i->getIntegral() + v + w);
529  }
530  }
531  }
532 
533  private:
535  };
536 
537 
538  /**
539  * Template class for spline interpolation in 1D
540  *
541  * This class implements the JFunction1D interface.
542  */
543  template<class JElement_t,
544  template<class, class> class JCollection_t,
545  class JResult_t = typename JElement_t::ordinate_type,
548  public JHermiteSplineFunction<JElement_t, JCollection_t, JResult_t, JDistance_t>,
549  public virtual JFunction1D<typename JElement_t::abscissa_type, JResult_t>
550  {
551  public:
552 
554 
559 
564 
566 
570 
571 
572  /**
573  * Default contructor.
574  */
576  {}
577  };
578 
579 
580  /**
581  * Functional map with spline interpolation.
582  */
583  template<class JKey_t,
584  class JValue_t,
585  template<class, class, class> class JMap_t,
586  class JResult_t,
587  class JDistance_t = JDistance<JKey_t> >
589  public JMap_t<JKey_t, JValue_t, JDistance_t>,
590  public JFunction<JKey_t, JResult_t>
591  {
592  public:
593 
594  typedef JMap_t<JKey_t, JValue_t, JDistance_t> collection_type;
596 
597  typedef typename collection_type::abscissa_type abscissa_type;
598  typedef typename collection_type::ordinate_type ordinate_type;
599  typedef typename collection_type::value_type value_type;
600  typedef typename collection_type::distance_type distance_type;
601 
602  typedef typename collection_type::const_iterator const_iterator;
603  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
604  typedef typename collection_type::iterator iterator;
605  typedef typename collection_type::reverse_iterator reverse_iterator;
606 
609  typedef typename function_type::JExceptionHandler exceptionhandler_type;
610 
615 
616 
617  /**
618  * Default constructor.
619  */
621  {}
622 
623 
624  /**
625  * Recursive interpolation method implementation.
626  *
627  * \param pX pointer to abscissa values
628  * \return function value
629  */
630  virtual result_type evaluate(const argument_type* pX) const override
631  {
632  const argument_type x = *pX;
633 
634  ++pX; // next argument value
635 
636  const_iterator p = this->begin();
637 
638  for (typename JHermiteSplineFunction1D_t::iterator q = buffer.begin(); q != buffer.end(); ++q, ++p) {
639  q->getY() = JFunction<argument_type, data_type>::getValue(p->getY(), pX);
640  }
641 
642  buffer.compile();
643 
644  return buffer(x);
645  }
646 
647 
648  private:
649  /**
650  * Function compilation.
651  */
652  virtual void do_compile() override
653  {
654  buffer.clear();
655 
656  for (iterator i = this->begin(); i != this->end(); ++i) {
657  buffer.put(i->getX(), data_type());
658  }
659  }
660 
661 
663  };
664 
665 
666  /**
667  * Conversion of data points to integral values.
668  *
669  * The integration includes the use of 2nd derivatives of the data points of the input spline interpolating function.
670  *
671  * \param input collection
672  * \param output mappable collection
673  * \return integral
674  */
675  template<class JElement_t,
676  template<class, class> class JCollection_t,
677  class JResult_t,
678  class JDistance_t>
679  inline typename JElement_t::ordinate_type
681  typename JMappable<JElement_t>::map_type& output)
682  {
683  typedef typename JElement_t::ordinate_type ordinate_type;
685 
686  ordinate_type V(JMATH::zero);
687 
688  if (input.getSize() > 1) {
689 
690  output.put(input.begin()->getX(), V);
691 
692  for (const_iterator j = input.begin(), i = j++; j != input.end(); ++i, ++j) {
693 
694  const double dx = input.getDistance(i->getX(), j->getX());
695  const ordinate_type y = i->getY() + j->getY();
696  const ordinate_type z = i->getU() - j->getU();
697 
698  const ordinate_type v = dx * 0.50 * y;
699  const ordinate_type w = dx * 1.00 * z*dx/12.0;
700 
701  V += v + w;
702 
703  output.put(j->getX(), V);
704  }
705  }
706 
707  return V;
708  }
709 
710 
711  /**
712  * Conversion of data points to integral values.
713  *
714  * The integration directly uses the integral values of the input spline interpolating function.
715  *
716  * \param input collection
717  * \param output mappable collection
718  * \return integral
719  */
720  template<class JElement_t,
721  template<class, class> class JCollection_t,
722  class JDistance_t>
723  inline typename JElement_t::ordinate_type
724  integrate(const JHermiteSplineFunction1D<JElement_t, JCollection_t, JResultPDF<typename JElement_t::ordinate_type>, JDistance_t>& input,
725  typename JMappable<JElement_t>::map_type& output)
726  {
727  typedef typename JElement_t::ordinate_type ordinate_type;
730 
731  if (input.getSize() > 1) {
732 
733  for (const_iterator i = input.begin(); i != input.end(); ++i) {
734  output.put(i->getX(), i->getIntegral());
735  }
736 
737  return input.rbegin()->getIntegral();
738  }
739 
740  return JMATH::zero;
741  }
742 }
743 
744 #endif
static abscissa_type h01p(abscissa_type t)
function_type::argument_type argument_type
JHermiteSplineFunction1D_t buffer
collection_type::value_type value_type
data_type w[N+1][M+1]
Definition: JPolint.hh:823
collection_type::value_type value_type
Exceptions.
then fatal No hydrophone data file $HYDROPHONE_TXT fi sort gr k
collection_type::const_iterator const_iterator
void compile(const bool monotone)
Determination of derivatives.
collection_type::ordinate_type ordinate_type
std::vector< event_type > data_type
Definition: JPerth.cc:78
static abscissa_type h00(abscissa_type t)
Exception for a functional operation.
Definition: JException.hh:142
collection_type::const_reverse_iterator const_reverse_iterator
JHermiteSplineFunction1D()
Default contructor.
Template class for distance evaluation.
Definition: JDistance.hh:24
This include file containes various data structures that can be used as specific return types for the...
functional_type::argument_type argument_type
Definition: JFunctional.hh:323
JHermiteSplineMap()
Default constructor.
virtual void do_compile() override
Function compilation.
JCollection_t< JElement_t, JDistance_t > collection_type
collection_type::value_type value_type
static abscissa_type H01(abscissa_type t)
Template interface definition for associative collection of elements.
virtual void do_compile() override
Function compilation.
Definition: JPolint.hh:816
static const JZero zero
Function object to assign zero value.
Definition: JZero.hh:105
static abscissa_type h11(abscissa_type t)
V(JDAQEvent-JTriggerReprocessor)*1.0/(JDAQEvent+1.0e-10)
Definition of zero value for any class.
Template definition of function object interface in one dimension.
Definition: JFunctional.hh:317
collection_type::const_iterator const_iterator
double getDistance(const JFirst_t &first, const JSecond_t &second)
Get distance between objects.
JFunction1D< abscissa_type, JResult_t > function_type
JResultType< ordinate_type >::result_type data_type
Exception for missing value.
Definition: JException.hh:214
JMap_t< JKey_t, JValue_t, JDistance_t > collection_type
collection_type::abscissa_type abscissa_type
static result_type getValue(const JFunctional &function, const argument_type *pX)
Recursive function value evaluation.
Definition: JFunctional.hh:107
collection_type::reverse_iterator reverse_iterator
static abscissa_type h00p(abscissa_type t)
virtual void do_compile() override
Determination of derivatives.
static abscissa_type h11p(abscissa_type t)
collection_type::distance_type distance_type
virtual result_type evaluate(const argument_type *pX) const override
Recursive interpolation method implementation.
virtual result_type evaluate(const argument_type *pX) const override
Recursive interpolation method implementation.
collection_type::distance_type distance_type
collection_type::abscissa_type abscissa_type
static abscissa_type h10p(abscissa_type t)
function_type::result_type result_type
Data structure for result including value, first derivative and integrals of function.
Definition: JResult.hh:337
static abscissa_type H10(abscissa_type t)
Auxiliary data structure for handling std::ostream.
Template definition of function object interface.
Definition: JFunctional.hh:32
functional_type::result_type result_type
Definition: JFunctional.hh:324
Auxiliary class to evaluate result type.
Definition: JFunctional.hh:377
static abscissa_type h01(abscissa_type t)
collection_type::const_reverse_iterator const_reverse_iterator
void put(typename JClass< key_type >::argument_type key, typename JClass< mapped_type >::argument_type value)
Put pair-wise element (key,value) into collection.
collection_type::iterator iterator
Template class to define the corresponding JCollection for a given template JMap. ...
collection_type::abscissa_type abscissa_type
collection_type::reverse_iterator reverse_iterator
Template base class spline interpolations.
Data structure for result including value and first derivative of function.
Definition: JResult.hh:43
function_type::argument_type argument_type
virtual result_type evaluate(const argument_type *pX) const override
Recursive interpolation method implementation.
Exception handler for functional object.
Definition: JFunctional.hh:131
Template class for spline interpolation in 1D.
Template definition of function object interface in multidimensions.
Definition: JFunctional.hh:303
static abscissa_type H00(abscissa_type t)
function_type::JExceptionHandler exceptionhandler_type
collection_type::ordinate_type ordinate_type
functional_type::result_type result_type
Definition: JFunctional.hh:308
Functional map with spline interpolation.
collection_type::iterator iterator
const JExceptionHandler & getExceptionHandler() const
Get exception handler.
Definition: JFunctional.hh:277
Exception for division by zero.
Definition: JException.hh:286
JHermiteSplineCollection()
Default constructor.
collection_type::iterator iterator
collection_type::reverse_iterator reverse_iterator
JFunction< JKey_t, JResult_t > function_type
JHermiteSplineFunction1D< JSplineElement2D< argument_type, data_type >, JMapCollection< JMap_t >::template collection_type, result_type > JHermiteSplineFunction1D_t
int j
Definition: JPolint.hh:748
JHermiteSplineCollection< JElement_t, JCollection_t, JDistance_t > collection_type
collection_type::distance_type distance_type
Exception for accessing a value in a collection that is outside of its range.
Definition: JException.hh:178
Template definition for functional collection with spline interpolation.
functional_type::argument_type argument_type
Definition: JFunctional.hh:307
static abscissa_type H11(abscissa_type t)
collection_type::ordinate_type ordinate_type
data_type v[N+1][M+1]
Definition: JPolint.hh:822
double u[N+1]
Definition: JPolint.hh:821
static abscissa_type h10(abscissa_type t)
function_type::result_type result_type
collection_type::const_reverse_iterator const_reverse_iterator
JElement_t::ordinate_type integrate(const JCollection< JElement_t, JDistance_t > &input, typename JMappable< JElement_t >::map_type &output)
Conversion of data points to integral values.
Definition: JCollection.hh:812
collection_type::const_iterator const_iterator
function_type::JExceptionHandler exceptionhandler_type