1 /* _______________________________________________________________________
2
3 DAKOTA: Design Analysis Kit for Optimization and Terascale Applications
4 Copyright 2014-2020 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
5 This software is distributed under the GNU Lesser General Public License.
6 For more information, see the README file in the top Dakota directory.
7 _______________________________________________________________________ */
8
9 //- Class: NonDSparseGrid
10 //- Description: Wrapper class for C++ code from Pecos/packages/VPISparseGrid
11 //- Owner: Mike Eldred
12 //- Revised by:
13 //- Version:
14
15 #ifndef NOND_SPARSE_GRID_H
16 #define NOND_SPARSE_GRID_H
17
18 #include "dakota_data_types.hpp"
19 #include "DataMethod.hpp"
20 #include "NonDIntegration.hpp"
21 #include "SparseGridDriver.hpp"
22
23 namespace Dakota {
24
25
26 /// Derived nondeterministic class that generates N-dimensional
27 /// Smolyak sparse grids for numerical evaluation of expectation
28 /// integrals over independent standard random variables.
29
30 /** This class is used by NonDPolynomialChaos and
31 NonDStochCollocation, but could also be used for general numerical
32 integration of moments. It employs 1-D Clenshaw-Curtis and Gaussian
33 quadrature rules within Smolyak sparse grids. */
34
35 class NonDSparseGrid: public NonDIntegration
36 {
37 public:
38
39 //
40 //- Heading: Constructors and destructor
41 //
42
43 // alternate constructor for instantiations "on the fly"
44 NonDSparseGrid(Model& model, unsigned short ssg_level,
45 const RealVector& dim_pref,
46 short exp_coeffs_soln_approach, short driver_mode,
47 short growth_rate = Pecos::MODERATE_RESTRICTED_GROWTH,
48 //short refine_type = Pecos::NO_REFINEMENT,
49 short refine_control = Pecos::NO_CONTROL,
50 bool track_uniq_prod_wts = true);
51
52 ~NonDSparseGrid(); ///< destructor
53
54 //
55 //- Heading: Virtual function redefinitions
56 //
57
58 /// update the sparse grid level (e.g., from a level sequence)
59 void sparse_grid_level(unsigned short ssg_level);
60
61 /// increment ssgDriver::ssgLevel
62 void increment_grid();
63 /// update ssgDriver::ssgAnisoLevelWts and increment ssgDriver::ssgLevel
64 /// based on specified anisotropic weighting
65 void increment_grid_weights(const RealVector& aniso_wts);
66 /// increment ssgDriver::ssgLevel based on existing anisotropic weighting
67 void increment_grid_weights();
68 /// decrement ssgDriver::ssgLevel
69 void decrement_grid();
70
71 void evaluate_grid_increment();
72 void push_grid_increment();
73 void pop_grid_increment();
74 void merge_grid_increment();
75
76 /// reset ssgDriver level and dimension preference back to
77 /// {ssgLevel,dimPref}Spec for the active key, following refinement
78 /// or sequence advancement
79 void reset();
80 /// blow away all data for all keys
81 void reset_all();
82
83 /// returns SparseGridDriver::active_multi_index()
84 const std::set<UShortArray>& active_multi_index() const;
85
86 /// invokes SparseGridDriver::print_smolyak_multi_index()
87 void print_smolyak_multi_index() const;
88
89 /// invokes SparseGridDriver::initialize_sets()
90 void initialize_sets();
91 /// invokes SparseGridDriver::update_reference()
92 void update_reference();
93 /// invokes SparseGridDriver::increment_smolyak_multi_index()
94 void increment_set(const UShortArray& set);
95 /// invokes SparseGridDriver::unique_trial_points()
96 int increment_size() const;
97 /// invokes SparseGridDriver::push_set()
98 void push_set();
99 /// invokes SparseGridDriver::compute_trial_grid()
100 void evaluate_set();
101 /// invokes SparseGridDriver::pop_set()
102 void decrement_set();
103 /// invokes SparseGridDriver::update_sets()
104 void update_sets(const UShortArray& set_star);
105 /// invokes SparseGridDriver::finalize_sets()
106 void finalize_sets(bool output_sets, bool converged_within_tol,bool reverted);
107
108 int num_samples() const;
109
110 protected:
111
112 //
113 //- Heading: Constructors and destructor
114 //
115
116 NonDSparseGrid(ProblemDescDB& problem_db, Model& model); ///< constructor
117
118 //
119 //- Heading: Virtual function redefinitions
120 //
121
122 void initialize_grid(const std::vector<Pecos::BasisPolynomial>& poly_basis);
123
124 void get_parameter_sets(Model& model);
125
126 //void check_variables(const Pecos::ShortArray& x_types);
127
128 void sampling_reset(int min_samples, bool all_data_flag, bool stats_flag);
129
130 //
131 //- Heading: Member functions
132 //
133
134 const RealVector& anisotropic_weights() const;
135
136 private:
137
138 //
139 //- Heading: Data
140 //
141
142 /// type of sparse grid driver: combined, incremental, hierarchical, ...
143 short ssgDriverType;
144 /// convenience pointer to the numIntDriver representation
145 std::shared_ptr<Pecos::SparseGridDriver> ssgDriver;
146
147 /// the user specification for the Smolyak sparse grid level, rendered
148 /// anisotropic via dimPrefSpec
149 unsigned short ssgLevelSpec;
150 /// value of ssgDriver->level() prior to increment_grid(), for restoration
151 /// in decrement_grid() since increment must induce a change in grid size
152 /// and this adaptive increment in not reversible
153 unsigned short ssgLevelPrev;
154 };
155
156
increment_grid_weights()157 inline void NonDSparseGrid::increment_grid_weights()
158 { increment_grid_weights(ssgDriver->anisotropic_weights()); }
159
160
sparse_grid_level(unsigned short ssg_level)161 inline void NonDSparseGrid::sparse_grid_level(unsigned short ssg_level)
162 { ssgLevelSpec = ssg_level; reset(); }
163
164
reset()165 inline void NonDSparseGrid::reset()
166 {
167 // reset the grid for the current active key to its original user spec,
168 // prior to any grid refinement
169 // > also invokes SparseGridDriver::clear_size() if change is induced
170 // > updates to other keys are managed by {assign,increment}_specification_
171 // sequence() in multilevel expansion methods
172 ssgDriver->level(ssgLevelSpec);
173 ssgDriver->dimension_preference(dimPrefSpec);
174
175 // Clear grid size (may vary with either dist param change or grid
176 // level/anisotropy) and clear dist param update trackers
177 ssgDriver->reset();
178
179 // This fn does not clear history, such as accumulated 1D pts/wts -->
180 // reset_all() or reset_1d_collocation_points_weights() should be used
181 // when distribution param updates invalidate this history
182 }
183
184
reset_all()185 inline void NonDSparseGrid::reset_all()
186 {
187 // This "nuclear option" is not currently used
188
189 ssgDriver->clear_keys();
190 ssgDriver->update_active_iterators();
191 reset();
192 }
193
194
active_multi_index() const195 inline const std::set<UShortArray>& NonDSparseGrid::active_multi_index() const
196 { return ssgDriver->active_multi_index(); }
197
198
print_smolyak_multi_index() const199 inline void NonDSparseGrid::print_smolyak_multi_index() const
200 { return ssgDriver->print_smolyak_multi_index(); }
201
202
initialize_sets()203 inline void NonDSparseGrid::initialize_sets()
204 { ssgDriver->initialize_sets(); }
205
206
update_reference()207 inline void NonDSparseGrid::update_reference()
208 { ssgDriver->update_reference(); }
209
210
increment_set(const UShortArray & set)211 inline void NonDSparseGrid::increment_set(const UShortArray& set)
212 { ssgDriver->increment_smolyak_multi_index(set); }
213
214
increment_size() const215 inline int NonDSparseGrid::increment_size() const
216 { return ssgDriver->unique_trial_points(); }
217
218
push_set()219 inline void NonDSparseGrid::push_set()
220 { ssgDriver->push_set(); }
221
222
evaluate_set()223 inline void NonDSparseGrid::evaluate_set()
224 {
225 ssgDriver->compute_trial_grid(allSamples);
226 evaluate_parameter_sets(iteratedModel, true, false);
227 ++numIntegrations;
228 }
229
230
decrement_set()231 inline void NonDSparseGrid::decrement_set()
232 { ssgDriver->pop_set(); }
233
234
update_sets(const UShortArray & set_star)235 inline void NonDSparseGrid::update_sets(const UShortArray& set_star)
236 { ssgDriver->update_sets(set_star); }
237
238
239 inline void NonDSparseGrid::
finalize_sets(bool output_sets,bool converged_within_tol,bool reverted)240 finalize_sets(bool output_sets, bool converged_within_tol, bool reverted)
241 { ssgDriver->finalize_sets(output_sets, converged_within_tol, reverted); }
242
243
evaluate_grid_increment()244 inline void NonDSparseGrid::evaluate_grid_increment()
245 {
246 ssgDriver->compute_increment(allSamples);
247 evaluate_parameter_sets(iteratedModel, true, false);
248 ++numIntegrations;
249 }
250
251
push_grid_increment()252 inline void NonDSparseGrid::push_grid_increment()
253 { ssgDriver->push_increment(); }
254
255
pop_grid_increment()256 inline void NonDSparseGrid::pop_grid_increment()
257 { ssgDriver->pop_increment(); }
258
259
merge_grid_increment()260 inline void NonDSparseGrid::merge_grid_increment()
261 { ssgDriver->merge_unique(); }
262
263
num_samples() const264 inline int NonDSparseGrid::num_samples() const
265 { return ssgDriver->grid_size(); }
266
267
anisotropic_weights() const268 inline const RealVector& NonDSparseGrid::anisotropic_weights() const
269 { return ssgDriver->anisotropic_weights(); }
270
271 } // namespace Dakota
272
273 #endif
274