44 #include "Teuchos_UnitTestHarness.hpp" 45 #include "Teuchos_TestingHelpers.hpp" 46 #include "Teuchos_UnitTestRepository.hpp" 47 #include "Teuchos_GlobalMPISession.hpp" 55 template <
typename OrdinalType,
typename ValueType>
60 Teuchos::RCP<const Stokhos::CompletePolynomialBasis<OrdinalType,ValueType> >
basis;
61 Teuchos::RCP<Teuchos::SerialDenseMatrix<int,double> >
Bij;
62 Teuchos::RCP<Stokhos::Sparse3Tensor<int,double> >
Cijk;
63 Teuchos::RCP<Stokhos::Dense3Tensor<int,double> >
Dijk;
64 Teuchos::RCP<const Stokhos::Quadrature<OrdinalType,ValueType> >
quad;
65 Teuchos::RCP< Stokhos::DerivOrthogPolyExpansion<OrdinalType,ValueType> >
exp;
66 Stokhos::OrthogPolyApprox<OrdinalType,ValueType> x,
y,
u,
u2,
cx,
cu,
cu2,
sx,
su,
su2;
75 const OrdinalType d = 1;
76 const OrdinalType p = 7;
79 Teuchos::Array< Teuchos::RCP<const Stokhos::OneDOrthogPolyBasis<OrdinalType,ValueType> > > bases(d);
80 for (OrdinalType i=0; i<d; i++)
91 Bij =
basis->computeDerivDoubleProductTensor();
92 Cijk =
basis->computeTripleProductTensor();
113 for (OrdinalType i=0; i<d; i++) {
118 for (OrdinalType i=0; i<d; i++)
122 template <
class Func>
127 const Teuchos::Array<ValueType>& weights =
quad->getQuadWeights();
128 const Teuchos::Array< Teuchos::Array<ValueType> >& points =
129 quad->getQuadPoints();
130 const Teuchos::Array< Teuchos::Array<ValueType> >& values =
131 quad->getBasisAtQuadPoints();
132 OrdinalType nqp = weights.size();
135 for (OrdinalType i=0; i<c.
size(); i++)
140 for (OrdinalType k=0; k<nqp; k++) {
141 ValueType
val =
a.evaluate(points[k], values[k]);
143 for (
int i=0; i<c.
size(); i++)
144 c[i] += weights[k]*
val*values[k][i] /
basis->norm_squared(i);
148 template <
class Func>
154 const Teuchos::Array<ValueType>& weights =
quad->getQuadWeights();
155 const Teuchos::Array< Teuchos::Array<ValueType> >& points =
156 quad->getQuadPoints();
157 const Teuchos::Array< Teuchos::Array<ValueType> >& values =
158 quad->getBasisAtQuadPoints();
159 OrdinalType nqp = weights.size();
162 for (OrdinalType i=0; i<c.
size(); i++)
167 for (OrdinalType k=0; k<nqp; k++) {
168 ValueType val1 =
a.evaluate(points[k], values[k]);
169 ValueType val2 = b.
evaluate(points[k], values[k]);
170 ValueType
val = func(val1, val2);
171 for (
int i=0; i<c.
size(); i++)
172 c[i] += weights[k]*
val*values[k][i] /
basis->norm_squared(i);
176 template <
class Func>
183 const Teuchos::Array<ValueType>& weights =
quad->getQuadWeights();
184 const Teuchos::Array< Teuchos::Array<ValueType> >& points =
185 quad->getQuadPoints();
186 const Teuchos::Array< Teuchos::Array<ValueType> >& values =
187 quad->getBasisAtQuadPoints();
188 OrdinalType nqp = weights.size();
191 for (OrdinalType i=0; i<c.
size(); i++)
196 for (OrdinalType k=0; k<nqp; k++) {
197 ValueType val2 = b.
evaluate(points[k], values[k]);
198 ValueType
val = func(
a, val2);
199 for (
int i=0; i<c.
size(); i++)
200 c[i] += weights[k]*
val*values[k][i] /
basis->norm_squared(i);
204 template <
class Func>
211 const Teuchos::Array<ValueType>& weights =
quad->getQuadWeights();
212 const Teuchos::Array< Teuchos::Array<ValueType> >& points =
213 quad->getQuadPoints();
214 const Teuchos::Array< Teuchos::Array<ValueType> >& values =
215 quad->getBasisAtQuadPoints();
216 OrdinalType nqp = weights.size();
219 for (OrdinalType i=0; i<c.
size(); i++)
224 for (OrdinalType k=0; k<nqp; k++) {
225 ValueType val1 =
a.evaluate(points[k], values[k]);
226 ValueType
val = func(val1, b);
227 for (
int i=0; i<c.
size(); i++)
228 c[i] += weights[k]*
val*values[k][i] /
basis->norm_squared(i);
293 return 0.5*
std::log((1.0+a)/(1.0-a));
298 double operator() (
double a,
double b)
const {
return a + b; }
301 double operator() (
double a,
double b)
const {
return a - b; }
304 double operator() (
double a,
double b)
const {
return a * b; }
307 double operator() (
double a,
double b)
const {
return a / b; }
700 setup.exp->plus(ru, v, w);
816 setup.exp->minus(ru, v, w);
932 setup.exp->times(ru, v, w);
1048 setup.exp->divide(ru, v, w);
1233 setup.exp->plusEqual(ru, v);
1286 setup.exp->minusEqual(ru, v);
1287 setup.exp->unaryMinus(v, v);
1427 Teuchos::GlobalMPISession mpiSession(&argc, &
argv);
1428 return Teuchos::UnitTestRepository::runUnitTestsFromMain(argc,
argv);
double operator()(double a) const
KOKKOS_INLINE_FUNCTION PCE< Storage > sqrt(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > tan(const PCE< Storage > &a)
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > u2
double operator()(double a) const
double operator()(double a) const
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > cx
double operator()(double a) const
double operator()(double a) const
double operator()(double a, double b) const
KOKKOS_INLINE_FUNCTION PCE< Storage > sinh(const PCE< Storage > &a)
double operator()(double a) const
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > su2
Teuchos::RCP< const Stokhos::Quadrature< OrdinalType, ValueType > > quad
Teuchos::RCP< Stokhos::Dense3Tensor< int, double > > Dijk
KOKKOS_INLINE_FUNCTION PCE< Storage > pow(const PCE< Storage > &a, const PCE< Storage > &b)
void computePCE2(Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &c, const Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &a, const Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &b)
void computePCE2LC(Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &c, ValueType a, const Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &b)
bool comparePCEs(const PCEType &a1, const std::string &a1_name, const Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &a2, const std::string &a2_name, const ValueType &rel_tol, const ValueType &abs_tol, Teuchos::FancyOStream &out)
double operator()(double a) const
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > u
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > sx
KOKKOS_INLINE_FUNCTION PCE< Storage > tanh(const PCE< Storage > &a)
int main(int argc, char *argv[])
KOKKOS_INLINE_FUNCTION PCE< Storage > cbrt(const PCE< Storage > &a)
ordinal_type size() const
Return size.
double operator()(double a) const
double operator()(double a) const
KOKKOS_INLINE_FUNCTION PCE< Storage > acos(const PCE< Storage > &a)
double operator()(double a, double b) const
double operator()(double a) const
Teuchos::RCP< const Stokhos::CompletePolynomialBasis< OrdinalType, ValueType > > basis
Teuchos::RCP< Stokhos::DerivOrthogPolyExpansion< OrdinalType, ValueType > > exp
double operator()(double a) const
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > su
Teuchos::RCP< Stokhos::Sparse3Tensor< int, double > > Cijk
void reset(const Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > &new_basis, ordinal_type sz=0)
Reset to a new basis.
double operator()(double a, double b) const
double operator()(double a) const
KOKKOS_INLINE_FUNCTION PCE< Storage > cosh(const PCE< Storage > &a)
double operator()(double a, double b) const
Othogonal polynomial expansions based on derivative calculations.
TEUCHOS_UNIT_TEST(Stokhos_DerivExpansion, UMinus)
Multivariate orthogonal polynomial basis generated from a total-order complete-polynomial tensor prod...
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, double > > Bij
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > cu2
double operator()(double a) const
KOKKOS_INLINE_FUNCTION PCE< Storage > atan(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > exp(const PCE< Storage > &a)
double operator()(double a) const
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > y
UnitTestSetup< int, double > setup
double operator()(double a) const
ScalarType f(const Teuchos::Array< ScalarType > &x, double a, double b)
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > x
double operator()(double a) const
KOKKOS_INLINE_FUNCTION PCE< Storage > sin(const PCE< Storage > &a)
Stokhos::OrthogPolyApprox< OrdinalType, ValueType > cu
double operator()(double a, double b) const
double operator()(double a) const
double operator()(double a) const
value_type evaluate(const Teuchos::Array< value_type > &point) const
Evaluate polynomial approximation at a point.
KOKKOS_INLINE_FUNCTION PCE< Storage > log(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > log10(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > asin(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > cos(const PCE< Storage > &a)
Defines quadrature for a tensor product basis by tensor products of 1-D quadrature rules...
void computePCE2RC(Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &c, const Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &a, ValueType b)
reference term(ordinal_type dimension, ordinal_type order)
Get coefficient term for given dimension and order.
void computePCE1(Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &c, const Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &a)