Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_OrthogPolyApproxImp.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 #include "Stokhos_ProductBasis.hpp"
11 
12 template <typename ordinal_type, typename value_type, typename storage_type>
13 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
14 OrthogPolyApprox(
15  const Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >& theBasis,
16  ordinal_type sz,
17  const value_type* vals) :
18  basis_(theBasis),
19  coeff_(0)
20 {
21  if (sz != 0)
22  coeff_.resize(sz);
23  else if (theBasis != Teuchos::null)
24  coeff_.resize(theBasis->size());
25  else
26  coeff_.resize(ordinal_type(1));
27  if (vals != NULL)
28  coeff_.init(vals);
29 }
30 
31 template <typename ordinal_type, typename value_type, typename storage_type>
32 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
33 OrthogPolyApprox(const Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>& x) :
34  basis_(x.basis_),
35  coeff_(x.coeff_)
36 {
37 }
38 
39 template <typename ordinal_type, typename value_type, typename storage_type>
40 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
41 ~OrthogPolyApprox()
42 {
43 }
44 
45 template <typename ordinal_type, typename value_type, typename storage_type>
46 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>&
47 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
48 operator=(const Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>& x)
49 {
50  if (this != &x) {
51  basis_ = x.basis_;
52  coeff_ = x.coeff_;
53  }
54 
55  return *this;
56 }
57 
58 template <typename ordinal_type, typename value_type, typename storage_type>
59 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>&
60 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
61 operator=(const value_type& v)
62 {
63  coeff_.init(value_type(0));
64  coeff_.init(&v, 1);
65 
66  return *this;
67 }
68 
69 template <typename ordinal_type, typename value_type, typename storage_type>
70 void
71 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
72 init(const value_type& v)
73 {
74  coeff_.init(v);
75 }
76 
77 template <typename ordinal_type, typename value_type, typename storage_type>
78 void
79 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
80 init(const value_type* v)
81 {
82  coeff_.init(v);
83 }
84 
85 template <typename ordinal_type, typename value_type, typename storage_type>
86 void
87 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
88 load(value_type* v)
89 {
90  coeff_.load(v);
91 }
92 
93 template <typename ordinal_type, typename value_type, typename storage_type>
94 Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >
95 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
96 basis() const
97 {
98  return basis_;
99 }
100 
101 template <typename ordinal_type, typename value_type, typename storage_type>
102 void
103 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
104 reset(const Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >& new_basis, ordinal_type sz)
105 {
106  basis_ = new_basis;
107  if (sz != 0)
108  resize(sz);
109  else if (basis_ != Teuchos::null)
110  resize(basis_->size());
111  else
112  resize(1);
113 }
114 
115 template <typename ordinal_type, typename value_type, typename storage_type>
116 void
117 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
118 resize(ordinal_type sz)
119 {
120  coeff_.resize(sz);
121 }
122 
123 template <typename ordinal_type, typename value_type, typename storage_type>
124 ordinal_type
125 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
126 size() const
127 {
128  return coeff_.size();
129 }
130 
131 template <typename ordinal_type, typename value_type, typename storage_type>
132 typename Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::pointer
133 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
134 coeff()
135 {
136 #ifdef STOKHOS_DEBUG
137  TEUCHOS_TEST_FOR_EXCEPTION(coeff_.size() == 0, std::logic_error,
138  "Stokhos::OrthogPolyApprox::coeff(): " <<
139  "Coefficient array is empty!");
140 #endif
141  return coeff_.coeff();
142 }
143 
144 template <typename ordinal_type, typename value_type, typename storage_type>
145 typename Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::const_pointer
146 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
147 coeff() const
148 {
149 #ifdef STOKHOS_DEBUG
150  TEUCHOS_TEST_FOR_EXCEPTION(coeff_.size() == 0, std::logic_error,
151  "Stokhos::OrthogPolyApprox::coeff(): " <<
152  "Coefficient array is empty!");
153 #endif
154  return coeff_.coeff();
155 }
156 
157 template <typename ordinal_type, typename value_type, typename storage_type>
158 typename Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::reference
159 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
160 operator[](ordinal_type i)
161 {
162  return coeff_[i];
163 }
164 
165 template <typename ordinal_type, typename value_type, typename storage_type>
166 typename Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::const_reference
167 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
168 operator[](ordinal_type i) const
169 {
170  return coeff_[i];
171 }
172 
173 template <typename ordinal_type, typename value_type, typename storage_type>
174 typename Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::reference
175 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
176 term(ordinal_type dimension, ordinal_type theOrder)
177 {
178  Teuchos::RCP< const Stokhos::ProductBasis<ordinal_type, value_type> >
179  product_basis = Teuchos::rcp_dynamic_cast< const Stokhos::ProductBasis<ordinal_type, value_type> >(basis_, true);
180  ordinal_type d = basis_->dimension();
181  MultiIndex<ordinal_type> theTerm(d);
182  theTerm[dimension] = theOrder;
183  ordinal_type index = product_basis->index(theTerm);
184  return coeff_[index];
185 }
186 
187 template <typename ordinal_type, typename value_type, typename storage_type>
188 typename Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::const_reference
189 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
190 term(ordinal_type dimension, ordinal_type theOrder) const
191 {
192  Teuchos::RCP< const Stokhos::ProductBasis<ordinal_type, value_type> >
193  product_basis = Teuchos::rcp_dynamic_cast< const Stokhos::ProductBasis<ordinal_type, value_type> >(basis_, true);
194  ordinal_type d = basis_->dimension();
195  MultiIndex<ordinal_type> theTerm(d);
196  theTerm[dimension] = theOrder;
197  ordinal_type index = product_basis->index(theTerm);
198  return coeff_[index];
199 }
200 
201 template <typename ordinal_type, typename value_type, typename storage_type>
202 const Stokhos::MultiIndex<ordinal_type>&
203 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
204 order(ordinal_type theTerm) const
205 {
206  Teuchos::RCP< const Stokhos::ProductBasis<ordinal_type, value_type> >
207  product_basis = Teuchos::rcp_dynamic_cast< const Stokhos::ProductBasis<ordinal_type, value_type> >(basis_, true);
208  return product_basis->term(theTerm);
209 }
210 
211 template <typename ordinal_type, typename value_type, typename storage_type>
212 value_type
213 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
214 evaluate(const Teuchos::Array<value_type>& point) const
215 {
216  Teuchos::Array<value_type> basis_vals(basis_->size());
217  basis_->evaluateBases(point, basis_vals);
218  return evaluate(point, basis_vals);
219 }
220 
221 template <typename ordinal_type, typename value_type, typename storage_type>
222 value_type
223 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
224 evaluate(const Teuchos::Array<value_type>& point,
225  const Teuchos::Array<value_type>& basis_vals) const
226 {
227  value_type val = value_type(0.0);
228  for (ordinal_type i=0; i<static_cast<ordinal_type>(coeff_.size()); i++)
229  val += coeff_[i]*basis_vals[i];
230 
231  return val;
232 }
233 
234 template <typename ordinal_type, typename value_type, typename storage_type>
235 value_type
236 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
237 mean() const
238 {
239  return coeff_[0];
240 }
241 
242 template <typename ordinal_type, typename value_type, typename storage_type>
243 value_type
244 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
245 standard_deviation() const
246 {
247  value_type std_dev = 0.0;
248  for (ordinal_type i=1; i<static_cast<ordinal_type>(coeff_.size()); i++) {
249  std_dev += coeff_[i]*coeff_[i]*basis_->norm_squared(i);
250  }
251  std_dev = std::sqrt(std_dev);
252  return std_dev;
253 }
254 
255 template <typename ordinal_type, typename value_type, typename storage_type>
256 value_type
257 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
258 two_norm() const
259 {
260  return std::sqrt(this->two_norm_squared());
261 }
262 
263 template <typename ordinal_type, typename value_type, typename storage_type>
264 value_type
265 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
266 two_norm_squared() const
267 {
268  value_type nrm = 0.0;
269  if (basis_ == Teuchos::null) { // Check for special case of constants
270  TEUCHOS_TEST_FOR_EXCEPTION(
271  coeff_.size() != 1, std::logic_error,
272  "basis_ == null && coeff_.size() > 1");
273  nrm = coeff_[0]*coeff_[0];
274  }
275  else {
276  for (ordinal_type i=0; i<static_cast<ordinal_type>(coeff_.size()); i++)
277  nrm += coeff_[i]*coeff_[i]*basis_->norm_squared(i);
278  }
279  return nrm;
280 }
281 
282 template <typename ordinal_type, typename value_type, typename storage_type>
283 value_type
284 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
285 inner_product(const Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>& b) const
286 {
287  // Check a and b are compatible
288  TEUCHOS_TEST_FOR_EXCEPTION(
289  basis_ == Teuchos::null && coeff_.size() != 1, std::logic_error,
290  "basis_ == null && coeff_.size() > 1");
291  TEUCHOS_TEST_FOR_EXCEPTION(
292  b.basis_ == Teuchos::null && b.coeff_.size() != 1, std::logic_error,
293  "b.basis_ == null && b.coeff_.size() > 1");
294  TEUCHOS_TEST_FOR_EXCEPTION(
295  coeff_.size() != b.coeff_.size() &&
296  coeff_.size() != 1 && b.coeff_.size() != 1, std::logic_error,
297  "Coefficient array sizes do not match");
298 
299  value_type v = 0.0;
300  if (coeff_.size() == 1 || b.coeff_.size() == 1)
301  v = coeff_[0]*b.coeff_[0];
302  else
303  for (ordinal_type i=0; i<static_cast<ordinal_type>(coeff_.size()); i++)
304  v += coeff_[i]*b.coeff_[i]*basis_->norm_squared(i);
305 
306  return v;
307 }
308 
309 template <typename ordinal_type, typename value_type, typename storage_type>
310 std::ostream&
311 Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>::
312 print(std::ostream& os) const
313 {
314  os << "Stokhos::OrthogPolyApprox of size " << coeff_.size() << " in basis "
315  << "\n\t" << basis_->getName() << ":" << std::endl;
316 
317  Teuchos::RCP< const Stokhos::ProductBasis<ordinal_type, value_type> >
318  product_basis = Teuchos::rcp_dynamic_cast< const Stokhos::ProductBasis<ordinal_type, value_type> >(basis_);
319 
320  if (product_basis != Teuchos::null) {
321  for (ordinal_type i=0; i<static_cast<ordinal_type>(coeff_.size()); i++) {
322  const Stokhos::MultiIndex<ordinal_type>& trm = product_basis->term(i);
323  os << "\t\t(";
324  for (ordinal_type j=0; j<static_cast<ordinal_type>(trm.size())-1; j++)
325  os << trm[j] << ", ";
326  os << trm[trm.size()-1] << ") = " << coeff_[i] << std::endl;
327  }
328  }
329  else {
330  os << "[ ";
331 
332  for (ordinal_type i=0; i<static_cast<ordinal_type>(coeff_.size()); i++) {
333  os << coeff_[i] << " ";
334  }
335 
336  os << "]\n";
337  }
338 
339  return os;
340 }
341 
342 template <typename ordinal_type, typename value_type, typename storage_type>
343 std::ostream&
344 Stokhos::operator << (
345  std::ostream& os,
346  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, storage_type>& a)
347 {
348  os << "[ ";
349 
350  for (ordinal_type i=0; i<static_cast<ordinal_type>(a.size()); i++) {
351  os << a[i] << " ";
352  }
353 
354  os << "]\n";
355  return os;
356 }
Sacado::UQ::sqrt
KOKKOS_INLINE_FUNCTION PCE< Storage > sqrt(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:847
Stokhos::OrthogPolyApprox::two_norm_squared
value_type two_norm_squared() const
Compute the squared two-norm of expansion.
Definition: Stokhos_OrthogPolyApproxImp.hpp:266
Stokhos::OrthogPolyApprox::resize
void resize(ordinal_type sz)
Resize coefficient array (coefficients are preserved)
Definition: Stokhos_OrthogPolyApproxImp.hpp:118
Stokhos::OrthogPolyApprox::evaluate
value_type evaluate(const Teuchos::Array< value_type > &point) const
Evaluate polynomial approximation at a point.
Definition: Stokhos_OrthogPolyApproxImp.hpp:214
Stokhos::OrthogPolyApprox::basis_
Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > basis_
Basis expansion is relative to.
Definition: Stokhos_OrthogPolyApprox.hpp:147
Stokhos::OrthogPolyApprox::operator=
OrthogPolyApprox & operator=(const OrthogPolyApprox &x)
Assignment operator (deep copy)
TEUCHOS_TEST_FOR_EXCEPTION
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Stokhos::MultiIndex::size
ordinal_type size() const
Size.
Definition: Stokhos_ProductBasisUtils.hpp:67
Stokhos::OrthogPolyApprox::init
void init(const value_type &v)
Initialize coefficients to value.
Definition: Stokhos_OrthogPolyApproxImp.hpp:72
Stokhos::OrthogPolyApprox::OrthogPolyApprox
OrthogPolyApprox(const Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > &basis=Teuchos::null, ordinal_type sz=0, const value_type *vals=NULL)
Constructor with supplied size sz.
Definition: Stokhos_OrthogPolyApproxImp.hpp:14
Stokhos::OrthogPolyApprox::coeff
pointer coeff()
Return coefficient array.
Definition: Stokhos_OrthogPolyApproxImp.hpp:134
Stokhos_ProductBasis.hpp
Stokhos::OrthogPolyApprox::two_norm
value_type two_norm() const
Compute the two-norm of expansion.
Definition: Stokhos_OrthogPolyApproxImp.hpp:258
j
j
Definition: Sacado_Fad_Exp_MP_Vector.hpp:495
Stokhos::OrthogPolyApprox::inner_product
value_type inner_product(const OrthogPolyApprox &b) const
Compute the L2 inner product of 2 PCEs.
Definition: Stokhos_OrthogPolyApproxImp.hpp:285
Stokhos::OrthogPolyApprox::standard_deviation
value_type standard_deviation() const
Compute standard deviation of expansion.
Definition: Stokhos_OrthogPolyApproxImp.hpp:245
Stokhos::OrthogPolyBasis
Abstract base class for multivariate orthogonal polynomials.
Definition: Stokhos_OrthogPolyBasis.hpp:40
Stokhos::OrthogPolyApprox::operator[]
reference operator[](ordinal_type i)
Array access.
Definition: Stokhos_OrthogPolyApproxImp.hpp:160
Stokhos::OrthogPolyApprox::reference
storage_type::reference reference
Definition: Stokhos_OrthogPolyApprox.hpp:32
Stokhos::operator<<
std::ostream & operator<<(std::ostream &os, const ProductContainer< coeff_type > &vec)
Definition: Stokhos_ProductContainer.hpp:152
Stokhos::OrthogPolyApprox::coeff_
storage_type coeff_
OrthogPolyApprox coefficients.
Definition: Stokhos_OrthogPolyApprox.hpp:150
Stokhos::OrthogPolyApprox::reset
void reset(const Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > &new_basis, ordinal_type sz=0)
Reset to a new basis.
Definition: Stokhos_OrthogPolyApproxImp.hpp:104
Stokhos::OrthogPolyApprox::basis
Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > basis() const
Return basis.
Definition: Stokhos_OrthogPolyApproxImp.hpp:96
Stokhos::ProductBasis
Abstract base class for multivariate orthogonal polynomials generated from tensor products of univari...
Definition: Stokhos_ProductBasis.hpp:31
Stokhos::OrthogPolyApprox::load
void load(value_type *v)
Load coefficients to an array of values.
Definition: Stokhos_OrthogPolyApproxImp.hpp:88
Stokhos::OrthogPolyApprox::mean
value_type mean() const
Compute mean of expansion.
Definition: Stokhos_OrthogPolyApproxImp.hpp:237
Stokhos::OrthogPolyApprox::const_reference
storage_type::const_reference const_reference
Definition: Stokhos_OrthogPolyApprox.hpp:33
Stokhos::OrthogPolyApprox::const_pointer
storage_type::const_pointer const_pointer
Definition: Stokhos_OrthogPolyApprox.hpp:35
Stokhos::OrthogPolyApprox
Class to store coefficients of a projection onto an orthogonal polynomial basis.
Definition: Stokhos_OrthogPolyApprox.hpp:29
val
expr val()
Stokhos::StandardStorage::init
void init(const_reference v)
Initialize values to a constant value.
Definition: Stokhos_StandardStorage.hpp:42
Stokhos::OrthogPolyApprox::print
std::ostream & print(std::ostream &os) const
Print approximation in basis.
Definition: Stokhos_OrthogPolyApproxImp.hpp:312
Stokhos::OrthogPolyApprox::size
ordinal_type size() const
Return size.
Definition: Stokhos_OrthogPolyApproxImp.hpp:126
Stokhos::StandardStorage::resize
void resize(const ordinal_type &sz)
Resize to new size (values are preserved)
Definition: Stokhos_StandardStorage.hpp:60
Stokhos::OrthogPolyApprox::order
const MultiIndex< ordinal_type > & order(ordinal_type term) const
Get orders for a given term.
Definition: Stokhos_OrthogPolyApproxImp.hpp:204
Stokhos::StandardStorage::size
ordinal_type size() const
Return size.
Definition: Stokhos_StandardStorage.hpp:63
Stokhos::OrthogPolyApprox::term
reference term(ordinal_type dimension, ordinal_type order)
Get coefficient term for given dimension and order.
Definition: Stokhos_OrthogPolyApproxImp.hpp:176
Generated on Fri Nov 22 2024 09:18:13 for Stokhos Package Browser (Single Doxygen Collection) by   doxygen 1.8.5