42 #ifndef SACADO_PCE_ORTHOGPOLY_HPP 43 #define SACADO_PCE_ORTHOGPOLY_HPP 47 #ifdef HAVE_STOKHOS_SACADO 49 #include "Teuchos_RCP.hpp" 51 #include "Sacado_Traits.hpp" 52 #include "Sacado_Handle.hpp" 53 #include "Sacado_mpl_apply.hpp" 73 template <
typename T,
typename Storage >
98 typedef typename approx_type::pointer pointer;
99 typedef typename approx_type::const_pointer const_pointer;
100 typedef typename approx_type::reference reference;
101 typedef typename approx_type::const_reference const_reference;
104 template <
typename S>
106 typedef typename Sacado::mpl::apply<Storage,ordinal_type,S>::type
storage_type;
107 typedef OrthogPoly<S,storage_type> type;
126 OrthogPoly(
const Teuchos::RCP<expansion_type>& expansion);
132 OrthogPoly(
const Teuchos::RCP<expansion_type>& expansion,
136 OrthogPoly(
const OrthogPoly&
x);
142 void init(
const T& v) { th->init(v); }
145 void init(
const T* v) { th->init(v); }
148 template <
typename S>
149 void init(
const OrthogPoly<T,S>& v) { th->init(v.getOrthogPolyApprox()); }
152 void load(T* v) { th->load(v); }
155 template <
typename S>
156 void load(OrthogPoly<T,S>& v) { th->load(v.getOrthogPolyApprox()); }
162 void reset(
const Teuchos::RCP<expansion_type>& expansion);
168 void reset(
const Teuchos::RCP<expansion_type>& expansion,
181 void copyForWrite() { th.makeOwnCopy(); }
184 value_type evaluate(
const Teuchos::Array<value_type>& point)
const;
187 value_type evaluate(
const Teuchos::Array<value_type>& point,
188 const Teuchos::Array<value_type>& bvals)
const;
194 value_type standard_deviation()
const {
return th->standard_deviation(); }
197 value_type two_norm()
const {
return th->two_norm(); }
200 value_type two_norm_squared()
const {
return th->two_norm_squared(); }
203 value_type inner_product(
const OrthogPoly& b)
const {
204 return th->inner_product(b.getOrthogPolyApprox()); }
207 std::ostream& print(std::ostream& os)
const {
return th->print(os); }
210 bool isEqualTo(
const OrthogPoly&
x)
const;
221 OrthogPoly<T,Storage>& operator=(
const OrthogPoly<T,Storage>&
x);
231 Teuchos::RCP<const basis_type> basis()
const {
return th->basis(); }
234 Teuchos::RCP<expansion_type> expansion()
const {
return expansion_; }
244 const_reference
val()
const {
return (*th)[0]; }
247 reference
val() {
return (*th)[0]; }
260 bool hasFastAccess(
ordinal_type sz)
const {
return th->size()>=sz;}
263 const_pointer coeff()
const {
return th->coeff();}
266 pointer coeff() {
return th->coeff();}
280 return th->term(dimension, order); }
284 return th->term(dimension, order); }
287 Teuchos::Array<ordinal_type> order(
ordinal_type term)
const {
288 return th->order(term); }
304 OrthogPoly<T,Storage>& operator += (
const value_type&
x);
307 OrthogPoly<T,Storage>& operator -= (
const value_type&
x);
310 OrthogPoly<T,Storage>& operator *= (
const value_type&
x);
313 OrthogPoly<T,Storage>& operator /= (
const value_type&
x);
316 OrthogPoly<T,Storage>& operator += (
const OrthogPoly<T,Storage>&
x);
319 OrthogPoly<T,Storage>& operator -= (
const OrthogPoly<T,Storage>&
x);
322 OrthogPoly<T,Storage>& operator *= (
const OrthogPoly<T,Storage>&
x);
325 OrthogPoly<T,Storage>& operator /= (
const OrthogPoly<T,Storage>&
x);
330 const approx_type& getOrthogPolyApprox()
const {
return *th; }
333 approx_type& getOrthogPolyApprox() {
return *th; }
338 Teuchos::RCP<expansion_type> expansion_;
341 static Teuchos::RCP<expansion_type> const_expansion_;
343 Sacado::Handle< Stokhos::OrthogPolyApprox<int,value_type,Storage> > th;
347 template <
typename T,
typename Storage>
348 Teuchos::RCP<typename OrthogPoly<T, Storage>::expansion_type>
349 OrthogPoly<T,Storage>::const_expansion_ =
353 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
354 operator+(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
356 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
358 const OrthogPoly<T,Storage>& b);
360 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
361 operator+(
const OrthogPoly<T,Storage>& a,
364 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
365 operator-(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
367 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
369 const OrthogPoly<T,Storage>& b);
371 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
372 operator-(
const OrthogPoly<T,Storage>& a,
375 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
376 operator*(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
378 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
380 const OrthogPoly<T,Storage>& b);
382 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
383 operator*(
const OrthogPoly<T,Storage>& a,
386 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
387 operator/(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
389 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
391 const OrthogPoly<T,Storage>& b);
393 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
394 operator/(
const OrthogPoly<T,Storage>& a,
397 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
398 exp(
const OrthogPoly<T,Storage>& a);
400 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
401 log(
const OrthogPoly<T,Storage>& a);
403 template <
typename T,
typename Storage>
void 404 log(OrthogPoly<T,Storage>& c,
const OrthogPoly<T,Storage>& a);
406 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
407 log10(
const OrthogPoly<T,Storage>& a);
409 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
410 sqrt(
const OrthogPoly<T,Storage>& a);
412 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
413 cbrt(
const OrthogPoly<T,Storage>& a);
415 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
416 pow(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
418 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
420 const OrthogPoly<T,Storage>& b);
422 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
423 pow(
const OrthogPoly<T,Storage>& a,
426 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
427 cos(
const OrthogPoly<T,Storage>& a);
429 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
430 sin(
const OrthogPoly<T,Storage>& a);
432 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
433 tan(
const OrthogPoly<T,Storage>& a);
435 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
436 cosh(
const OrthogPoly<T,Storage>& a);
438 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
439 sinh(
const OrthogPoly<T,Storage>& a);
441 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
442 tanh(
const OrthogPoly<T,Storage>& a);
444 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
445 acos(
const OrthogPoly<T,Storage>& a);
447 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
448 asin(
const OrthogPoly<T,Storage>& a);
450 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
451 atan(
const OrthogPoly<T,Storage>& a);
453 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
454 atan2(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
456 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
458 const OrthogPoly<T,Storage>& b);
460 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
461 atan2(
const OrthogPoly<T,Storage>& a,
464 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
465 acosh(
const OrthogPoly<T,Storage>& a);
467 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
468 asinh(
const OrthogPoly<T,Storage>& a);
470 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
471 atanh(
const OrthogPoly<T,Storage>& a);
473 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
474 abs(
const OrthogPoly<T,Storage>& a);
476 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
477 fabs(
const OrthogPoly<T,Storage>& a);
479 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
480 max(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
482 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
484 const OrthogPoly<T,Storage>& b);
486 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
487 max(
const OrthogPoly<T,Storage>& a,
490 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
491 min(
const OrthogPoly<T,Storage>& a,
const OrthogPoly<T,Storage>& b);
493 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
495 const OrthogPoly<T,Storage>& b);
497 template <
typename T,
typename Storage> OrthogPoly<T,Storage>
498 min(
const OrthogPoly<T,Storage>& a,
501 template <
typename T,
typename Storage>
bool 503 const OrthogPoly<T,Storage>& b);
505 template <
typename T,
typename Storage>
bool 507 const OrthogPoly<T,Storage>& b);
509 template <
typename T,
typename Storage>
bool 513 template <
typename T,
typename Storage>
bool 515 const OrthogPoly<T,Storage>& b);
517 template <
typename T,
typename Storage>
bool 519 const OrthogPoly<T,Storage>& b);
521 template <
typename T,
typename Storage>
bool 525 template <
typename T,
typename Storage>
bool 526 operator<=(const OrthogPoly<T,Storage>& a,
527 const OrthogPoly<T,Storage>& b);
529 template <
typename T,
typename Storage>
bool 531 const OrthogPoly<T,Storage>& b);
533 template <
typename T,
typename Storage>
bool 534 operator<=(const OrthogPoly<T,Storage>& a,
537 template <
typename T,
typename Storage>
bool 539 const OrthogPoly<T,Storage>& b);
541 template <
typename T,
typename Storage>
bool 543 const OrthogPoly<T,Storage>& b);
545 template <
typename T,
typename Storage>
bool 549 template <
typename T,
typename Storage>
bool 550 operator<(const OrthogPoly<T,Storage>& a,
551 const OrthogPoly<T,Storage>& b);
553 template <
typename T,
typename Storage>
bool 555 const OrthogPoly<T,Storage>& b);
557 template <
typename T,
typename Storage>
bool 558 operator<(const OrthogPoly<T,Storage>& a,
561 template <
typename T,
typename Storage>
bool 562 operator>(
const OrthogPoly<T,Storage>& a,
563 const OrthogPoly<T,Storage>& b);
565 template <
typename T,
typename Storage>
bool 567 const OrthogPoly<T,Storage>& b);
569 template <
typename T,
typename Storage>
bool 570 operator>(
const OrthogPoly<T,Storage>& a,
573 template <
typename T,
typename Storage> std::ostream&
574 operator << (std::ostream& os, const OrthogPoly<T,Storage>& a);
576 template <
typename T,
typename Storage> std::istream&
577 operator >> (std::istream& os, OrthogPoly<T,Storage>& a);
586 #endif // HAVE_STOKHOS_SACADO 588 #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)
expr1 expr1 expr1 expr2 expr1 expr1 c expr2 expr1 c fastAccessCoeff(j) - expr2.val(j)
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.
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType * x
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)
Orthogonal polynomial expansion class for constant (size 1) expansions.
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.
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)