42#ifndef STOKHOS_PRODUCT_BASIS_UTILS_HPP
43#define STOKHOS_PRODUCT_BASIS_UTILS_HPP
45#include "Teuchos_Array.hpp"
46#include "Teuchos_RCP.hpp"
47#include "Teuchos_SerialDenseVector.hpp"
48#include "Teuchos_TimeMonitor.hpp"
59 template <
typename ordinal_type>
60 ordinal_type
n_choose_k(
const ordinal_type& n,
const ordinal_type& k) {
69 ordinal_type l = std::min(n-k,k);
70 for (ordinal_type i=0; i<l; i++) {
78 template <
typename ordinal_t>
155 std::ostream&
print(std::ostream& os)
const {
158 os <<
index[i] <<
" ";
198 template <
typename ordinal_type>
200 const MultiIndex<ordinal_type>& m) {
214 template <
typename ordinal_t>
281 class Iterator :
public std::iterator<std::input_iterator_tag,
285 typedef std::iterator<std::input_iterator_tag,multiindex_type>
base_type;
377 template <
typename ordinal_t>
452 class Iterator :
public std::iterator<std::input_iterator_tag,
456 typedef std::iterator<std::input_iterator_tag,multiindex_type>
base_type;
557 template <
typename ordinal_t>
639 class Iterator :
public std::iterator<std::input_iterator_tag,
643 typedef std::iterator<std::input_iterator_tag,multiindex_type>
base_type;
716 template <
typename ordinal_t,
typename element_t>
747 const Teuchos::Array<element_type>&
getTerm()
const {
return term; }
753 operator Teuchos::ArrayView<element_type>() {
return term; }
756 operator Teuchos::ArrayView<const element_type>()
const {
return term; }
766 std::ostream&
print(std::ostream& os)
const {
769 os <<
term[i] <<
" ";
777 Teuchos::Array<element_type>
term;
781 template <
typename ordinal_type,
typename element_type>
784 const TensorProductElement<ordinal_type,element_type>& m) {
797 template <
typename term_type,
798 typename compare_type = std::less<typename term_type::element_type> >
810 bool operator()(
const term_type& a,
const term_type& b)
const {
812 while(i < a.dimension() && i < b.dimension() &&
equal(a[i],b[i])) i++;
816 if (i == a.dimension())
return i != b.dimension();
820 if (i == b.dimension())
return false;
823 return cmp(a[i],b[i]);
833 return !(
cmp(a,b) ||
cmp(b,a));
847 template <
typename term_type,
848 typename compare_type = std::less<typename term_type::element_type> >
859 bool operator()(
const term_type& a,
const term_type& b)
const {
863 while (i < a.dimension() && i < b.dimension() &&
equal(a_order,b_order)) {
868 return cmp(a_order,b_order);
878 return !(
cmp(a,b) ||
cmp(b,a));
901 template <
typename term_type>
912 bool operator()(
const term_type& a,
const term_type& b)
const {
919 if ( (x < y) && (x < (x^y)) ) {
934 template <
typename value_type>
945 bool operator() (
const value_type& a,
const value_type& b)
const {
957 template <
typename ordinal_type>
971 template <
typename ordinal_type>
999 template <
typename ordinal_type>
1003 ordinal_type idx = ordinal_type(0);
1004 ordinal_type p = ordinal_type(0);
1005 ordinal_type den = ordinal_type(1);
1006 for (ordinal_type d=1; d<=dim; ++d) {
1012 ordinal_type num = 1;
1013 for (ordinal_type i=p; i<p+d; ++i)
1030 template <
typename ordinal_type>
1033 ordinal_type max_order) {
1035 ordinal_type idx = ordinal_type(0);
1036 for (ordinal_type d=1; d<=dim; ++d) {
1037 ordinal_type p = index[d-1];
1039 ordinal_type offset = ordinal_type(0);
1040 for (ordinal_type pp=0; pp<p; ++pp)
1063 template <
typename index_set_type,
1064 typename growth_rule_type,
1065 typename basis_set_type,
1066 typename basis_map_type>
1069 const growth_rule_type& growth_rule,
1070 basis_set_type& basis_set,
1071 basis_map_type& basis_map) {
1073 typedef typename index_set_type::ordinal_type ordinal_type;
1074 typedef typename index_set_type::iterator index_iterator_type;
1075 typedef typename basis_set_type::iterator basis_iterator_type;
1076 typedef typename basis_set_type::key_type coeff_type;
1078 ordinal_type dim = index_set.dimension();
1081 index_iterator_type index_iterator = index_set.begin();
1082 index_iterator_type index_iterator_end = index_set.end();
1083 for (; index_iterator != index_iterator_end; ++index_iterator) {
1086 coeff_type coeff(dim);
1087 for (ordinal_type i=0; i<dim; i++)
1088 coeff[i] = growth_rule[i]((*index_iterator)[i]);
1091 basis_set[coeff] = ordinal_type(0);
1095 basis_map.resize(basis_set.size());
1096 ordinal_type idx = 0;
1097 basis_iterator_type basis_iterator = basis_set.begin();
1098 basis_iterator_type basis_iterator_end = basis_set.end();
1099 for (; basis_iterator != basis_iterator_end; ++basis_iterator) {
1100 basis_iterator->second = idx;
1101 basis_map[idx] = basis_iterator->first;
1110 template <
typename index_set_type,
1111 typename basis_set_type,
1112 typename basis_map_type>
1115 basis_set_type& basis_set,
1116 basis_map_type& basis_map) {
1117 typedef typename index_set_type::ordinal_type ordinal_type;
1118 ordinal_type dim = index_set.dimension();
1119 Teuchos::Array< IdentityGrowthRule<ordinal_type> > growth_rule(dim);
1123 template <
typename ordinal_type,
1124 typename value_type,
1125 typename basis_set_type,
1126 typename basis_map_type,
1127 typename coeff_predicate_type,
1128 typename k_coeff_predicate_type>
1129 static Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> >
1131 const Teuchos::Array< Teuchos::RCP<
const OneDOrthogPolyBasis<ordinal_type, value_type> > >& bases,
1132 const basis_set_type& basis_set,
1133 const basis_map_type& basis_map,
1134 const coeff_predicate_type& coeff_pred,
1135 const k_coeff_predicate_type& k_coeff_pred,
1136 const value_type sparse_tol = 1.0e-12)
1138#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
1139 TEUCHOS_FUNC_TIME_MONITOR(
"Stokhos: Total Triple-Product Tensor Time");
1141 typedef typename basis_map_type::value_type coeff_type;
1142 ordinal_type d = bases.size();
1154 Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> > Cijk =
1155 Teuchos::rcp(
new Sparse3Tensor<ordinal_type, value_type>);
1158 Teuchos::Array< Teuchos::RCP<Sparse3Tensor<ordinal_type,value_type> > > Cijk_1d(d);
1159 for (ordinal_type i=0; i<d; i++) {
1161 bases[i]->computeSparseTripleProductTensor(bases[i]->order()+1);
1167 typedef Sparse3Tensor<ordinal_type,value_type> Cijk_type;
1168 typedef typename Cijk_type::k_iterator k_iterator;
1169 typedef typename Cijk_type::kj_iterator kj_iterator;
1170 typedef typename Cijk_type::kji_iterator kji_iterator;
1171 Teuchos::Array<k_iterator> k_iterators(d, Cijk_1d[0]->k_begin());
1172 Teuchos::Array<kj_iterator > j_iterators(d, Cijk_1d[0]->j_begin(k_iterators[0]));
1173 Teuchos::Array<kji_iterator > i_iterators(d, Cijk_1d[0]->i_begin(j_iterators[0]));
1174 coeff_type terms_i(d), terms_j(d), terms_k(d);
1175 for (ordinal_type dim=0; dim<d; dim++) {
1176 k_iterators[dim] = Cijk_1d[dim]->k_begin();
1177 j_iterators[dim] = Cijk_1d[dim]->j_begin(k_iterators[dim]);
1178 i_iterators[dim] = Cijk_1d[dim]->i_begin(j_iterators[dim]);
1179 terms_i[dim] = index(i_iterators[dim]);
1180 terms_j[dim] = index(j_iterators[dim]);
1181 terms_k[dim] = index(k_iterators[dim]);
1187 bool valid_i = coeff_pred(terms_i);
1188 bool valid_j = coeff_pred(terms_j);
1189 bool valid_k = k_coeff_pred(terms_k);
1194 ordinal_type cnt = 0;
1198 if (valid_i && valid_j && valid_k) {
1200 typename basis_set_type::const_iterator k =
1201 basis_set.find(terms_k);
1205 typename basis_set_type::const_iterator i =
1206 basis_set.find(terms_i);
1210 typename basis_set_type::const_iterator
j =
1211 basis_set.find(terms_j);
1214 value_type c = value_type(1.0);
1215 value_type nrm = value_type(1.0);
1216 for (ordinal_type dim=0; dim<d; dim++) {
1217 c *= value(i_iterators[dim]);
1218 nrm *= bases[dim]->norm_squared(terms_i[dim]);
1220 if (std::abs(c/nrm) > sparse_tol)
1221 Cijk->add_term(I,J,K,c);
1225 ordinal_type cdim = 0;
1230 while (inc && cdim < d) {
1231 ordinal_type cur_dim = cdim;
1232 ++i_iterators[cdim];
1234 if (i_iterators[cdim] != Cijk_1d[cdim]->i_end(j_iterators[cdim])) {
1235 terms_i[cdim] = index(i_iterators[cdim]);
1236 valid_i = coeff_pred(terms_i);
1238 if (i_iterators[cdim] == Cijk_1d[cdim]->i_end(j_iterators[cdim]) ||
1240 ++j_iterators[cdim];
1242 if (j_iterators[cdim] != Cijk_1d[cdim]->j_end(k_iterators[cdim])) {
1243 terms_j[cdim] = index(j_iterators[cdim]);
1244 valid_j = coeff_pred(terms_j);
1246 if (j_iterators[cdim] == Cijk_1d[cdim]->j_end(k_iterators[cdim]) ||
1248 ++k_iterators[cdim];
1250 if (k_iterators[cdim] != Cijk_1d[cdim]->k_end()) {
1251 terms_k[cdim] = index(k_iterators[cdim]);
1252 valid_k = k_coeff_pred(terms_k);
1254 if (k_iterators[cdim] == Cijk_1d[cdim]->k_end() || !valid_k) {
1255 k_iterators[cdim] = Cijk_1d[cdim]->k_begin();
1257 terms_k[cur_dim] = index(k_iterators[cur_dim]);
1258 valid_k = k_coeff_pred(terms_k);
1262 j_iterators[cur_dim] =
1263 Cijk_1d[cur_dim]->j_begin(k_iterators[cur_dim]);
1264 terms_j[cur_dim] = index(j_iterators[cur_dim]);
1265 valid_j = coeff_pred(terms_j);
1269 i_iterators[cur_dim] =
1270 Cijk_1d[cur_dim]->i_begin(j_iterators[cur_dim]);
1271 terms_i[cur_dim] = index(i_iterators[cur_dim]);
1272 valid_i = coeff_pred(terms_i);
1277 if (!valid_i || !valid_j || !valid_k)
1287 Cijk->fillComplete();
1292 template <
typename ordinal_type>
1320 if (!more)
return false;
1330 ordinal_type& delta_j,
1331 ordinal_type& delta_k) {
1334 ordinal_type i0 =
i;
1335 ordinal_type j0 =
j;
1336 ordinal_type k0 =
k;
1339 while (more && zero) {
1380 if (
k+
j <
i ||
i+
j <
k ||
i+
k <
j)
return true;
1386 template <
typename ordinal_type,
1387 typename value_type,
1388 typename basis_set_type,
1389 typename basis_map_type,
1390 typename coeff_predicate_type,
1391 typename k_coeff_predicate_type>
1392 static Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> >
1394 const Teuchos::Array< Teuchos::RCP<
const OneDOrthogPolyBasis<ordinal_type, value_type> > >& bases,
1395 const basis_set_type& basis_set,
1396 const basis_map_type& basis_map,
1397 const coeff_predicate_type& coeff_pred,
1398 const k_coeff_predicate_type& k_coeff_pred,
1399 bool symmetric =
false,
1400 const value_type sparse_tol = 1.0e-12) {
1401#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
1402 TEUCHOS_FUNC_TIME_MONITOR(
"Stokhos: Total Triple-Product Tensor Time");
1404 typedef typename basis_map_type::value_type coeff_type;
1405 ordinal_type d = bases.size();
1417 Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> > Cijk =
1418 Teuchos::rcp(
new Sparse3Tensor<ordinal_type, value_type>);
1421 Teuchos::Array< Teuchos::RCP<Sparse3Tensor<ordinal_type,value_type> > > Cijk_1d(d);
1422 for (ordinal_type i=0; i<d; i++) {
1424 bases[i]->computeSparseTripleProductTensor(bases[i]->order()+1);
1428 typedef Cijk_1D_Iterator<ordinal_type> Cijk_Iterator;
1429 Teuchos::Array< Cijk_1D_Iterator<ordinal_type> > Cijk_1d_iterators(d);
1430 coeff_type terms_i(d), terms_j(d), terms_k(d);
1431 for (ordinal_type dim=0; dim<d; dim++) {
1432 Cijk_1d_iterators[dim] = Cijk_Iterator(bases[dim]->order(), symmetric);
1441 bool valid_i = coeff_pred(terms_i);
1442 bool valid_j = coeff_pred(terms_j);
1443 bool valid_k = k_coeff_pred(terms_k);
1444 bool stop = !valid_i || !valid_j || !valid_k;
1445 ordinal_type cnt = 0;
1448 typename basis_set_type::const_iterator k = basis_set.find(terms_k);
1449 typename basis_set_type::const_iterator i = basis_set.find(terms_i);
1450 typename basis_set_type::const_iterator
j = basis_set.find(terms_j);
1455 value_type c = value_type(1.0);
1456 for (ordinal_type dim=0; dim<d; dim++) {
1457 c *= Cijk_1d[dim]->getValue(terms_i[dim], terms_j[dim], terms_k[dim]);
1459 if (std::abs(c) > sparse_tol) {
1460 Cijk->add_term(I,J,K,c);
1470 ordinal_type cdim = 0;
1472 while (inc && cdim < d) {
1473 bool more = Cijk_1d_iterators[cdim].increment();
1474 terms_i[cdim] = Cijk_1d_iterators[cdim].i;
1475 terms_j[cdim] = Cijk_1d_iterators[cdim].j;
1476 terms_k[cdim] = Cijk_1d_iterators[cdim].k;
1481 valid_i = coeff_pred(terms_i);
1482 valid_j = coeff_pred(terms_j);
1483 valid_k = k_coeff_pred(terms_k);
1485 if (valid_i && valid_j && valid_k)
1496 Cijk->fillComplete();
1505 template <
typename ordinal_type,
typename value_type>
1521 Teuchos::Array< MultiIndex<ordinal_type> >& terms,
1522 Teuchos::Array<ordinal_type>& num_terms);
1536 Teuchos::Array< MultiIndex<ordinal_type> >& terms,
1537 Teuchos::Array<ordinal_type>& num_terms);
1545 const Teuchos::Array< MultiIndex<ordinal_type> >& terms,
1546 const Teuchos::Array<ordinal_type>& num_terms,
1547 ordinal_type max_p);
1553template <
typename ordinal_type,
typename value_type>
1559 Teuchos::Array<ordinal_type>& num_terms)
1561 Teuchos::Array<ordinal_type> basis_orders(d, p);
1562 return compute_terms(basis_orders, sz, terms, num_terms);
1565template <
typename ordinal_type,
typename value_type>
1568compute_terms(
const Teuchos::Array<ordinal_type>& basis_orders,
1571 Teuchos::Array<ordinal_type>& num_terms)
1609 for (ordinal_type i=0; i<d; i++) {
1610 if (basis_orders[i] > p)
1611 p = basis_orders[i];
1615 Teuchos::Array< Teuchos::Array< MultiIndex<ordinal_type> > > terms_order(p+1);
1621 terms_order[0].resize(1);
1627 Teuchos::Array<ordinal_type> cnt(d), cnt_next(d);
1628 MultiIndex<ordinal_type> term(d);
1629 for (ordinal_type
j=0;
j<d;
j++) {
1630 if (basis_orders[
j] >= 1)
1639 for (ordinal_type k=1; k<=p; k++) {
1641 num_terms[k] = num_terms[k-1];
1647 for (ordinal_type
j=0;
j<d;
j++) {
1650 for (ordinal_type i=0; i<cnt[
j]; i++) {
1651 if (terms_order[k-1][prev+i][
j] < basis_orders[
j]) {
1652 term = terms_order[k-1][prev+i];
1654 terms_order[k].push_back(term);
1657 for (ordinal_type l=0; l<=
j; l++)
1664 prev += cnt[
j] - cnt[
j+1];
1669 for (ordinal_type
j=0;
j<d;
j++) {
1670 cnt[
j] = cnt_next[
j];
1676 num_terms[p+1] = sz;
1681 for (ordinal_type k=0; k<=p; k++) {
1683 for (ordinal_type
j=0;
j<num_k;
j++)
1684 terms[i++] = terms_order[k][
j];
1690template <
typename ordinal_type,
typename value_type>
1695 const Teuchos::Array<ordinal_type>& num_terms,
1704 ordinal_type ord = 0;
1705 for (ordinal_type i=0; i<d; i++)
1707 TEUCHOS_TEST_FOR_EXCEPTION(ord < 0 || ord > max_p, std::logic_error,
1708 "Stokhos::CompletePolynomialBasis::compute_index(): " <<
1709 "Term has invalid order " << ord);
1716 k = num_terms[ord-1];
1717 ordinal_type k_max=num_terms[ord];
1719 while (k < k_max && !found) {
1720 bool found_term =
true;
1721 for (ordinal_type
j=0;
j<d;
j++) {
1722 found_term = found_term && (term[
j] == terms[k][
j]);
1729 TEUCHOS_TEST_FOR_EXCEPTION(k >= k_max && !found, std::logic_error,
1730 "Stokhos::CompletePolynomialBasis::compute_index(): " <<
1731 "Could not find specified term.");
Iterator class for iterating over elements of the index set.
const multiindex_type * const_pointer
base_type::difference_type difference_type
base_type::value_type value_type
multiindex_type index
Current value of iterator.
ordinal_type max_order
Maximum order of iterator.
ordinal_type dim
Dimension.
Iterator(ordinal_type max_order_, const multiindex_type &component_max_order_, const multiindex_type &index_)
Constructor.
multiindex_type component_max_order
Maximum order for each component.
const_pointer operator->() const
Dereference.
base_type::pointer pointer
base_type::reference reference
base_type::iterator_category iterator_category
Iterator & operator++()
Prefix increment, i.e., ++iterator.
std::iterator< std::input_iterator_tag, multiindex_type > base_type
const_reference operator*() const
Dereference.
bool operator==(const Iterator &it) const
Compare equality of iterators.
const multiindex_type & const_reference
Teuchos::Array< ordinal_type > orders
Maximum orders for each term to determine how to increment.
bool operator!=(const Iterator &it) const
Compare inequality of iterators.
Iterator & operator++(int)
Postfix increment, i.e., iterator++.
An anisotropic total order index set.
multiindex_type max_orders() const
Return maximum order for each dimension.
const_iterator end() const
Return iterator for end of the index set.
ordinal_type lower_order
Lower order of index set.
multiindex_type lower
Component-wise lower bounds.
ordinal_type dim
Dimension.
MultiIndex< ordinal_type > multiindex_type
ordinal_type dimension() const
Return dimension.
AnisotropicTotalOrderIndexSet(ordinal_type upper_order_, const multiindex_type &lower_, const multiindex_type &upper_)
Constructor.
AnisotropicTotalOrderIndexSet(ordinal_type upper_order_, const multiindex_type &upper_)
Constructor.
multiindex_type upper
Component-wise upper bounds.
const_iterator begin() const
Return iterator for first element in the set.
ordinal_type upper_order
Upper order of index set.
Utilities for indexing a multi-variate complete polynomial basis.
static ordinal_type compute_terms(ordinal_type p, ordinal_type d, ordinal_type &sz, Teuchos::Array< MultiIndex< ordinal_type > > &terms, Teuchos::Array< ordinal_type > &num_terms)
Compute the 2-D array of basis terms which maps a basis index into the orders for each basis dimensio...
static ordinal_type compute_terms(const Teuchos::Array< ordinal_type > &basis_orders, ordinal_type &sz, Teuchos::Array< MultiIndex< ordinal_type > > &terms, Teuchos::Array< ordinal_type > &num_terms)
Compute the 2-D array of basis terms which maps a basis index into the orders for each basis dimensio...
static ordinal_type compute_index(const MultiIndex< ordinal_type > &term, const Teuchos::Array< MultiIndex< ordinal_type > > &terms, const Teuchos::Array< ordinal_type > &num_terms, ordinal_type max_p)
Compute basis index given the orders for each basis dimension.
A functor for comparing floating-point numbers to some tolerance.
FloatingPointLess(const value_type &tol_=1.0e-12)
Constructor.
bool operator()(const value_type &a, const value_type &b) const
Compare if a < b.
~FloatingPointLess()
Destructor.
A comparison functor implementing a strict weak ordering based lexographic ordering.
bool operator()(const term_type &a, const term_type &b) const
Determine if a is less than b.
LexographicLess(const compare_type &cmp_=compare_type())
Constructor.
compare_type cmp
Element comparison functor.
bool equal(const element_type &a, const element_type &b) const
Determine if two elements a and b are equal.
term_type product_element_type
term_type::ordinal_type ordinal_type
term_type::element_type element_type
A comparison functor implementing a strict weak ordering based Morton Z-ordering.
term_type product_element_type
term_type::element_type element_type
MortonZLess()
Constructor.
term_type::ordinal_type ordinal_type
bool operator()(const term_type &a, const term_type &b) const
A multidimensional index.
void resize(ordinal_type d, ordinal_type v=ordinal_type(0))
Resize.
Teuchos::Array< ordinal_type > index
index terms
Teuchos::Array< element_type > & getTerm()
Term access.
ordinal_type dimension() const
Dimension.
const Teuchos::Array< element_type > & getTerm() const
Term access.
bool operator==(const MultiIndex &idx) const
Compare equality.
MultiIndex & termWiseMin(const MultiIndex &idx)
Replace multiindex with min of this and other multiindex.
bool operator!=(const MultiIndex &idx) const
Compare equality.
MultiIndex(ordinal_type dim, ordinal_type v=ordinal_type(0))
Constructor.
ordinal_type order() const
Compute total order of index.
std::ostream & print(std::ostream &os) const
Print multiindex.
MultiIndex & termWiseMin(const ordinal_type idx)
Replace multiindex with min of this and given value.
ordinal_type size() const
Size.
void init(ordinal_type v)
Initialize.
MultiIndex & termWiseMax(const MultiIndex &idx)
Replace multiindex with max of this and other multiindex.
bool termWiseLEQ(const MultiIndex &idx) const
Compare term-wise less-than or equal-to.
MultiIndex & termWiseMax(const ordinal_type idx)
Replace multiindex with max of this and given value.
const ordinal_type & operator[](ordinal_type i) const
Term access.
Utilities for indexing a multi-variate complete polynomial basis.
static void buildProductBasis(const index_set_type &index_set, basis_set_type &basis_set, basis_map_type &basis_map)
Generate a product basis from an index set.
static Teuchos::RCP< Stokhos::Sparse3Tensor< ordinal_type, value_type > > computeTripleProductTensorNew(const Teuchos::Array< Teuchos::RCP< const OneDOrthogPolyBasis< ordinal_type, value_type > > > &bases, const basis_set_type &basis_set, const basis_map_type &basis_map, const coeff_predicate_type &coeff_pred, const k_coeff_predicate_type &k_coeff_pred, bool symmetric=false, const value_type sparse_tol=1.0e-12)
static void buildProductBasis(const index_set_type &index_set, const growth_rule_type &growth_rule, basis_set_type &basis_set, basis_map_type &basis_map)
Generate a product basis from an index set.
static Teuchos::RCP< Stokhos::Sparse3Tensor< ordinal_type, value_type > > computeTripleProductTensor(const Teuchos::Array< Teuchos::RCP< const OneDOrthogPolyBasis< ordinal_type, value_type > > > &bases, const basis_set_type &basis_set, const basis_map_type &basis_map, const coeff_predicate_type &coeff_pred, const k_coeff_predicate_type &k_coeff_pred, const value_type sparse_tol=1.0e-12)
Container storing a term in a generalized tensor product.
element_type order() const
Compute total order of tensor product element.
std::ostream & print(std::ostream &os) const
Print multiindex.
Teuchos::Array< element_type > & getTerm()
Term access.
TensorProductElement()
Default constructor.
ordinal_type dimension() const
Return dimension.
ordinal_type size() const
Return dimension.
~TensorProductElement()
Destructor.
const Teuchos::Array< element_type > & getTerm() const
Term access.
const element_type & operator[](ordinal_type i) const
Term access.
TensorProductElement(ordinal_type dim, const element_type &val=element_type(0))
Constructor.
Teuchos::Array< element_type > term
Array storing term elements.
Iterator class for iterating over elements of the index set.
const multiindex_type & const_reference
multiindex_type index
Current value of iterator.
base_type::value_type value_type
const_pointer operator->() const
Dereference.
base_type::difference_type difference_type
std::iterator< std::input_iterator_tag, multiindex_type > base_type
Iterator & operator++(int)
Postfix increment, i.e., iterator++.
ordinal_type dim
Dimension.
Iterator(const multiindex_type &upper_, const multiindex_type &index_)
Constructor.
multiindex_type upper
Upper bound of iterator.
const_reference operator*() const
Dereference.
bool operator!=(const Iterator &it) const
Compare inequality of iterators.
base_type::iterator_category iterator_category
base_type::pointer pointer
Iterator & operator++()
Prefix increment, i.e., ++iterator.
bool operator==(const Iterator &it) const
Compare equality of iterators.
base_type::reference reference
const multiindex_type * const_pointer
A tensor product index set.
multiindex_type upper
Upper bound of index set.
TensorProductIndexSet(const multiindex_type &lower_, const multiindex_type &upper_)
Constructor.
ordinal_type dimension() const
Return dimension.
TensorProductIndexSet(ordinal_type dim_, ordinal_type upper_)
Constructor.
TensorProductIndexSet(ordinal_type dim_, ordinal_type lower_, ordinal_type upper_)
Constructor.
const_iterator begin() const
Return iterator for first element in the set.
TensorProductIndexSet(const multiindex_type &upper_)
Constructor.
const_iterator end() const
Return iterator for end of the index set.
MultiIndex< ordinal_type > multiindex_type
ordinal_type dim
Dimension.
multiindex_type lower
Lower bound of index set.
multiindex_type max_orders() const
Return maximum order for each dimension.
Iterator class for iterating over elements of the index set.
bool operator==(const Iterator &it) const
Compare equality of iterators.
std::iterator< std::input_iterator_tag, multiindex_type > base_type
Iterator(ordinal_type max_order_, const multiindex_type &index_)
Constructor.
Iterator & operator++()
Prefix increment, i.e., ++iterator.
base_type::value_type value_type
const_pointer operator->() const
Dereference.
multiindex_type index
Current value of iterator.
const multiindex_type * const_pointer
base_type::difference_type difference_type
Teuchos::Array< ordinal_type > orders
Maximum orders for each term to determine how to increment.
bool operator!=(const Iterator &it) const
Compare inequality of iterators.
const_reference operator*() const
Dereference.
ordinal_type dim
Dimension.
base_type::pointer pointer
const multiindex_type & const_reference
base_type::reference reference
base_type::iterator_category iterator_category
ordinal_type max_order
Maximum order of iterator.
Iterator & operator++(int)
Postfix increment, i.e., iterator++.
An isotropic total order index set.
const_iterator end() const
Return iterator for end of the index set.
const_iterator begin() const
Return iterator for first element in the set.
MultiIndex< ordinal_type > multiindex_type
ordinal_type upper
Upper order of index set.
multiindex_type max_orders() const
Return maximum order for each dimension.
ordinal_type dim
Dimension.
ordinal_type lower
Lower order of index set.
ordinal_type dimension() const
Return dimension.
TotalOrderIndexSet(ordinal_type dim_, ordinal_type lower_, ordinal_type upper_)
Constructor.
TotalOrderIndexSet(ordinal_type dim_, ordinal_type upper_)
Constructor.
A comparison functor implementing a strict weak ordering based total-order ordering,...
term_type::element_type element_type
term_type::ordinal_type ordinal_type
TotalOrderLess(const compare_type &cmp_=compare_type())
Constructor.
term_type product_element_type
bool operator()(const term_type &a, const term_type &b) const
bool equal(const element_type &a, const element_type &b) const
Determine if two elements a and b are equal.
compare_type cmp
Element comparison functor.
Top-level namespace for Stokhos classes and functions.
ordinal_type totalOrderMapping(const Stokhos::MultiIndex< ordinal_type > &index)
ordinal_type lexicographicMapping(const Stokhos::MultiIndex< ordinal_type > &index, ordinal_type max_order)
std::ostream & operator<<(std::ostream &os, const ProductContainer< coeff_type > &vec)
ordinal_type n_choose_k(const ordinal_type &n, const ordinal_type &k)
Compute bionomial coefficient (n ; k) = n!/( k! (n-k)! )
bool increment(ordinal_type &delta_i, ordinal_type &delta_j, ordinal_type &delta_k)
Cijk_1D_Iterator(ordinal_type p_i, ordinal_type p_j, ordinal_type p_k, bool sym)
Cijk_1D_Iterator(ordinal_type p=0, bool sym=false)
Predicate functor for building sparse triple products.
MultiIndex< ordinal_type > coeff_type
bool operator()(const coeff_type &term) const
TensorProductPredicate(const coeff_type &orders_)
Predicate functor for building sparse triple products based on total order.
MultiIndex< ordinal_type > coeff_type
bool operator()(const coeff_type &term) const
TotalOrderPredicate(ordinal_type max_order_, const coeff_type &orders_)