Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Sacado_PCE_OrthogPoly.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Stokhos Package
4 //
5 // Copyright 2009 NTESS and the Stokhos contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef SACADO_PCE_ORTHOGPOLY_HPP
11 #define SACADO_PCE_ORTHOGPOLY_HPP
12 
13 #include "Stokhos_ConfigDefs.h"
14 
15 #ifdef HAVE_STOKHOS_SACADO
16 
17 #include "Teuchos_RCP.hpp"
18 
19 #include "Sacado_Traits.hpp"
20 #include "Sacado_Handle.hpp"
21 #include "Sacado_mpl_apply.hpp"
22 
26 
27 #include <cmath>
28 #include <algorithm> // for std::min and std::max
29 #include <ostream> // for std::ostream
30 
31 namespace Sacado {
32 
34  namespace PCE {
35 
37 
41  template <typename T, typename Storage >
42  class OrthogPoly {
43  public:
44 
46  typedef T value_type;
47 
49  typedef typename ScalarType<T>::type scalar_type;
50 
52  typedef int ordinal_type;
53 
55  typedef Storage storage_type;
56 
59 
62 
65 
66  typedef typename approx_type::pointer pointer;
67  typedef typename approx_type::const_pointer const_pointer;
68  typedef typename approx_type::reference reference;
69  typedef typename approx_type::const_reference const_reference;
70 
72  template <typename S>
73  struct apply {
74  typedef typename Sacado::mpl::apply<Storage,ordinal_type,S>::type storage_type;
75  typedef OrthogPoly<S,storage_type> type;
76  };
77 
79 
82  OrthogPoly();
83 
85 
88  OrthogPoly(const value_type& x);
89 
91 
94  OrthogPoly(const Teuchos::RCP<expansion_type>& expansion);
95 
97 
100  OrthogPoly(const Teuchos::RCP<expansion_type>& expansion,
101  ordinal_type sz);
102 
104  OrthogPoly(const OrthogPoly& x);
105 
107  ~OrthogPoly();
108 
110  void init(const T& v) { th->init(v); }
111 
113  void init(const T* v) { th->init(v); }
114 
116  template <typename S>
117  void init(const OrthogPoly<T,S>& v) { th->init(v.getOrthogPolyApprox()); }
118 
120  void load(T* v) { th->load(v); }
121 
123  template <typename S>
124  void load(OrthogPoly<T,S>& v) { th->load(v.getOrthogPolyApprox()); }
125 
127 
130  void reset(const Teuchos::RCP<expansion_type>& expansion);
131 
133 
136  void reset(const Teuchos::RCP<expansion_type>& expansion,
137  ordinal_type sz);
138 
140 
149  void copyForWrite() { th.makeOwnCopy(); }
150 
152  value_type evaluate(const Teuchos::Array<value_type>& point) const;
153 
155  value_type evaluate(const Teuchos::Array<value_type>& point,
156  const Teuchos::Array<value_type>& bvals) const;
157 
159  value_type mean() const {return th->mean(); }
160 
162  value_type standard_deviation() const { return th->standard_deviation(); }
163 
165  value_type two_norm() const { return th->two_norm(); }
166 
168  value_type two_norm_squared() const { return th->two_norm_squared(); }
169 
171  value_type inner_product(const OrthogPoly& b) const {
172  return th->inner_product(b.getOrthogPolyApprox()); }
173 
175  std::ostream& print(std::ostream& os) const { return th->print(os); }
176 
178  bool isEqualTo(const OrthogPoly& x) const;
179 
184 
186  OrthogPoly<T,Storage>& operator=(const value_type& val);
187 
189  OrthogPoly<T,Storage>& operator=(const OrthogPoly<T,Storage>& x);
190 
192 
197 
199  Teuchos::RCP<const basis_type> basis() const { return th->basis(); }
200 
202  Teuchos::RCP<expansion_type> expansion() const { return expansion_; }
203 
205 
210 
212  const_reference val() const { return (*th)[0]; }
213 
215  reference val() { return (*th)[0]; }
216 
218 
223 
225  ordinal_type size() const { return th->size();}
226 
228  bool hasFastAccess(ordinal_type sz) const { return th->size()>=sz;}
229 
231  const_pointer coeff() const { return th->coeff();}
232 
234  pointer coeff() { return th->coeff();}
235 
237  value_type coeff(ordinal_type i) const {
238  value_type tmp= i<th->size() ? (*th)[i]:value_type(0.); return tmp;}
239 
241  reference fastAccessCoeff(ordinal_type i) { return (*th)[i];}
242 
244  value_type fastAccessCoeff(ordinal_type i) const { return (*th)[i];}
245 
247  reference term(ordinal_type dimension, ordinal_type order) {
248  return th->term(dimension, order); }
249 
251  const_reference term(ordinal_type dimension, ordinal_type order) const {
252  return th->term(dimension, order); }
253 
255  Teuchos::Array<ordinal_type> order(ordinal_type term) const {
256  return th->order(term); }
257 
259 
264 
266  OrthogPoly<T,Storage> operator + () const;
267 
269  OrthogPoly<T,Storage> operator - () const;
270 
272  OrthogPoly<T,Storage>& operator += (const value_type& x);
273 
275  OrthogPoly<T,Storage>& operator -= (const value_type& x);
276 
278  OrthogPoly<T,Storage>& operator *= (const value_type& x);
279 
281  OrthogPoly<T,Storage>& operator /= (const value_type& x);
282 
284  OrthogPoly<T,Storage>& operator += (const OrthogPoly<T,Storage>& x);
285 
287  OrthogPoly<T,Storage>& operator -= (const OrthogPoly<T,Storage>& x);
288 
290  OrthogPoly<T,Storage>& operator *= (const OrthogPoly<T,Storage>& x);
291 
293  OrthogPoly<T,Storage>& operator /= (const OrthogPoly<T,Storage>& x);
294 
296 
298  const approx_type& getOrthogPolyApprox() const { return *th; }
299 
301  approx_type& getOrthogPolyApprox() { return *th; }
302 
303  protected:
304 
306  Teuchos::RCP<expansion_type> expansion_;
307 
309  Teuchos::RCP<expansion_type> const_expansion_;
310 
311  Sacado::Handle< Stokhos::OrthogPolyApprox<int,value_type,Storage> > th;
312 
313  }; // class Hermite
314 
315  // Operations
316  template <typename T, typename Storage> OrthogPoly<T,Storage>
317  operator+(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
318 
319  template <typename T, typename Storage> OrthogPoly<T,Storage>
320  operator+(const typename OrthogPoly<T,Storage>::value_type& a,
321  const OrthogPoly<T,Storage>& b);
322 
323  template <typename T, typename Storage> OrthogPoly<T,Storage>
324  operator+(const OrthogPoly<T,Storage>& a,
325  const typename OrthogPoly<T,Storage>::value_type& b);
326 
327  template <typename T, typename Storage> OrthogPoly<T,Storage>
328  operator-(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
329 
330  template <typename T, typename Storage> OrthogPoly<T,Storage>
331  operator-(const typename OrthogPoly<T,Storage>::value_type& a,
332  const OrthogPoly<T,Storage>& b);
333 
334  template <typename T, typename Storage> OrthogPoly<T,Storage>
335  operator-(const OrthogPoly<T,Storage>& a,
336  const typename OrthogPoly<T,Storage>::value_type& b);
337 
338  template <typename T, typename Storage> OrthogPoly<T,Storage>
339  operator*(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
340 
341  template <typename T, typename Storage> OrthogPoly<T,Storage>
342  operator*(const typename OrthogPoly<T,Storage>::value_type& a,
343  const OrthogPoly<T,Storage>& b);
344 
345  template <typename T, typename Storage> OrthogPoly<T,Storage>
346  operator*(const OrthogPoly<T,Storage>& a,
347  const typename OrthogPoly<T,Storage>::value_type& b);
348 
349  template <typename T, typename Storage> OrthogPoly<T,Storage>
350  operator/(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
351 
352  template <typename T, typename Storage> OrthogPoly<T,Storage>
353  operator/(const typename OrthogPoly<T,Storage>::value_type& a,
354  const OrthogPoly<T,Storage>& b);
355 
356  template <typename T, typename Storage> OrthogPoly<T,Storage>
357  operator/(const OrthogPoly<T,Storage>& a,
358  const typename OrthogPoly<T,Storage>::value_type& b);
359 
360  template <typename T, typename Storage> OrthogPoly<T,Storage>
361  exp(const OrthogPoly<T,Storage>& a);
362 
363  template <typename T, typename Storage> OrthogPoly<T,Storage>
364  log(const OrthogPoly<T,Storage>& a);
365 
366  template <typename T, typename Storage> void
367  log(OrthogPoly<T,Storage>& c, const OrthogPoly<T,Storage>& a);
368 
369  template <typename T, typename Storage> OrthogPoly<T,Storage>
370  log10(const OrthogPoly<T,Storage>& a);
371 
372  template <typename T, typename Storage> OrthogPoly<T,Storage>
373  sqrt(const OrthogPoly<T,Storage>& a);
374 
375  template <typename T, typename Storage> OrthogPoly<T,Storage>
376  cbrt(const OrthogPoly<T,Storage>& a);
377 
378  template <typename T, typename Storage> OrthogPoly<T,Storage>
379  pow(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
380 
381  template <typename T, typename Storage> OrthogPoly<T,Storage>
382  pow(const T& a,
383  const OrthogPoly<T,Storage>& b);
384 
385  template <typename T, typename Storage> OrthogPoly<T,Storage>
386  pow(const OrthogPoly<T,Storage>& a,
387  const T& b);
388 
389  template <typename T, typename Storage> OrthogPoly<T,Storage>
390  cos(const OrthogPoly<T,Storage>& a);
391 
392  template <typename T, typename Storage> OrthogPoly<T,Storage>
393  sin(const OrthogPoly<T,Storage>& a);
394 
395  template <typename T, typename Storage> OrthogPoly<T,Storage>
396  tan(const OrthogPoly<T,Storage>& a);
397 
398  template <typename T, typename Storage> OrthogPoly<T,Storage>
399  cosh(const OrthogPoly<T,Storage>& a);
400 
401  template <typename T, typename Storage> OrthogPoly<T,Storage>
402  sinh(const OrthogPoly<T,Storage>& a);
403 
404  template <typename T, typename Storage> OrthogPoly<T,Storage>
405  tanh(const OrthogPoly<T,Storage>& a);
406 
407  template <typename T, typename Storage> OrthogPoly<T,Storage>
408  acos(const OrthogPoly<T,Storage>& a);
409 
410  template <typename T, typename Storage> OrthogPoly<T,Storage>
411  asin(const OrthogPoly<T,Storage>& a);
412 
413  template <typename T, typename Storage> OrthogPoly<T,Storage>
414  atan(const OrthogPoly<T,Storage>& a);
415 
416  template <typename T, typename Storage> OrthogPoly<T,Storage>
417  atan2(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
418 
419  template <typename T, typename Storage> OrthogPoly<T,Storage>
420  atan2(const typename OrthogPoly<T,Storage>::value_type& a,
421  const OrthogPoly<T,Storage>& b);
422 
423  template <typename T, typename Storage> OrthogPoly<T,Storage>
424  atan2(const OrthogPoly<T,Storage>& a,
425  const typename OrthogPoly<T,Storage>::value_type& b);
426 
427  template <typename T, typename Storage> OrthogPoly<T,Storage>
428  acosh(const OrthogPoly<T,Storage>& a);
429 
430  template <typename T, typename Storage> OrthogPoly<T,Storage>
431  asinh(const OrthogPoly<T,Storage>& a);
432 
433  template <typename T, typename Storage> OrthogPoly<T,Storage>
434  atanh(const OrthogPoly<T,Storage>& a);
435 
436  template <typename T, typename Storage> OrthogPoly<T,Storage>
437  abs(const OrthogPoly<T,Storage>& a);
438 
439  template <typename T, typename Storage> OrthogPoly<T,Storage>
440  fabs(const OrthogPoly<T,Storage>& a);
441 
442  template <typename T, typename Storage> OrthogPoly<T,Storage>
443  max(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
444 
445  template <typename T, typename Storage> OrthogPoly<T,Storage>
446  max(const typename OrthogPoly<T,Storage>::value_type& a,
447  const OrthogPoly<T,Storage>& b);
448 
449  template <typename T, typename Storage> OrthogPoly<T,Storage>
450  max(const OrthogPoly<T,Storage>& a,
451  const typename OrthogPoly<T,Storage>::value_type& b);
452 
453  template <typename T, typename Storage> OrthogPoly<T,Storage>
454  min(const OrthogPoly<T,Storage>& a, const OrthogPoly<T,Storage>& b);
455 
456  template <typename T, typename Storage> OrthogPoly<T,Storage>
457  min(const typename OrthogPoly<T,Storage>::value_type& a,
458  const OrthogPoly<T,Storage>& b);
459 
460  template <typename T, typename Storage> OrthogPoly<T,Storage>
461  min(const OrthogPoly<T,Storage>& a,
462  const typename OrthogPoly<T,Storage>::value_type& b);
463 
464  template <typename T, typename Storage> bool
465  operator==(const OrthogPoly<T,Storage>& a,
466  const OrthogPoly<T,Storage>& b);
467 
468  template <typename T, typename Storage> bool
470  const OrthogPoly<T,Storage>& b);
471 
472  template <typename T, typename Storage> bool
473  operator==(const OrthogPoly<T,Storage>& a,
474  const typename OrthogPoly<T,Storage>::value_type& b);
475 
476  template <typename T, typename Storage> bool
477  operator!=(const OrthogPoly<T,Storage>& a,
478  const OrthogPoly<T,Storage>& b);
479 
480  template <typename T, typename Storage> bool
482  const OrthogPoly<T,Storage>& b);
483 
484  template <typename T, typename Storage> bool
485  operator!=(const OrthogPoly<T,Storage>& a,
486  const typename OrthogPoly<T,Storage>::value_type& b);
487 
488  template <typename T, typename Storage> bool
489  operator<=(const OrthogPoly<T,Storage>& a,
490  const OrthogPoly<T,Storage>& b);
491 
492  template <typename T, typename Storage> bool
494  const OrthogPoly<T,Storage>& b);
495 
496  template <typename T, typename Storage> bool
497  operator<=(const OrthogPoly<T,Storage>& a,
498  const typename OrthogPoly<T,Storage>::value_type& b);
499 
500  template <typename T, typename Storage> bool
501  operator>=(const OrthogPoly<T,Storage>& a,
502  const OrthogPoly<T,Storage>& b);
503 
504  template <typename T, typename Storage> bool
506  const OrthogPoly<T,Storage>& b);
507 
508  template <typename T, typename Storage> bool
509  operator>=(const OrthogPoly<T,Storage>& a,
510  const typename OrthogPoly<T,Storage>::value_type& b);
511 
512  template <typename T, typename Storage> bool
513  operator<(const OrthogPoly<T,Storage>& a,
514  const OrthogPoly<T,Storage>& b);
515 
516  template <typename T, typename Storage> bool
518  const OrthogPoly<T,Storage>& b);
519 
520  template <typename T, typename Storage> bool
521  operator<(const OrthogPoly<T,Storage>& a,
522  const typename OrthogPoly<T,Storage>::value_type& b);
523 
524  template <typename T, typename Storage> bool
525  operator>(const OrthogPoly<T,Storage>& a,
526  const OrthogPoly<T,Storage>& b);
527 
528  template <typename T, typename Storage> bool
530  const OrthogPoly<T,Storage>& b);
531 
532  template <typename T, typename Storage> bool
533  operator>(const OrthogPoly<T,Storage>& a,
534  const typename OrthogPoly<T,Storage>::value_type& b);
535 
536  template <typename T, typename Storage> std::ostream&
537  operator << (std::ostream& os, const OrthogPoly<T,Storage>& a);
538 
539  template <typename T, typename Storage> std::istream&
540  operator >> (std::istream& os, OrthogPoly<T,Storage>& a);
541 
542  } // namespace PCE
543 
544 } // namespace Sacado
545 
548 
549 #endif // HAVE_STOKHOS_SACADO
550 
551 #endif // SACADO_PCE_ORTHOGPOLY_HPP
OrthogPoly< T, Storage > exp(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > log(const OrthogPoly< T, Storage > &a)
Stokhos::StandardStorage< int, double > storage_type
OrthogPoly< T, Storage > sin(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > sqrt(const OrthogPoly< T, Storage > &a)
bool operator>=(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
OrthogPoly< T, Storage > operator-(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
OrthogPoly< T, Storage > pow(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
OrthogPoly< T, Storage > atan(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > cbrt(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > acos(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > atanh(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > cosh(const OrthogPoly< T, Storage > &a)
std::istream & operator>>(std::istream &is, OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > sinh(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > tan(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > asin(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > operator+(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
OrthogPoly< T, Storage > operator/(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
OrthogPoly< T, Storage > max(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
atan2(expr1.val(), expr2.val())
Stokhos::LegendreBasis< int, double > basis_type
Abstract base class for orthogonal polynomial-based expansions.
bool operator!=(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
Abstract base class for multivariate orthogonal polynomials.
OrthogPoly< T, Storage > cos(const OrthogPoly< T, Storage > &a)
bool operator==(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
OrthogPoly< T, Storage > acosh(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > min(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
expr1 expr1 expr1 expr2 expr1 expr1 c expr2 expr1 c fastAccessCoeff(j)-expr2.val(j)
OrthogPoly< T, Storage > log10(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > abs(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > fabs(const OrthogPoly< T, Storage > &a)
Class to store coefficients of a projection onto an orthogonal polynomial basis.
expr val()
OrthogPoly< T, Storage > asinh(const OrthogPoly< T, Storage > &a)
OrthogPoly< T, Storage > tanh(const OrthogPoly< T, Storage > &a)
bool operator>(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)
OrthogPoly< T, Storage > operator*(const OrthogPoly< T, Storage > &a, const OrthogPoly< T, Storage > &b)