queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
TensorProductQuadrature.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // QUESO - a library to support the Quantification of Uncertainty
5 // for Estimation, Simulation and Optimization
6 //
7 // Copyright (C) 2008-2017 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 #include <queso/TensorProductQuadrature.h>
26 #include <queso/asserts.h>
27 #include <queso/1DQuadrature.h>
28 #include <queso/MultiDimensionalIndexing.h>
29 #include <queso/VectorSpace.h>
30 #include <queso/GslVector.h>
31 #include <queso/GslMatrix.h>
32 
33 namespace QUESO
34 {
35  template<class V,class M>
37  const std::vector<QUESO::SharedPtr<Base1DQuadrature>::Type> & q_rules )
38  : MultiDQuadratureBase<V,M>(domain)
39  {
40  const unsigned int dim = domain.vectorSpace().dimGlobal();
41 
42  queso_require_equal_to_msg(dim, q_rules.size(), "Mismatched quadrature rule size and vector space dimension!");
43 
44  // We figure the number of q_points in each dimension and the total
45  unsigned int total_n_q_points = 1;
46  std::vector<unsigned int> n_q_points(dim);
47  for( unsigned int i = 0; i < dim; i++ )
48  {
49  n_q_points[i] = q_rules[i]->positions().size();
50  total_n_q_points *= q_rules[i]->positions().size();
51  }
52 
53  this->m_positions.resize(total_n_q_points,SharedPtr<GslVector>::Type());
54  this->m_weights.resize(total_n_q_points);
55 
56  // The positions are just a tensor product of the positions for each of the 1-D rules.
57  // Thus, we leverage the MultiDimensionalIndexing to move from "global" quadrature point
58  // number to the index in each of the dimensions. Each of these dimensional indices tells
59  // us which component of the corresponding 1-D rule to extract.
60  // The final quadrature weight will be the product of each of the corresponding weights
61  // for the 1-D rule, with the correspondence dictated by the "local" index in that dimension
62  // (as with the position).
63  std::vector<unsigned int> indices(dim);
64  for( unsigned int q = 0; q < total_n_q_points; q++ )
65  {
66  MultiDimensionalIndexing::globalToCoord( q, n_q_points, indices );
67 
68  // Make sure indices size isn't something we don't expect after population
69  queso_assert_equal_to(indices.size(),dim);
70 
71  // Create new vector from our vector space
72  // We are taking ownership of this pointer.
73  typename QUESO::SharedPtr<V>::Type domain_vec(domain.vectorSpace().newVector());
74 
75  // Now populate the weights and positions
76  double weight = 1.0;
77  for( unsigned int i = 0; i < dim; i++ )
78  {
79  unsigned int idx = indices[i];
80  (*domain_vec)[i] = q_rules[i]->positions()[idx];
81  weight *= q_rules[i]->weights()[idx];
82  }
83 
84  // m_positions is taking ownership of the V*
85  this->m_positions[q] = domain_vec;
86  this->m_weights[q] = weight;
87  }
88  }
89 
90  // Instantiate
92 
93 } // end namespace QUESO
std::shared_ptr< T > Type
Definition: SharedPtr.h:69
A templated class for handling subsets.
Definition: VectorSubset.h:46
Numerical quadrature using a tensor product of Base1DQuadrature rules.
static void globalToCoord(unsigned int global, const std::vector< unsigned int > &n_points, std::vector< unsigned int > &coord_indices)
Inverse of coordToGlobal map.
std::vector< double > m_weights
int dim
Definition: ann_test.cpp:472
MonteCarloSGOptions::MonteCarloSGOptions(const BaseEnvironment &env, const char *prefix, const McOptionsValues &alternativeOptionsValues queso_require_equal_to_msg)(m_env.optionsInputFileName(), std::string(""), std::string("this constructor is incompatible with the existence of an options input file"))
Base class for multi-dimensional quadrature rules.
TensorProductQuadrature(const VectorSubset< V, M > &domain, const std::vector< QUESO::SharedPtr< Base1DQuadrature >::Type > &q_rules)
std::vector< typename QUESO::SharedPtr< V >::Type > m_positions
Locations of quadrature points.
const VectorSpace< V, M > & vectorSpace() const
Vector space to which this set belongs to. See template specialization.
Definition: VectorSubset.C:68

Generated on Tue Jun 5 2018 19:48:55 for queso-0.57.1 by  doxygen 1.8.5