42#include "Teuchos_Assert.hpp"
44#include "Teuchos_ConfigDefs.hpp"
46#define UQ_PREP_F77 F77_FUNC_(uq_prep,UQ_PREP)
47#define UQ_PROD2_F77 F77_FUNC_(uq_prod2,UQ_PROD2)
48#define UQ_DIV_F77 F77_FUNC_(uq_div,UQ_DIV)
49#define UQ_EXP_F77 F77_FUNC_(uq_exp,UQ_EXP)
50#define UQ_LOG_F77 F77_FUNC_(uq_log,UQ_LOG)
51#define UQ_SQRT_F77 F77_FUNC_(uq_sqrt,UQ_SQRT)
52#define UQ_EXP_INT_F77 F77_FUNC_(uq_exp_int,UQ_EXP)
53#define UQ_LOG_INT_F77 F77_FUNC_(uq_log_int,UQ_LOG)
58 void UQ_DIV_F77(
const double*,
const double*,
double*,
int*);
59 void UQ_EXP_F77(
const double*,
double*,
int*,
int*,
double*,
int*);
60 void UQ_LOG_F77(
const double*,
double*,
int*,
int*,
double*,
int*);
66template <
typename ordinal_type,
typename value_type>
67Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
68ForUQTKOrthogPolyExpansion(
71 EXPANSION_METHOD method_,
77 order = this->basis->order();
78 dim = this->basis->dimension();
84template <
typename ordinal_type,
typename value_type>
86Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
88 const value_type&
val)
90 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::timesEqual(c,
val);
93template <
typename ordinal_type,
typename value_type>
95Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
97 const value_type&
val)
99 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::divideEqual(c,
val);
102template <
typename ordinal_type,
typename value_type>
104Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
116 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
117 "Stokhos::ForUQTKOrthogPolyExpansion::timesEqual()" <<
118 ": Expansion size (" << sz <<
119 ") is too small for computation.");
126 if (p > 1 && xp > 1) {
127 TEUCHOS_TEST_FOR_EXCEPTION(pc != xp, std::logic_error,
128 "Stokhos::ForUQTKOrthogPolyExpansion::timesEqual()"
129 <<
": Arguments have incompatible sizes: "
130 <<
"x.size() = " << xp <<
", c.size() = " << pc <<
".");
140 for (ordinal_type i=0; i<p; i++)
144 for (ordinal_type i=1; i<xp; i++)
153template <
typename ordinal_type,
typename value_type>
155Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
157 const OrthogPolyApprox<ordinal_type, value_type, node_type>& x)
166 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
167 "Stokhos::ForUQTKOrthogPolyExpansion::divideEqual()" <<
168 ": Expansion size (" << sz <<
169 ") is too small for computation.");
177 TEUCHOS_TEST_FOR_EXCEPTION(pc != xp, std::logic_error,
178 "Stokhos::ForUQTKOrthogPolyExpansion::divideEqual()"
179 <<
": Arguments have incompatible sizes: "
180 <<
"x.size() = " << xp <<
", c.size() = " << pc <<
".");
190 for (ordinal_type i=0; i<pc; i++)
195template <
typename ordinal_type,
typename value_type>
197Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
205 if (pa > 1 && pb > 1)
209 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
210 "Stokhos::ForUQTKOrthogPolyExpansion::times()" <<
211 ": Expansion size (" << sz <<
212 ") is too small for computation.");
220 if (pa > 1 && pb > 1) {
221 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc || pb != pc, std::logic_error,
222 "Stokhos::ForUQTKOrthogPolyExpansion::times()"
223 <<
": Arguments have incompatible sizes: "
224 <<
"a.size() = " << pa <<
", b.size() = " << pb
225 <<
", required size = " << pc <<
".");
231 for (ordinal_type i=0; i<pc; i++)
235 for (ordinal_type i=0; i<pc; i++)
243template <
typename ordinal_type,
typename value_type>
245Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
250 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::times(c,a,b);
253template <
typename ordinal_type,
typename value_type>
255Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
260 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::times(c,a,b);
263template <
typename ordinal_type,
typename value_type>
265Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
277 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
278 "Stokhos::ForUQTKOrthogPolyExpansion::divide()" <<
279 ": Expansion size (" << sz <<
280 ") is too small for computation.");
289 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc || pb != pc, std::logic_error,
290 "Stokhos::ForUQTKOrthogPolyExpansion::divide()"
291 <<
": Arguments have incompatible sizes: "
292 <<
"a.size() = " << pa <<
", b.size() = " << pb
293 <<
", required size = " << pc <<
".");
299 for (ordinal_type i=0; i<pa; i++)
304template <
typename ordinal_type,
typename value_type>
306Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
324 TEUCHOS_TEST_FOR_EXCEPTION(pb != pc, std::logic_error,
325 "Stokhos::ForUQTKOrthogPolyExpansion::divide()"
326 <<
": Arguments have incompatible sizes: "
327 <<
"b.size() = " << pb
328 <<
", required size = " << pc <<
".");
339template <
typename ordinal_type,
typename value_type>
341Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
346 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::divide(c,a,b);
349template <
typename ordinal_type,
typename value_type>
351Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
368 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
369 "Stokhos::ForUQTKOrthogPolyExpansion::exp()"
370 <<
": Arguments have incompatible sizes: "
371 <<
"a.size() = " << pa <<
", c.size() = " << pc
374 if (method == TAYLOR) {
382 cc[0] = std::exp(ca[0]);
385template <
typename ordinal_type,
typename value_type>
387Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
404 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
405 "Stokhos::ForUQTKOrthogPolyExpansion::log()"
406 <<
": Arguments have incompatible sizes: "
407 <<
"a.size() = " << pa <<
", c.size() = " << pc
410 if (method == TAYLOR) {
418 cc[0] = std::log(ca[0]);
421template <
typename ordinal_type,
typename value_type>
423Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
429 divide(c,c,std::log(10.0));
434 c[0] = std::log10(a[0]);
438template <
typename ordinal_type,
typename value_type>
440Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
457 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
458 "Stokhos::ForUQTKOrthogPolyExpansion::sqrt()"
459 <<
": Arguments have incompatible sizes: "
460 <<
"a.size() = " << pa <<
", c.size() = " << pc
467 cc[0] = std::sqrt(ca[0]);
470template <
typename ordinal_type,
typename value_type>
472Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
484 c[0] = std::cbrt(a[0]);
488template <
typename ordinal_type,
typename value_type>
490Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
503 c[0] = std::pow(a[0], b[0]);
507template <
typename ordinal_type,
typename value_type>
509Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
515 times(c,std::log(a),b);
521 c[0] = std::pow(a, b[0]);
525template <
typename ordinal_type,
typename value_type>
527Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
540 c[0] = std::pow(a[0], b);
544template <
typename ordinal_type,
typename value_type>
546Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
553 s[0] = std::sin(a[0]);
556 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
557 "Stokhos::ForUQTKOrthogPolyExpansion::sin()"
558 <<
": Method not implemented!");
561template <
typename ordinal_type,
typename value_type>
563Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
570 c[0] = std::cos(a[0]);
573 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
574 "Stokhos::ForUQTKOrthogPolyExpansion::cos()"
575 <<
": Method not implemented!");
578template <
typename ordinal_type,
typename value_type>
580Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
587 t[0] = std::tan(a[0]);
590 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
591 "Stokhos::ForUQTKOrthogPolyExpansion::tan()"
592 <<
": Method not implemented!");
595template <
typename ordinal_type,
typename value_type>
597Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
607 this->minusEqual(s, t);
613 s[0] = std::sinh(a[0]);
617template <
typename ordinal_type,
typename value_type>
619Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
629 this->plusEqual(c, t);
635 c[0] = std::cosh(a[0]);
639template <
typename ordinal_type,
typename value_type>
641Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
652 this->minus(t, c, s);
653 this->plusEqual(c, s);
659 t[0] = std::tanh(a[0]);
663template <
typename ordinal_type,
typename value_type>
665Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
672 c[0] = std::acos(a[0]);
675 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
676 "Stokhos::ForUQTKOrthogPolyExpansion::acos()"
677 <<
": Method not implemented!");
680template <
typename ordinal_type,
typename value_type>
682Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
689 c[0] = std::asin(a[0]);
692 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
693 "Stokhos::ForUQTKOrthogPolyExpansion::asin()"
694 <<
": Method not implemented!");
697template <
typename ordinal_type,
typename value_type>
699Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
706 c[0] = std::atan(a[0]);
709 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
710 "Stokhos::ForUQTKOrthogPolyExpansion::atan()"
711 <<
": Method not implemented!");
714template <
typename ordinal_type,
typename value_type>
716Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
721 if (a.
size() == 1 && b.
size() == 1) {
724 c[0] = std::atan2(a[0], b[0]);
727 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
728 "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
729 <<
": Method not implemented!");
732template <
typename ordinal_type,
typename value_type>
734Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
742 c[0] = std::atan2(a, b[0]);
745 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
746 "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
747 <<
": Method not implemented!");
750template <
typename ordinal_type,
typename value_type>
752Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
760 c[0] = std::atan2(a[0], b);
763 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
764 "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
765 <<
": Method not implemented!");
768template <
typename ordinal_type,
typename value_type>
770Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
777 c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]-
value_type(1.0)));
780 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
781 "Stokhos::ForUQTKOrthogPolyExpansion::acosh()"
782 <<
": Method not implemented!");
785template <
typename ordinal_type,
typename value_type>
787Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
794 c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]+
value_type(1.0)));
797 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
798 "Stokhos::ForUQTKOrthogPolyExpansion::asinh()"
799 <<
": Method not implemented!");
802template <
typename ordinal_type,
typename value_type>
804Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
814 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
815 "Stokhos::ForUQTKOrthogPolyExpansion::atanh()"
816 <<
": Method not implemented!");
Class to store coefficients of a projection onto an orthogonal polynomial basis.
void resize(ordinal_type sz)
Resize coefficient array (coefficients are preserved)
pointer coeff()
Return coefficient array.
ordinal_type size() const
Return size.
Abstract base class for multivariate orthogonal polynomials.
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType * x
static void copy(const T *src, T *dest, int sz)
Copy array from src to dest of length sz.
static T * get_and_fill(int sz)
Get memory for new array of length sz and fill with zeros.