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:       NestedModel
10 //- Description: Implementation code for the NestedModel class
11 //- Owner:       Mike Eldred
12 //- Checked by:
13 
14 #include "NestedModel.hpp"
15 #include "ProblemDescDB.hpp"
16 #include "MarginalsCorrDistribution.hpp"
17 #include "dakota_system_defs.hpp"
18 #include "pecos_global_defs.hpp"
19 #include "EvaluationStore.hpp"
20 
21 static const char rcsId[]="@(#) $Id: NestedModel.cpp 7024 2010-10-16 01:24:42Z mseldre $";
22 
23 //#define DEBUG
24 
25 
26 namespace Dakota {
27 
NestedModel(ProblemDescDB & problem_db)28 NestedModel::NestedModel(ProblemDescDB& problem_db):
29   Model(BaseConstructor(), problem_db),
30   nestedModelEvalCntr(0), firstUpdate(true), outerMIPLIndex(0),
31   subIteratorSched(parallelLib,
32 		   true, // peer 1 must assign jobs to peers 2-n
33 		   problem_db.get_int("model.nested.iterator_servers"),
34 		   problem_db.get_int("model.nested.processors_per_iterator"),
35 		   problem_db.get_short("model.nested.iterator_scheduling")),
36   subMethodPointer(problem_db.get_string("model.nested.sub_method_pointer")),
37   subIteratorJobCntr(0),
38   optInterfacePointer(problem_db.get_string("model.interface_pointer"))
39 {
40   ignoreBounds = problem_db.get_bool("responses.ignore_bounds");
41   centralHess  = problem_db.get_bool("responses.central_hess");
42 
43   // Retrieve the variable mapping inputs
44   const StringArray& primary_var_mapping
45     = problem_db.get_sa("model.nested.primary_variable_mapping");
46   const StringArray& secondary_var_mapping
47     = problem_db.get_sa("model.nested.secondary_variable_mapping");
48 
49   // NestedModel may set the DB list nodes, but is required to restore them
50   // to their previous setting in order to remove the need to continuously
51   // reset at the environment level (which would be wasteful since the type
52   // of derived model may not be known at the environment level).
53   size_t method_index = problem_db.get_db_method_node(); // for restoration
54   size_t model_index  = problem_db.get_db_model_node();  // for restoration
55 
56   // interface for non-nested data is optional
57   if (optInterfacePointer.empty())
58     numOptInterfPrimary = numOptInterfIneqCon = numOptInterfEqCon = 0;
59   else {
60     const String& oi_resp_ptr
61       = problem_db.get_string("model.optional_interface_responses_pointer");
62     bool oi_resp_ptr_defined = !oi_resp_ptr.empty();
63     if (oi_resp_ptr_defined)
64       problem_db.set_db_responses_node(oi_resp_ptr);
65     // JAS: We need to work a little harder here to make sure that the
66     // optional interface responses have a gradient and hessian spec that is
67     // compatible with the gradient and hessian spec for the nested model's
68     // responses. Currently, if the nested model's responses specify analytic
69     // or mixed gradients or hessians but the optional interface responses just
70     // have numerical or no, Dakota dies with a segfault and no error message.
71     // An even better solution might be to wrap the optional interface in a
72     // SimulationModel, which would allow gradient/hessian requests to be
73     // satisfied however the user wants, and would make it easier for us to
74     // honor a request for scaling.
75 
76     optInterfGradientType = problem_db.get_string("responses.gradient_type");
77     optInterfHessianType = problem_db.get_string("responses.hessian_type");
78     optInterfGradIdAnalytic
79       = problem_db.get_is("responses.gradients.mixed.id_analytic");
80     optInterfHessIdAnalytic
81       = problem_db.get_is("responses.hessians.mixed.id_analytic");
82 
83     numOptInterfIneqCon
84       = problem_db.get_sizet("responses.num_nonlinear_inequality_constraints");
85     numOptInterfEqCon
86       = problem_db.get_sizet("responses.num_nonlinear_equality_constraints");
87     optInterfaceResponse
88       = problem_db.get_response(SIMULATION_RESPONSE, currentVariables);
89     optionalInterface = problem_db.get_interface();
90     size_t num_fns = optInterfaceResponse.num_functions();
91     numOptInterfPrimary = num_fns - numOptInterfIneqCon - numOptInterfEqCon;
92 
93     if (oi_resp_ptr_defined) {
94       // Echo a warning if there is a user specification of constraint bounds/
95       // targets that will be superceded by top level constraint bounds/targets.
96       const RealVector& interf_ineq_l_bnds
97 	= problem_db.get_rv("responses.nonlinear_inequality_lower_bounds");
98       const RealVector& interf_ineq_u_bnds
99 	= problem_db.get_rv("responses.nonlinear_inequality_upper_bounds");
100       const RealVector& interf_eq_targets
101 	= problem_db.get_rv("responses.nonlinear_equality_targets");
102       bool warning_flag = false;
103       size_t i;
104       Real dbl_inf = std::numeric_limits<Real>::infinity();
105       for (i=0; i<numOptInterfIneqCon; ++i)
106 	if ( interf_ineq_l_bnds[i] > -dbl_inf || interf_ineq_u_bnds[i] != 0. )
107 	  warning_flag = true;
108       for (i=0; i<numOptInterfEqCon; ++i)
109 	if ( interf_eq_targets[i] != 0. )
110 	  warning_flag = true;
111       if (warning_flag)
112 	Cerr << "Warning: nonlinear constraint bounds and targets in nested "
113 	     << "model optional interfaces\n         are superceded by "
114 	     << "composite response constraint bounds and targets."
115 	     << std::endl;
116     }
117 
118     // db_responses restore not needed since set_db_list_nodes below will reset
119   }
120 
121   problem_db.set_db_list_nodes(subMethodPointer); // even if empty
122 
123   subModel = problem_db.get_model();
124   //check_submodel_compatibility(subModel); // sanity checks performed below
125   // if outer level output is verbose/debug, request fine-grained evaluation
126   // reporting for purposes of the final output summary.  This allows verbose
127   // final summaries without verbose output on every sub-iterator completion.
128   if (outputLevel > NORMAL_OUTPUT)
129     subModel.fine_grained_evaluation_counters();
130 
131   problem_db.set_db_method_node(method_index); // restore method only
132   problem_db.set_db_model_nodes(model_index);  // restore all model nodes
133 
134   // Perform error checks on variable mapping inputs and convert from
135   // strings to indices for efficiency at run time.
136   size_t i, num_var_map_1 = primary_var_mapping.size(),
137     num_var_map_2 = secondary_var_mapping.size(),
138     num_curr_cv   = currentVariables.cv(),
139     num_curr_div  = currentVariables.div(),
140     num_curr_dsv  = currentVariables.dsv(),
141     num_curr_drv  = currentVariables.drv(),
142     num_curr_vars = num_curr_cv + num_curr_div + num_curr_dsv + num_curr_drv;
143   // Error checks: maps can be empty strings, but must be present to assure
144   // correct association.
145   if ( ( num_var_map_1 && num_var_map_1 != num_curr_vars ) ||
146        ( num_var_map_2 && num_var_map_2 != num_var_map_1 ) ) {
147     Cerr << "\nError: length of variable mapping specification(s) does not "
148 	 << "match number of active variables." << std::endl;
149     abort_handler(MODEL_ERROR);
150   }
151   // active are sized based on totals due to different mapping options
152   active1ACVarMapIndices.resize(num_curr_vars);
153   active1ADIVarMapIndices.resize(num_curr_vars);
154   active1ADSVarMapIndices.resize(num_curr_vars);
155   active1ADRVarMapIndices.resize(num_curr_vars);
156   extraCVarsData.resize(num_curr_cv);
157   extraDIVarsData.resize(num_curr_div);
158   extraDSVarsData.resize(num_curr_dsv);
159   extraDRVarsData.resize(num_curr_drv);
160   if (num_var_map_2) {
161     active2ACVarMapTargets.resize(num_curr_vars);
162     active2ADIVarMapTargets.resize(num_curr_vars);
163     active2ADSVarMapTargets.resize(num_curr_vars);
164     active2ADRVarMapTargets.resize(num_curr_vars);
165   }
166   short inactive_sm_view = EMPTY_VIEW;
167   UShortMultiArrayConstView curr_c_types
168     = currentVariables.continuous_variable_types();
169   UShortMultiArrayConstView curr_di_types
170     = currentVariables.discrete_int_variable_types();
171   UShortMultiArrayConstView curr_ds_types
172     = currentVariables.discrete_string_variable_types();
173   UShortMultiArrayConstView curr_dr_types
174     = currentVariables.discrete_real_variable_types();
175   UShortMultiArrayConstView submodel_a_c_types
176     = subModel.all_continuous_variable_types();
177   UShortMultiArrayConstView submodel_a_di_types
178     = subModel.all_discrete_int_variable_types();
179   UShortMultiArrayConstView submodel_a_ds_types
180     = subModel.all_discrete_string_variable_types();
181   UShortMultiArrayConstView submodel_a_dr_types
182     = subModel.all_discrete_real_variable_types();
183   size_t curr_i, sm_acv_cntr = 0, sm_adiv_cntr = 0, sm_adsv_cntr = 0,
184     sm_adrv_cntr = 0, sm_acv_avail = 0, sm_adiv_avail = 0, sm_adsv_avail = 0,
185     sm_adrv_avail = 0;
186   unsigned short prev_c_type = USHRT_MAX, prev_di_type = USHRT_MAX,
187     prev_ds_type = USHRT_MAX, prev_dr_type = USHRT_MAX;
188   String empty_str;
189 
190   const SharedVariablesData& svd = currentVariables.shared_data();
191 
192   // Map ACTIVE CONTINUOUS VARIABLES from currentVariables
193   for (i=0; i<num_curr_cv; ++i) {
194     curr_i = svd.cv_index_to_active_index(i);
195     const String& map1
196       = (num_var_map_1) ? primary_var_mapping[curr_i] : empty_str;
197     if (map1.empty()) {
198       // for default mappings between consistent types, propagate bounds/labels
199       extraCVarsData.set(i);
200       // default mapping: inactive subModel vars = active currentVariables
201       update_inactive_view(currentVariables.view().first, inactive_sm_view);
202       // Can't use label matching, since subModel labels may not be updated
203       // until runtime.  index() returns the _first_ instance of the type.
204       unsigned short curr_c_type = curr_c_types[i];
205       size_t sm_acv_offset = find_index(submodel_a_c_types, curr_c_type);
206       if (sm_acv_offset == _NPOS) {
207 	Cerr << "\nError: active continuous variable type '" << curr_c_type
208 	     << "' could not be matched within all sub-model continuous "
209 	     << "variable types." << std::endl;
210 	abort_handler(MODEL_ERROR);
211       }
212       // For multiple types, sm_acv_cntr must be reset to 0 at the type change
213       if (curr_c_type != prev_c_type) {
214 	sm_acv_cntr  = 0;
215 	sm_acv_avail = std::count(submodel_a_c_types.begin(),
216 				  submodel_a_c_types.end(), curr_c_type);
217       }
218       active1ACVarMapIndices[curr_i]  = sm_acv_offset + sm_acv_cntr++;
219       active1ADIVarMapIndices[curr_i] = active1ADSVarMapIndices[curr_i]
220 	= active1ADRVarMapIndices[curr_i] = _NPOS;
221       prev_c_type = curr_c_type;
222       if (sm_acv_cntr > sm_acv_avail) {
223 	Cerr << "\nError: default insertions of type '" << curr_c_type
224 	     << "' exceed sub-model allocation." << std::endl;
225 	abort_handler(MODEL_ERROR);
226       }
227       if (num_var_map_2)
228 	active2ACVarMapTargets[curr_i] = active2ADIVarMapTargets[curr_i]
229 	  = active2ADSVarMapTargets[curr_i] = active2ADRVarMapTargets[curr_i]
230 	  = Pecos::NO_TARGET;
231     }
232     else {
233       extraCVarsData.reset(i); // not a default mapping based on types
234       const String& map2
235 	= (num_var_map_2) ? secondary_var_mapping[curr_i] : empty_str;
236       resolve_real_variable_mapping(map1, map2, curr_i, inactive_sm_view);
237     }
238   }
239 
240   // Map ACTIVE DISCRETE INTEGER VARIABLES from currentVariables
241   for (i=0; i<num_curr_div; ++i) {
242     curr_i = svd.div_index_to_active_index(i);
243     const String& map1
244       = (num_var_map_1) ? primary_var_mapping[curr_i] : empty_str;
245     if (map1.empty()) {
246       // for default mappings between consistent types, propagate bounds/labels
247       extraDIVarsData.set(i);
248       // default mapping: inactive subModel vars = active currentVariables
249       update_inactive_view(currentVariables.view().first, inactive_sm_view);
250       // Can't use label matching, since subModel labels may not be updated
251       // until runtime.  index() returns the _first_ instance of the type.
252       unsigned short curr_di_type = curr_di_types[i];
253       size_t sm_adiv_offset = find_index(submodel_a_di_types, curr_di_type);
254       if (sm_adiv_offset == _NPOS) {
255 	Cerr << "\nError: active discrete integer variable type '"
256 	     << curr_di_type << "' could not be matched within all sub-model "
257 	     << "discrete integer variable types." << std::endl;
258 	abort_handler(MODEL_ERROR);
259       }
260       // For multiple types, sm_adiv_cntr must be reset to 0 at the type change
261       if (curr_di_type != prev_di_type) {
262 	sm_adiv_cntr  = 0;
263 	sm_adiv_avail = std::count(submodel_a_di_types.begin(),
264 				   submodel_a_di_types.end(), curr_di_type);
265       }
266       active1ADIVarMapIndices[curr_i] = sm_adiv_offset + sm_adiv_cntr++;
267       active1ACVarMapIndices[curr_i]  = active1ADSVarMapIndices[curr_i]
268 	= active1ADRVarMapIndices[curr_i] = _NPOS;
269       prev_di_type = curr_di_type;
270       if (sm_adiv_cntr > sm_adiv_avail) {
271 	Cerr << "\nError: default insertions of type '" << curr_di_type
272 	     << "' exceed sub-model allocation." << std::endl;
273 	abort_handler(MODEL_ERROR);
274       }
275       if (num_var_map_2)
276 	active2ACVarMapTargets[curr_i] = active2ADIVarMapTargets[curr_i]
277 	  = active2ADSVarMapTargets[curr_i] = active2ADRVarMapTargets[curr_i]
278 	  = Pecos::NO_TARGET;
279     }
280     else {
281       extraDIVarsData.reset(i); // not a default mapping based on types
282       const String& map2
283 	= (num_var_map_2) ? secondary_var_mapping[curr_i] : empty_str;
284       resolve_integer_variable_mapping(map1, map2, curr_i, inactive_sm_view);
285     }
286   }
287 
288   // Map ACTIVE DISCRETE STRING VARIABLES from currentVariables
289   for (i=0; i<num_curr_dsv; ++i) {
290     curr_i = svd.dsv_index_to_active_index(i);
291     const String& map1
292       = (num_var_map_1) ? primary_var_mapping[curr_i] : empty_str;
293     if (map1.empty()) {
294       // for default mappings between consistent types, propagate bounds/labels
295       extraDSVarsData.set(i);
296       // default mapping: inactive subModel vars = active currentVariables
297       update_inactive_view(currentVariables.view().first, inactive_sm_view);
298       // Can't use label matching, since subModel labels may not be updated
299       // until runtime.  index() returns the _first_ instance of the type.
300       unsigned short curr_ds_type = curr_ds_types[i];
301       size_t sm_adsv_offset = find_index(submodel_a_ds_types, curr_ds_type);
302       if (sm_adsv_offset == _NPOS) {
303 	Cerr << "\nError: active discrete string variable type '"
304 	     << curr_ds_type << "' could not be matched within all sub-model "
305 	     << "discrete string variable types." << std::endl;
306 	abort_handler(MODEL_ERROR);
307       }
308       // For multiple types, sm_adsv_cntr must be reset to 0 at the type change
309       if (curr_ds_type != prev_ds_type) {
310 	sm_adsv_cntr  = 0;
311 	sm_adsv_avail = std::count(submodel_a_ds_types.begin(),
312 				   submodel_a_ds_types.end(), curr_ds_type);
313       }
314       active1ADSVarMapIndices[curr_i] = sm_adsv_offset + sm_adsv_cntr++;
315       active1ACVarMapIndices[curr_i]  = active1ADIVarMapIndices[curr_i]
316 	= active1ADRVarMapIndices[curr_i] = _NPOS;
317       prev_ds_type = curr_ds_type;
318       if (sm_adsv_cntr > sm_adsv_avail) {
319 	Cerr << "\nError: default insertions of type '" << curr_ds_type
320 	     << "' exceed sub-model allocation." << std::endl;
321 	abort_handler(MODEL_ERROR);
322       }
323       if (num_var_map_2)
324 	active2ACVarMapTargets[curr_i] = active2ADIVarMapTargets[curr_i]
325 	  = active2ADSVarMapTargets[curr_i] = active2ADRVarMapTargets[curr_i]
326 	  = Pecos::NO_TARGET;
327     }
328     else {
329       extraDSVarsData.reset(i); // not a default mapping based on types
330       const String& map2
331 	= (num_var_map_2) ? secondary_var_mapping[curr_i] : empty_str;
332       resolve_string_variable_mapping(map1, map2, curr_i, inactive_sm_view);
333     }
334   }
335 
336   // Map ACTIVE DISCRETE REAL VARIABLES from currentVariables
337   for (i=0; i<num_curr_drv; ++i) {
338     curr_i = svd.drv_index_to_active_index(i);
339     const String& map1
340       = (num_var_map_1) ? primary_var_mapping[curr_i] : empty_str;
341     if (map1.empty()) {
342       // for default mappings between consistent types, propagate bounds/labels
343       extraDRVarsData.set(i);
344       // default mapping: inactive subModel vars = active currentVariables
345       update_inactive_view(currentVariables.view().first, inactive_sm_view);
346       // Can't use label matching, since subModel labels may not be updated
347       // until runtime.  index() returns the _first_ instance of the type.
348       unsigned short curr_dr_type = curr_dr_types[i];
349       size_t sm_adrv_offset = find_index(submodel_a_dr_types, curr_dr_type);
350       if (sm_adrv_offset == _NPOS) {
351 	Cerr << "\nError: active discrete real variable type '" << curr_dr_type
352 	     << "' could not be matched within all sub-model discrete real "
353 	     << "variable types." << std::endl;
354 	abort_handler(MODEL_ERROR);
355       }
356       // For multiple types, sm_adrv_cntr must be reset to 0 at the type change
357       if (curr_dr_type != prev_dr_type) {
358 	sm_adrv_cntr  = 0;
359 	sm_adrv_avail = std::count(submodel_a_dr_types.begin(),
360 				   submodel_a_dr_types.end(), curr_dr_type);
361       }
362       active1ADRVarMapIndices[curr_i] = sm_adrv_offset + sm_adrv_cntr++;
363       active1ACVarMapIndices[curr_i]  = active1ADIVarMapIndices[curr_i]
364 	= active1ADSVarMapIndices[curr_i] = _NPOS;
365       prev_dr_type = curr_dr_type;
366       if (sm_adrv_cntr > sm_adrv_avail) {
367 	Cerr << "\nError: default insertions of type '" << curr_dr_type
368 	     << "' exceed sub-model allocation." << std::endl;
369 	abort_handler(MODEL_ERROR);
370       }
371       if (num_var_map_2)
372 	active2ACVarMapTargets[curr_i] = active2ADIVarMapTargets[curr_i]
373 	  = active2ADSVarMapTargets[curr_i] = active2ADRVarMapTargets[curr_i]
374 	  = Pecos::NO_TARGET;
375     }
376     else {
377       extraDRVarsData.reset(i); // not a default mapping based on types
378       const String& map2
379 	= (num_var_map_2) ? secondary_var_mapping[curr_i] : empty_str;
380       resolve_real_variable_mapping(map1, map2, curr_i, inactive_sm_view);
381     }
382   }
383 
384   size_t num_curr_ccv = currentVariables.acv() - num_curr_cv,//active complement
385     num_curr_cdiv = currentVariables.adiv() - num_curr_div, // active complement
386     num_curr_cdsv = currentVariables.adsv() - num_curr_dsv, // active complement
387     num_curr_cdrv = currentVariables.adrv() - num_curr_drv; // active complement
388   // complement can be sized based on corresponding currentVariables sizes
389   // due to restriction to default mappings
390   if (num_curr_ccv)
391     complement1ACVarMapIndices.resize(num_curr_ccv);
392   if (num_curr_cdiv)
393     complement1ADIVarMapIndices.resize(num_curr_cdiv);
394   if (num_curr_cdsv)
395     complement1ADSVarMapIndices.resize(num_curr_cdsv);
396   if (num_curr_cdrv)
397     complement1ADRVarMapIndices.resize(num_curr_cdrv);
398   UShortMultiArrayConstView curr_ac_types
399     = currentVariables.all_continuous_variable_types();
400   UShortMultiArrayConstView curr_adi_types
401     = currentVariables.all_discrete_int_variable_types();
402   UShortMultiArrayConstView curr_ads_types
403     = currentVariables.all_discrete_string_variable_types();
404   UShortMultiArrayConstView curr_adr_types
405     = currentVariables.all_discrete_real_variable_types();
406   unsigned short prev_ac_type = USHRT_MAX, prev_adi_type = USHRT_MAX,
407     prev_ads_type = USHRT_MAX, prev_adr_type = USHRT_MAX;
408 
409   // Map COMPLEMENT CONTINUOUS VARIABLES from currentVariables
410   for (i=0; i<num_curr_ccv; ++i) {
411     curr_i = svd.ccv_index_to_acv_index(i);
412     // Can't use label matching, since subModel labels may not be updated
413     // until runtime.  index() returns the _first_ instance of the type.
414     unsigned short curr_ac_type = curr_ac_types[curr_i];
415     size_t sm_acv_offset = find_index(submodel_a_c_types, curr_ac_type);
416     if (sm_acv_offset == _NPOS) {
417       Cerr << "\nError: complement continuous variable type '" << curr_ac_type
418 	   << "' could not be matched within all sub-model continuous "
419 	   << "variable types." << std::endl;
420       abort_handler(MODEL_ERROR);
421     }
422     // For multiple types, sm_acv_cntr must be reset to 0 at the type change
423     if (curr_ac_type != prev_ac_type) {
424       sm_acv_cntr  = 0;
425       sm_acv_avail = std::count(submodel_a_c_types.begin(),
426 				submodel_a_c_types.end(), curr_ac_type);
427     }
428     complement1ACVarMapIndices[i] = sm_acv_offset + sm_acv_cntr++;
429     prev_ac_type = curr_ac_type;
430     if (sm_acv_cntr > sm_acv_avail) {
431       Cerr << "\nError: default insertions of type '" << curr_ac_type
432 	   << "' exceed sub-model allocation." << std::endl;
433       abort_handler(MODEL_ERROR);
434     }
435   }
436 
437   // Map COMPLEMENT DISCRETE INTEGER VARIABLES from currentVariables
438   for (i=0; i<num_curr_cdiv; ++i) {
439     curr_i = svd.cdiv_index_to_adiv_index(i);
440     // Can't use label matching, since subModel labels may not be updated
441     // until runtime.  index() returns the _first_ instance of the type.
442     unsigned short curr_adi_type = curr_adi_types[curr_i];
443     size_t sm_adiv_offset = find_index(submodel_a_di_types, curr_adi_type);
444     if (sm_adiv_offset == _NPOS) {
445       Cerr << "\nError: complement discrete integer variable type '"
446 	   << curr_adi_type << "' could not be matched within all sub-model "
447 	   << "discrete integer variable types." << std::endl;
448       abort_handler(MODEL_ERROR);
449     }
450     // For multiple types, sm_adiv_cntr must be reset to 0 at the type change
451     if (curr_adi_type != prev_adi_type) {
452       sm_adiv_cntr  = 0;
453       sm_adiv_avail = std::count(submodel_a_di_types.begin(),
454 				 submodel_a_di_types.end(), curr_adi_type);
455     }
456     complement1ADIVarMapIndices[i] = sm_adiv_offset + sm_adiv_cntr++;
457     prev_adi_type = curr_adi_type;
458     if (sm_adiv_cntr > sm_adiv_avail) {
459       Cerr << "\nError: default insertions of type '" << curr_adi_type
460 	   << "' exceed sub-model allocation." << std::endl;
461       abort_handler(MODEL_ERROR);
462     }
463   }
464 
465   // Map COMPLEMENT DISCRETE STRING VARIABLES from currentVariables
466   for (i=0; i<num_curr_cdsv; ++i) {
467     curr_i = svd.cdsv_index_to_adsv_index(i);
468     // Can't use label matching, since subModel labels may not be updated
469     // until runtime.  index() returns the _first_ instance of the type.
470     unsigned short curr_ads_type = curr_ads_types[curr_i];
471     size_t sm_adsv_offset = find_index(submodel_a_ds_types, curr_ads_type);
472     if (sm_adsv_offset == _NPOS) {
473       Cerr << "\nError: complement discrete string variable type '"
474 	   << curr_ads_type << "' could not be matched within all sub-model "
475 	   << "discrete string variable types." << std::endl;
476       abort_handler(MODEL_ERROR);
477     }
478     // For multiple types, sm_adsv_cntr must be reset to 0 at the type change
479     if (curr_ads_type != prev_ads_type) {
480       sm_adsv_cntr  = 0;
481       sm_adsv_avail = std::count(submodel_a_ds_types.begin(),
482 				 submodel_a_ds_types.end(), curr_ads_type);
483     }
484     complement1ADSVarMapIndices[i] = sm_adsv_offset + sm_adsv_cntr++;
485     prev_ads_type = curr_ads_type;
486     if (sm_adsv_cntr > sm_adsv_avail) {
487       Cerr << "\nError: default insertions of type '" << curr_ads_type
488 	   << "' exceed sub-model allocation." << std::endl;
489       abort_handler(MODEL_ERROR);
490     }
491   }
492 
493   // Map COMPLEMENT DISCRETE REAL VARIABLES from currentVariables
494   for (i=0; i<num_curr_cdrv; ++i) {
495     curr_i = svd.cdrv_index_to_adrv_index(i);
496     // Can't use label matching, since subModel labels may not be updated
497     // until runtime.  index() returns the _first_ instance of the type.
498     unsigned short curr_adr_type = curr_adr_types[curr_i];
499     size_t sm_adrv_offset = find_index(submodel_a_dr_types, curr_adr_type);
500     if (sm_adrv_offset == _NPOS) {
501       Cerr << "\nError: complement discrete real variable type '"
502 	   << curr_adr_type << "' could not be matched within all sub-model "
503 	   << "discrete real variable types." << std::endl;
504       abort_handler(MODEL_ERROR);
505     }
506     // For multiple types, sm_adrv_cntr must be reset to 0 at the type change
507     if (curr_adr_type != prev_adr_type) {
508       sm_adrv_cntr  = 0;
509       sm_adrv_avail = std::count(submodel_a_dr_types.begin(),
510 				 submodel_a_dr_types.end(), curr_adr_type);
511     }
512     complement1ADRVarMapIndices[i] = sm_adrv_offset + sm_adrv_cntr++;
513     prev_adr_type = curr_adr_type;
514     if (sm_adrv_cntr > sm_adrv_avail) {
515       Cerr << "\nError: default insertions of type '" << curr_adr_type
516 	   << "' exceed sub-model allocation." << std::endl;
517       abort_handler(MODEL_ERROR);
518     }
519   }
520 
521 #ifdef DEBUG
522   Cout << "\nactive primary variable mapping indices\nACV:\n"
523        << active1ACVarMapIndices << "ADIV:\n" << active1ADIVarMapIndices
524        << "ADSV:\n" << active1ADSVarMapIndices
525        << "ADRV:\n" << active1ADRVarMapIndices
526        << "\nactive secondary variable mapping targets:\nACV:\n"
527        << active2ACVarMapTargets << "ADIV:\n" << active2ADIVarMapTargets
528        << "ADSV:\n" << active2ADSVarMapTargets
529        << "ADRV:\n" << active2ADRVarMapTargets
530        << "\ncomplement primary variable mapping indices\nACV:\n"
531        << complement1ACVarMapIndices << "ADIV:\n" << complement1ADIVarMapIndices
532        << "ADSV:\n" << complement1ADSVarMapIndices
533        << "ADRV:\n" << complement1ADRVarMapIndices << '\n';
534 #endif // DEBUG
535 
536   // subModel view updating must be performed before subIterator instantiation
537   // since any model recastings will pick up the inactive view (and inactive
538   // view differences cause problems with recursion updating).
539   if (inactive_sm_view != EMPTY_VIEW)
540     subModel.inactive_view(inactive_sm_view); // recurse
541 }
542 
543 
declare_sources()544 void NestedModel::declare_sources() {
545   evaluationsDB.declare_source(modelId, modelType, subIterator.method_id(),
546       "iterator");
547   if(!optionalInterface.is_null())
548     evaluationsDB.declare_source(modelId, modelType, optionalInterface.interface_id(),
549       "interface");
550 }
551 
552 /** Asynchronous flags need to be initialized for the subModel.  In
553     addition, max_eval_concurrency is the outer level iterator
554     concurrency, not the subIterator concurrency that subModel will
555     see, and recomputing the message_lengths on the subModel is
556     probably not a bad idea either.  Therefore, recompute everything
557     on subModel using init_communicators(). */
558 void NestedModel::
derived_init_communicators(ParLevLIter pl_iter,int max_eval_concurrency,bool recurse_flag)559 derived_init_communicators(ParLevLIter pl_iter, int max_eval_concurrency,
560 			   bool recurse_flag)
561 {
562   // initialize optionalInterface for parallel operations
563   if (!optInterfacePointer.empty()) {
564     // allow recursion to progress - don't store/set/restore
565     parallelLib.parallel_configuration_iterator(modelPCIter);
566     optionalInterface.init_communicators(messageLengths, max_eval_concurrency);
567   }
568 
569   if (!recurse_flag)
570     return;
571 
572   // Due to Model::init_communicators(), we know that pl_iter is the lowest
573   // level within modelPCIter->miPLIters, from which we may further subdivide.
574 
575   // initializations for subIteratorSched:
576   // > incoming max_eval_concurrency is for concurrent execs of NestedModel,
577   //   which must be distinguished from eval concurrency within subIterator.
578   // > circular dependency between concurrent iterator partitioning and
579   //   subIterator instantiation is managed as described in
580   //   ConcurrentMetaIterator::init_communicators().
581   // > as for constructors, we recursively set and restore DB list nodes
582   //   (initiated from the restored starting point following construction).
583   size_t method_index = probDescDB.get_db_method_node(),
584          model_index  = probDescDB.get_db_model_node();  // for restoration
585   probDescDB.set_db_list_nodes(subMethodPointer);
586 
587   // > init_eval_concurrency instantiates subIterator on previous pl ranks
588   subIteratorSched.update(modelPCIter);
589   // > define min and max processors per iterator
590   IntIntPair ppi_pr
591     = subIteratorSched.configure(probDescDB, subIterator, subModel);
592   // > passed in max_eval_concurrency is the outer nested model concurrency
593   subIteratorSched.partition(max_eval_concurrency, ppi_pr);
594   // > now augment prev subIterator instantiations for additional mi_pl ranks
595   //   (new mi_pl is used via miPLIndex update in partition())
596   // > idle server is managed here; a dedicated master processor is managed
597   //   within IteratorScheduler::init_iterator().
598   if (subIteratorSched.iteratorServerId <= subIteratorSched.numIteratorServers)
599     subIteratorSched.init_iterator(probDescDB, subIterator, subModel);
600 
601   // > restore all DB nodes
602   probDescDB.set_db_method_node(method_index);
603   probDescDB.set_db_model_nodes(model_index);
604 
605   // > now that subIterator is constructed, perform downstream updates
606   if (!subIterator.is_null()) {
607     init_sub_iterator(); // follow DB restore: extracts data from nested spec
608     if (subIteratorSched.messagePass) {
609       // msg lengths: vars from this model, set & final results from subIterator
610       MPIPackBuffer buff; int eval_id = 0;
611       const Response& si_resp = subIterator.response_results();
612       buff << currentVariables << si_resp.active_set() << eval_id;
613       int params_buff_len = buff.size(); buff.reset();
614       buff << si_resp;
615       subIteratorSched.iterator_message_lengths(params_buff_len, buff.size());
616     }
617   }
618 }
619 
620 
621 void NestedModel::
derived_set_communicators(ParLevLIter pl_iter,int max_eval_concurrency,bool recurse_flag)622 derived_set_communicators(ParLevLIter pl_iter, int max_eval_concurrency,
623 			  bool recurse_flag)
624 {
625   // Outer context:
626   outerMIPLIndex = modelPCIter->mi_parallel_level_index(pl_iter);
627 
628   if (!optInterfacePointer.empty()) {
629     // allow recursion to progress - don't store/set/restore
630     parallelLib.parallel_configuration_iterator(modelPCIter);
631     optionalInterface.set_communicators(messageLengths, max_eval_concurrency);
632     // initial setting for asynchEvalFlag & evaluationCapacity based on
633     // optInterface (may be updated below)
634     set_ie_asynchronous_mode(max_eval_concurrency);
635   }
636   if (recurse_flag) {
637     // Inner context: set comms for subIterator
638     // > pl_iter is incoming context prior to subIterator partitioning
639     // > mi_pl_index reflects the miPL depth after subIterator partitioning
640     size_t mi_pl_index = outerMIPLIndex + 1;
641     subIteratorSched.update(modelPCIter, mi_pl_index);
642     if (subIteratorSched.iteratorServerId <=
643 	subIteratorSched.numIteratorServers) {
644       ParLevLIter si_pl_iter
645 	= modelPCIter->mi_parallel_level_iterator(mi_pl_index);
646       subIteratorSched.set_iterator(subIterator, si_pl_iter);
647     }
648 
649     // update asynchEvalFlag & evaluationCapacity based on subIteratorSched
650     if (subIteratorSched.messagePass)
651       asynchEvalFlag = true;
652     if (subIteratorSched.numIteratorServers > evaluationCapacity)
653       evaluationCapacity = subIteratorSched.numIteratorServers;
654   }
655 }
656 
657 
658 void NestedModel::
derived_free_communicators(ParLevLIter pl_iter,int max_eval_concurrency,bool recurse_flag)659 derived_free_communicators(ParLevLIter pl_iter, int max_eval_concurrency,
660 			   bool recurse_flag)
661 {
662   /*
663   if (!optInterfacePointer.empty()) {
664     // allow recursion to progress - don't store/set/restore
665     parallelLib.parallel_configuration_iterator(modelPCIter);
666     optionalInterface.free_communicators();
667   }
668   */
669   if (recurse_flag) {
670     // finalize comms for subIterator
671     // > pl_iter is incoming context prior to subIterator partitioning
672     // > mi_pl_index reflects the miPL depth after subIterator partitioning
673     size_t mi_pl_index = modelPCIter->mi_parallel_level_index(pl_iter) + 1;
674     subIteratorSched.update(modelPCIter, mi_pl_index);
675     if (subIteratorSched.iteratorServerId <=
676 	subIteratorSched.numIteratorServers) {
677       ParLevLIter si_pl_iter
678 	= modelPCIter->mi_parallel_level_iterator(mi_pl_index);
679       subIteratorSched.free_iterator(subIterator, si_pl_iter);
680     }
681     subIteratorSched.free_iterator_parallelism();
682   }
683 }
684 
685 
init_sub_iterator()686 void NestedModel::init_sub_iterator()
687 {
688   // activates additional data processing options and manages output
689   subIterator.sub_iterator_flag(true);
690   // passes data through to subModel via recursions layered by subIterator;
691   // VarMaps are used by ProbabilityTransformModel, which generally appears
692   // within these Model recursions (so can't call subModel here without
693   // requiring inclusion in the update_from_sub_model() chain)
694   subIterator.nested_variable_mappings(active1ACVarMapIndices,
695     active1ADIVarMapIndices, active1ADSVarMapIndices, active1ADRVarMapIndices,
696     active2ACVarMapTargets,  active2ADIVarMapTargets, active2ADSVarMapTargets,
697     active2ADRVarMapTargets);
698 
699   // Back out the number of eq/ineq constraints within secondaryRespCoeffs
700   // (subIterator constraints) from the total number of equality/inequality
701   // constraints and the number of interface equality/inequality constraints.
702   size_t num_mapped_ineq_con
703     = probDescDB.get_sizet("responses.num_nonlinear_inequality_constraints"),
704     num_mapped_eq_con
705     = probDescDB.get_sizet("responses.num_nonlinear_equality_constraints");
706   numSubIterMappedIneqCon = num_mapped_ineq_con - numOptInterfIneqCon;
707   numSubIterMappedEqCon   = num_mapped_eq_con   - numOptInterfEqCon;
708 
709   // For auto-generated identity map, we don't support an optional
710   // interface and require that the total number of NestedModel
711   // responses (pri+sec) equal the number of sub-iterator results
712   // functions.
713   //
714   // Rationale: Too complex to get right and communicate to
715   // user. Secondary responses from an optional interface augment any
716   // mapped secondary responses. Nested model primary responses come
717   // from an overlay of optional interface and mapped responses, with
718   // no restrictions on which is larger or smaller.
719   size_t num_mapped_total = currentResponse.num_functions(),
720     num_mapped_sec = num_mapped_ineq_con + num_mapped_eq_con,
721     num_mapped_pri = num_mapped_total - num_mapped_sec;
722   numSubIterFns = subIterator.response_results().num_functions();
723 
724   identityRespMap = probDescDB.get_bool("model.nested.identity_resp_map");
725   const RealVector& primary_resp_coeffs
726     = probDescDB.get_rv("model.nested.primary_response_mapping");
727   const RealVector& secondary_resp_coeffs
728     = probDescDB.get_rv("model.nested.secondary_response_mapping");
729 
730   if (identityRespMap) {
731     bool found_error = false;
732     if (!optInterfacePointer.empty()) {
733       Cerr << "\nError: identity_response_mapping not supported in conjunction"
734        << " with optional_interface_pointer; use explicit primary/secondary_"
735        << "response_mapping instead.\n";
736       found_error = true;
737     }
738     if (!primary_resp_coeffs.empty() || !secondary_resp_coeffs.empty()) {
739       Cerr << "\nError: Neither primary_response_mapping nor secondary_"
740         << "response_mapping may be specified in conjunction with identity_"
741         << "response_mapping.\n";
742       found_error = true;
743     }
744     if (num_mapped_total != numSubIterFns) {
745       Cerr << "\nError: For identity_response_mapping, number of nested model "
746         << "responses (primary + secondary functions) must equal the number of "
747         << "sub-method final results. Specified nested model has "
748         << num_mapped_total << " functions, while there are " << numSubIterFns
749         << " sub-method results.\n";
750       if (outputLevel >= VERBOSE_OUTPUT)
751 	Cerr << "Info: Sub-method returns these results:\n"
752 	     << subIterator.response_results().function_labels() << "\n";
753       else
754 	Cerr << "Info: Re-run with 'output verbose' to list the sub-method "
755 	     << "results.\n";
756       found_error = true;
757     }
758     if (found_error)
759       abort_handler(-1);
760 
761     if (outputLevel >= VERBOSE_OUTPUT)
762       Cout << "Info: NestedModel using identity response mapping." << std::endl;
763     subIterMappedPri = num_mapped_pri;
764     subIterMappedSec = num_mapped_sec;
765 
766     // BMA TODO: When using identity map, propagate labels when not
767     // present (there's no way to detect that the user didn't give
768     // them as defaults get generated at parse time)
769   }
770   else if (primary_resp_coeffs.empty() && secondary_resp_coeffs.empty()) {
771     Cerr << "\nError: no mappings provided for sub-iterator functions in "
772 	 << "NestedModel initialization." << std::endl;
773     abort_handler(MODEL_ERROR);
774   }
775   // Convert vectors to matrices using the number of subIterator response
776   // results (e.g., the number of UQ statistics) as the number of columns
777   // (rows are open ended since the number of subModel primary fns may
778   // be different from the total number of primary fns and subModel's
779   // constraints are a subset of the total constraints).
780   if (!primary_resp_coeffs.empty()) {
781     // this error would also be caught within copy_data(), but by checking here,
782     // we can provide a more helpful error message.
783     if (primary_resp_coeffs.length() % numSubIterFns) {
784       Cerr << "\nError: number of entries in primary_response_mapping ("
785 	   << primary_resp_coeffs.length() << ") not evenly divisible"
786 	   << "\n       by number of sub-iterator final results functions ("
787 	   << numSubIterFns << ") in NestedModel initialization." << std::endl;
788       Cerr << "\nInfo: The primary_response_mapping must have between 1 and "
789 	   << num_mapped_pri
790 	   << " (number of nested model primary response functions) row(s).\n"
791 	   << "It must have " << numSubIterFns
792            << " columns corresponding to the sub-method final results.\n";
793       if (outputLevel >= VERBOSE_OUTPUT)
794 	Cerr << "Info: Sub-method returns these results:\n"
795 	     << subIterator.response_results().function_labels() << "\n";
796       else
797 	Cerr << "Info: Re-run with 'output verbose' to list the sub-method "
798 	     << "results.\n";
799       abort_handler(MODEL_ERROR);
800     }
801     copy_data(primary_resp_coeffs, primaryRespCoeffs, 0, (int)numSubIterFns);
802     subIterMappedPri = primaryRespCoeffs.numRows();
803   }
804   if (!secondary_resp_coeffs.empty()) {
805     // this error would also be caught within copy_data(), but by checking here,
806     // we can provide a more helpful error message.
807     if (secondary_resp_coeffs.length() % numSubIterFns) {
808       Cerr << "\nError: number of entries in secondary_response_mapping ("
809 	   << secondary_resp_coeffs.length() << ") not evenly divisible"
810 	   << "\n       by number of sub-iterator final results functions ("
811 	   << numSubIterFns << ") in NestedModel initialization." << std::endl;
812       Cerr << "\nInfo: The secondary_response_mapping must have "
813 	   << numSubIterMappedIneqCon + numSubIterMappedEqCon
814 	   << " (number of nested model secondary response functions, less any "
815            << "optional interface secondary response functions) row(s).\n"
816 	   << "It must have " << numSubIterFns
817            << " columns corresponding to the sub-method final results.\n";
818       if (outputLevel >= VERBOSE_OUTPUT)
819 	Cerr << "Info: Sub-method returns these results:\n"
820 	     << subIterator.response_results().function_labels() << "\n";
821       else
822 	Cerr << "Info: Re-run with 'output verbose' to list the sub-method "
823 	     << "results.\n";
824       abort_handler(MODEL_ERROR);
825     }
826     copy_data(secondary_resp_coeffs, secondaryRespCoeffs, 0,(int)numSubIterFns);
827     subIterMappedSec = secondaryRespCoeffs.numRows();
828   }
829 }
830 
831 
832 void NestedModel::
resolve_map1(const String & map1,size_t & ac_index1,size_t & adi_index1,size_t & ads_index1,size_t & adr_index1,size_t curr_index,short & inactive_sm_view)833 resolve_map1(const String& map1, size_t& ac_index1, size_t& adi_index1,
834 	     size_t& ads_index1, size_t& adr_index1, size_t curr_index,
835 	     short& inactive_sm_view)
836 {
837   adi_index1 = ads_index1 = adr_index1 = _NPOS;
838   ac_index1 = find_index(subModel.all_continuous_variable_labels(), map1);
839   if (ac_index1 == _NPOS) {
840     adi_index1 = find_index(subModel.all_discrete_int_variable_labels(), map1);
841     if (adi_index1 == _NPOS) {
842       ads_index1
843 	= find_index(subModel.all_discrete_string_variable_labels(), map1);
844       if (ads_index1 == _NPOS) {
845 	adr_index1
846 	  = find_index(subModel.all_discrete_real_variable_labels(), map1);
847         if (adr_index1 == _NPOS) {
848 	  Cerr << "\nError: primary mapping " << map1 << " could not be "
849 	       << "matched within any sub-model variable labels." << std::endl;
850 	  abort_handler(MODEL_ERROR);
851 	}
852 	else if (find_index(subModel.discrete_real_variable_labels(), map1)
853 		 == _NPOS) // inactive DRV target
854 	  update_inactive_view(
855 	    subModel.all_discrete_real_variable_types()[adr_index1],
856 	    inactive_sm_view);
857       }
858       else if (find_index(subModel.discrete_string_variable_labels(), map1)
859 	       == _NPOS) // inactive DSV target
860 	update_inactive_view(
861 	  subModel.all_discrete_string_variable_types()[ads_index1],
862 	  inactive_sm_view);
863     }
864     else if (find_index(subModel.discrete_int_variable_labels(), map1)
865 	     == _NPOS) // inactive DIV target
866       update_inactive_view(
867 	subModel.all_discrete_int_variable_types()[adi_index1],
868 	inactive_sm_view);
869   }
870   else if (find_index(subModel.continuous_variable_labels(), map1)
871 	   == _NPOS) // inactive CV target
872     update_inactive_view(subModel.all_continuous_variable_types()[ac_index1],
873 			 inactive_sm_view);
874 
875   active1ACVarMapIndices[curr_index]  =  ac_index1;
876   active1ADIVarMapIndices[curr_index] = adi_index1;
877   active1ADSVarMapIndices[curr_index] = ads_index1;
878   active1ADRVarMapIndices[curr_index] = adr_index1;
879 }
880 
881 
882 void NestedModel::
resolve_real_variable_mapping(const String & map1,const String & map2,size_t curr_index,short & inactive_sm_view)883 resolve_real_variable_mapping(const String& map1, const String& map2,
884 			      size_t curr_index, short& inactive_sm_view)
885 {
886   size_t ac_index1, adi_index1, ads_index1, adr_index1;
887   resolve_map1(map1, ac_index1, adi_index1, ads_index1, adr_index1,
888 	       curr_index, inactive_sm_view);
889 
890   if (!active2ACVarMapTargets.empty()) { // indicates num_var_map_2
891     if (map2.empty())
892       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
893 	= active2ADSVarMapTargets[curr_index]
894 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
895     else if (ac_index1 != _NPOS) {
896       unsigned short type = subModel.all_continuous_variable_types()[ac_index1];
897       if (type == CONTINUOUS_DESIGN || type == CONTINUOUS_STATE) {
898 	if (map2 == "lower_bound")
899 	  active2ACVarMapTargets[curr_index] = Pecos::CR_LWR_BND;
900 	else if (map2 == "upper_bound")
901 	  active2ACVarMapTargets[curr_index] = Pecos::CR_UPR_BND;
902 	else {
903 	  Cerr << "\nError: " << map2 << " mapping not supported for "
904 	       << "continuous range variables." << std::endl;
905 	  abort_handler(MODEL_ERROR);
906 	}
907       }
908       else if (type == NORMAL_UNCERTAIN) {
909 	if (map2 == "mean")
910 	  active2ACVarMapTargets[curr_index] = Pecos::N_MEAN;
911 	else if (map2 == "std_deviation")
912 	  active2ACVarMapTargets[curr_index] = Pecos::N_STD_DEV;
913 	else if (map2 == "lower_bound")
914 	  active2ACVarMapTargets[curr_index] = Pecos::N_LWR_BND;
915 	else if (map2 == "upper_bound")
916 	  active2ACVarMapTargets[curr_index] = Pecos::N_UPR_BND;
917 	else if (map2 == "location")
918 	  active2ACVarMapTargets[curr_index] = Pecos::N_LOCATION;
919 	else if (map2 == "scale")
920 	  active2ACVarMapTargets[curr_index] = Pecos::N_SCALE;
921 	else {
922 	  Cerr << "\nError: " << map2 << " mapping not supported for "
923 	       << "normal distributions." << std::endl;
924 	  abort_handler(MODEL_ERROR);
925 	}
926       }
927       else if (type == LOGNORMAL_UNCERTAIN) {
928 	// TO DO: 3 parameter lognormal adds a location parameter
929 	if (map2 == "mean")
930 	  active2ACVarMapTargets[curr_index] = Pecos::LN_MEAN;
931 	else if (map2 == "std_deviation")
932 	  active2ACVarMapTargets[curr_index] = Pecos::LN_STD_DEV;
933 	else if (map2 == "lambda")
934 	    active2ACVarMapTargets[curr_index] = Pecos::LN_LAMBDA;
935 	else if (map2 == "zeta")
936 	  active2ACVarMapTargets[curr_index] = Pecos::LN_ZETA;
937 	else if (map2 == "error_factor")
938 	  active2ACVarMapTargets[curr_index] = Pecos::LN_ERR_FACT;
939 	else if (map2 == "lower_bound")
940 	  active2ACVarMapTargets[curr_index] = Pecos::LN_LWR_BND;
941 	else if (map2 == "upper_bound")
942 	  active2ACVarMapTargets[curr_index] = Pecos::LN_UPR_BND;
943 	else {
944 	  Cerr << "\nError: " << map2 << " mapping not supported for "
945 	       << "lognormal distributions." << std::endl;
946 	  abort_handler(MODEL_ERROR);
947 	}
948       }
949       else if (type == UNIFORM_UNCERTAIN) {
950 	if (map2 == "lower_bound")
951 	  active2ACVarMapTargets[curr_index] = Pecos::U_LWR_BND;
952 	else if (map2 == "upper_bound")
953 	  active2ACVarMapTargets[curr_index] = Pecos::U_UPR_BND;
954 	else if (map2 == "location")
955 	  active2ACVarMapTargets[curr_index] = Pecos::U_LOCATION;
956 	else if (map2 == "scale")
957 	  active2ACVarMapTargets[curr_index] = Pecos::U_SCALE;
958 	else {
959 	  Cerr << "\nError: " << map2 << " mapping not supported for "
960 	       << "uniform distributions." << std::endl;
961 	  abort_handler(MODEL_ERROR);
962 	}
963       }
964       else if (type == LOGUNIFORM_UNCERTAIN) {
965 	// TO DO: mean/std deviation or location/scale
966 	if (map2 == "lower_bound")
967 	  active2ACVarMapTargets[curr_index] = Pecos::LU_LWR_BND;
968 	else if (map2 == "upper_bound")
969 	  active2ACVarMapTargets[curr_index] = Pecos::LU_UPR_BND;
970 	else {
971 	  Cerr << "\nError: " << map2 << " mapping not supported for "
972 	       << "loguniform distributions." << std::endl;
973 	  abort_handler(MODEL_ERROR);
974 	}
975       }
976       else if (type == TRIANGULAR_UNCERTAIN) {
977 	if (map2 == "mode")
978 	  active2ACVarMapTargets[curr_index] = Pecos::T_MODE;
979 	else if (map2 == "lower_bound")
980 	  active2ACVarMapTargets[curr_index] = Pecos::T_LWR_BND;
981 	else if (map2 == "upper_bound")
982 	  active2ACVarMapTargets[curr_index] = Pecos::T_UPR_BND;
983 	else if (map2 == "location")
984 	  active2ACVarMapTargets[curr_index] = Pecos::T_LOCATION;
985 	else if (map2 == "scale")
986 	  active2ACVarMapTargets[curr_index] = Pecos::T_SCALE;
987 	else {
988 	  Cerr << "\nError: " << map2 << " mapping not supported for "
989 	       << "triangular distributions." << std::endl;
990 	  abort_handler(MODEL_ERROR);
991 	}
992       }
993       else if (type == EXPONENTIAL_UNCERTAIN) {
994 	// TO DO: mean/std deviation or location/scale
995 	if (map2 == "beta")
996 	  active2ACVarMapTargets[curr_index] = Pecos::E_BETA;
997 	else {
998 	  Cerr << "\nError: " << map2 << " mapping not supported for "
999 	       << "exponential distributions." << std::endl;
1000 	  abort_handler(MODEL_ERROR);
1001 	}
1002       }
1003       else if (type == BETA_UNCERTAIN) {
1004 	// TO DO: mean/std deviation or location/scale
1005 	if (map2 == "alpha")
1006 	  active2ACVarMapTargets[curr_index] = Pecos::BE_ALPHA;
1007 	else if (map2 == "beta")
1008 	  active2ACVarMapTargets[curr_index] = Pecos::BE_BETA;
1009 	else if (map2 == "lower_bound")
1010 	  active2ACVarMapTargets[curr_index] = Pecos::BE_LWR_BND;
1011 	else if (map2 == "upper_bound")
1012 	  active2ACVarMapTargets[curr_index] = Pecos::BE_UPR_BND;
1013 	else {
1014 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1015 	       << "beta distributions." << std::endl;
1016 	  abort_handler(MODEL_ERROR);
1017 	}
1018       }
1019       else if (type == GAMMA_UNCERTAIN) {
1020 	// TO DO: mean/std deviation or location/scale
1021 	if (map2 == "alpha")
1022 	  active2ACVarMapTargets[curr_index] = Pecos::GA_ALPHA;
1023 	else if (map2 == "beta")
1024 	  active2ACVarMapTargets[curr_index] = Pecos::GA_BETA;
1025 	else {
1026 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1027 	       << "gamma distributions." << std::endl;
1028 	  abort_handler(MODEL_ERROR);
1029 	}
1030       }
1031       else if (type == GUMBEL_UNCERTAIN) {
1032 	// TO DO: mean/std deviation or location/scale
1033 	if (map2 == "alpha")
1034 	  active2ACVarMapTargets[curr_index] = Pecos::GU_ALPHA;
1035 	else if (map2 == "beta")
1036 	  active2ACVarMapTargets[curr_index] = Pecos::GU_BETA;
1037 	else {
1038 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1039 	       << "gumbel distributions." << std::endl;
1040 	  abort_handler(MODEL_ERROR);
1041 	}
1042       }
1043       else if (type == FRECHET_UNCERTAIN) {
1044 	// TO DO: mean/std deviation or location/scale
1045 	if (map2 == "alpha")
1046 	  active2ACVarMapTargets[curr_index] = Pecos::F_ALPHA;
1047 	else if (map2 == "beta")
1048 	  active2ACVarMapTargets[curr_index] = Pecos::F_BETA;
1049 	else {
1050 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1051 	       << "frechet distributions." << std::endl;
1052 	  abort_handler(MODEL_ERROR);
1053 	}
1054       }
1055       else if (type == WEIBULL_UNCERTAIN) {
1056 	// TO DO: mean/std deviation or location/scale
1057 	if (map2 == "alpha")
1058 	  active2ACVarMapTargets[curr_index] = Pecos::W_ALPHA;
1059 	else if (map2 == "beta")
1060 	  active2ACVarMapTargets[curr_index] = Pecos::W_BETA;
1061 	else {
1062 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1063 	       << "weibull distributions." << std::endl;
1064 	  abort_handler(MODEL_ERROR);
1065 	}
1066       }
1067       else {
1068 	Cerr << "\nError: " << type << " variable type not supported in "
1069 	     << "secondary real mappings\n       for primary continuous "
1070 	     << "variable targets." << std::endl;
1071 	abort_handler(MODEL_ERROR);
1072       }
1073       active2ADIVarMapTargets[curr_index] = active2ADSVarMapTargets[curr_index]
1074 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1075     }
1076     else if (adi_index1 != _NPOS) {
1077       unsigned short type
1078 	= subModel.all_discrete_int_variable_types()[adi_index1];
1079       if (type == POISSON_UNCERTAIN) {
1080 	if (map2 == "lambda")
1081 	  active2ADIVarMapTargets[curr_index] = Pecos::P_LAMBDA;
1082 	else {
1083 	  Cerr << "\nError: " << map2 << " real mapping not supported for "
1084 	       << "poisson uncertain variables." << std::endl;
1085 	  abort_handler(MODEL_ERROR);
1086 	}
1087       }
1088       if (type == BINOMIAL_UNCERTAIN) {
1089 	if (map2 == "prob_per_trial")
1090 	  active2ADIVarMapTargets[curr_index] = Pecos::BI_P_PER_TRIAL;
1091 	else {
1092 	  Cerr << "\nError: " << map2 << " real mapping not supported for "
1093 	       << "binomial uncertain variables." << std::endl;
1094 	  abort_handler(MODEL_ERROR);
1095 	}
1096       }
1097       if (type == NEGATIVE_BINOMIAL_UNCERTAIN) {
1098 	if (map2 == "prob_per_trial")
1099 	  active2ADIVarMapTargets[curr_index] = Pecos::NBI_P_PER_TRIAL;
1100 	else {
1101 	  Cerr << "\nError: " << map2 << " real mapping not supported for "
1102 	       << "negative binomial uncertain variables." << std::endl;
1103 	  abort_handler(MODEL_ERROR);
1104 	}
1105       }
1106       if (type == GEOMETRIC_UNCERTAIN) {
1107 	if (map2 == "prob_per_trial")
1108 	  active2ADIVarMapTargets[curr_index] = Pecos::GE_P_PER_TRIAL;
1109 	else {
1110 	  Cerr << "\nError: " << map2 << " real mapping not supported for "
1111 	       << "geometric uncertain variables." << std::endl;
1112 	  abort_handler(MODEL_ERROR);
1113 	}
1114       }
1115       else {
1116 	Cerr << "\nError: " << type << " variable type not supported in "
1117 	     << "secondary real mappings\n       for primary discrete integer "
1118 	     << "variable targets." << std::endl;
1119 	abort_handler(MODEL_ERROR);
1120       }
1121       active2ACVarMapTargets[curr_index] = active2ADSVarMapTargets[curr_index]
1122 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1123     }
1124     else if (ads_index1 != _NPOS) {
1125       unsigned short type
1126 	= subModel.all_discrete_string_variable_types()[ads_index1];
1127       Cerr << "\nError: " << type << " variable type not supported in "
1128 	   << "secondary real mappings\n       for primary discrete string "
1129 	   << "variable targets." << std::endl;
1130       abort_handler(MODEL_ERROR);
1131       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1132 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1133     }
1134     else if (adr_index1 != _NPOS) {
1135       unsigned short type
1136 	= subModel.all_discrete_real_variable_types()[adr_index1];
1137       Cerr << "\nError: " << type << " variable type not supported in "
1138 	   << "secondary real mappings\n       for primary discrete real "
1139 	   << "variable targets." << std::endl;
1140       abort_handler(MODEL_ERROR);
1141       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1142 	= active2ADSVarMapTargets[curr_index] = Pecos::NO_TARGET;
1143     }
1144   }
1145 }
1146 
1147 
1148 void NestedModel::
resolve_integer_variable_mapping(const String & map1,const String & map2,size_t curr_index,short & inactive_sm_view)1149 resolve_integer_variable_mapping(const String& map1, const String& map2,
1150 				 size_t curr_index, short& inactive_sm_view)
1151 {
1152   size_t ac_index1, adi_index1, ads_index1, adr_index1;
1153   resolve_map1(map1, ac_index1, adi_index1, ads_index1, adr_index1,
1154 	       curr_index, inactive_sm_view);
1155 
1156   if (!active2ACVarMapTargets.empty()) { // indicates num_var_map_2
1157     if (map2.empty())
1158       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1159 	= active2ADSVarMapTargets[curr_index]
1160 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1161     else if (ac_index1 != _NPOS) {
1162       unsigned short type = subModel.all_continuous_variable_types()[ac_index1];
1163       Cerr << "\nError: " << type << " variable type not supported in "
1164 	   << "secondary integer mappings\n       for primary continuous "
1165 	   << "variable targets." << std::endl;
1166       abort_handler(MODEL_ERROR);
1167       active2ADIVarMapTargets[curr_index] = active2ADSVarMapTargets[curr_index]
1168 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1169     }
1170     else if (adi_index1 != _NPOS) {
1171       unsigned short type
1172 	= subModel.all_discrete_int_variable_types()[adi_index1];
1173       if (type == DISCRETE_DESIGN_RANGE || type == DISCRETE_STATE_RANGE) {
1174 	if (map2 == "lower_bound")
1175 	  active2ADIVarMapTargets[curr_index] = Pecos::DR_LWR_BND;
1176 	else if (map2 == "upper_bound")
1177 	  active2ADIVarMapTargets[curr_index] = Pecos::DR_UPR_BND;
1178 	else {
1179 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1180 	       << "discrete range variables." << std::endl;
1181 	  abort_handler(MODEL_ERROR);
1182 	}
1183       }
1184       else if (type == BINOMIAL_UNCERTAIN) {
1185 	if (map2 == "num_trials")
1186 	  active2ADIVarMapTargets[curr_index] = Pecos::BI_TRIALS;
1187 	else {
1188 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1189 	       << "binomial uncertain variables." << std::endl;
1190 	  abort_handler(MODEL_ERROR);
1191 	}
1192       }
1193       else if (type == NEGATIVE_BINOMIAL_UNCERTAIN) {
1194 	if (map2 == "num_trials")
1195 	  active2ADIVarMapTargets[curr_index] = Pecos::NBI_TRIALS;
1196 	else {
1197 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1198 	       << "negative binomial uncertain variables." << std::endl;
1199 	  abort_handler(MODEL_ERROR);
1200 	}
1201       }
1202       else if (type == HYPERGEOMETRIC_UNCERTAIN) {
1203 	if (map2 == "total_population")
1204 	  active2ADIVarMapTargets[curr_index] = Pecos::HGE_TOT_POP;
1205 	else if (map2 == "selected_population")
1206 	  active2ADIVarMapTargets[curr_index] = Pecos::HGE_SEL_POP;
1207 	else if (map2 == "num_drawn")
1208 	  active2ADIVarMapTargets[curr_index] = Pecos::HGE_DRAWN;
1209 	else {
1210 	  Cerr << "\nError: " << map2 << " mapping not supported for "
1211 	       << "hypergeometric uncertain variables." << std::endl;
1212 	  abort_handler(MODEL_ERROR);
1213 	}
1214       }
1215       else {
1216 	Cerr << "\nError: " << type << " variable type not supported in "
1217 	     << "secondary integer mappings\n       for primary discrete "
1218 	     << "integer variable targets." << std::endl;
1219 	abort_handler(MODEL_ERROR);
1220       }
1221       active2ACVarMapTargets[curr_index] = active2ADSVarMapTargets[curr_index]
1222 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1223     }
1224     else if (ads_index1 != _NPOS) {
1225       unsigned short type
1226 	= subModel.all_discrete_string_variable_types()[ads_index1];
1227       Cerr << "\nError: " << type << " variable type not supported in "
1228 	   << "secondary integer mappings\n       for primary discrete string "
1229 	   << "variable targets." << std::endl;
1230       abort_handler(MODEL_ERROR);
1231       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1232 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1233     }
1234     else if (adr_index1 != _NPOS) {
1235       unsigned short type
1236 	= subModel.all_discrete_real_variable_types()[adr_index1];
1237       Cerr << "\nError: " << type << " variable type not supported in "
1238 	   << "secondary integer mappings\n       for primary discrete real "
1239 	   << "variable targets." << std::endl;
1240       abort_handler(MODEL_ERROR);
1241       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1242 	= active2ADSVarMapTargets[curr_index] = Pecos::NO_TARGET;
1243     }
1244   }
1245 }
1246 
1247 
1248 void NestedModel::
resolve_string_variable_mapping(const String & map1,const String & map2,size_t curr_index,short & inactive_sm_view)1249 resolve_string_variable_mapping(const String& map1, const String& map2,
1250 				size_t curr_index, short& inactive_sm_view)
1251 {
1252   size_t ac_index1, adi_index1, ads_index1, adr_index1;
1253   resolve_map1(map1, ac_index1, adi_index1, ads_index1, adr_index1,
1254 	       curr_index, inactive_sm_view);
1255 
1256   if (!active2ACVarMapTargets.empty()) { // indicates num_var_map_2
1257     if (map2.empty())
1258       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1259 	= active2ADSVarMapTargets[curr_index]
1260 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1261     else if (ac_index1 != _NPOS) {
1262       unsigned short type = subModel.all_continuous_variable_types()[ac_index1];
1263       Cerr << "\nError: " << type << " variable type not supported in "
1264 	   << "secondary string mappings\n       for primary continuous "
1265 	   << "variable targets." << std::endl;
1266       abort_handler(MODEL_ERROR);
1267       active2ADIVarMapTargets[curr_index] = active2ADSVarMapTargets[curr_index]
1268 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1269     }
1270     else if (adi_index1 != _NPOS) {
1271       unsigned short type
1272 	= subModel.all_discrete_int_variable_types()[adi_index1];
1273       Cerr << "\nError: " << type << " variable type not supported in "
1274 	   << "secondary string mappings\n       for primary discrete integer "
1275 	   << "variable targets." << std::endl;
1276       abort_handler(MODEL_ERROR);
1277       active2ACVarMapTargets[curr_index] = active2ADSVarMapTargets[curr_index]
1278 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1279     }
1280     else if (ads_index1 != _NPOS) {
1281       unsigned short type
1282 	= subModel.all_discrete_string_variable_types()[ads_index1];
1283       Cerr << "\nError: " << type << " variable type not supported in "
1284 	   << "secondary string mappings\n       for primary discrete string "
1285 	   << "variable targets." << std::endl;
1286       abort_handler(MODEL_ERROR);
1287       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1288 	= active2ADRVarMapTargets[curr_index] = Pecos::NO_TARGET;
1289     }
1290     else if (adr_index1 != _NPOS) {
1291       unsigned short type
1292 	= subModel.all_discrete_real_variable_types()[adr_index1];
1293       Cerr << "\nError: " << type << " variable type not supported in "
1294 	   << "secondary string mappings\n       for primary discrete real "
1295 	   << "variable targets." << std::endl;
1296       abort_handler(MODEL_ERROR);
1297       active2ACVarMapTargets[curr_index] = active2ADIVarMapTargets[curr_index]
1298 	= active2ADSVarMapTargets[curr_index] = Pecos::NO_TARGET;
1299     }
1300   }
1301 }
1302 
1303 
1304 /** Update subModel's inactive variables with active variables from
1305     currentVariables, compute the optional interface and sub-iterator
1306     responses, and map these to the total model response. */
derived_evaluate(const ActiveSet & set)1307 void NestedModel::derived_evaluate(const ActiveSet& set)
1308 {
1309   ++nestedModelEvalCntr;
1310 
1311   // Set currentResponse asv and extract opt_interface_set/sub_iterator_set
1312   currentResponse.active_set(set); currentResponse.reset();
1313   ActiveSet opt_interface_set, sub_iterator_set;
1314   bool      opt_interface_map, sub_iterator_map;
1315   set_mapping(set, opt_interface_set, opt_interface_map,
1316 	           sub_iterator_set,  sub_iterator_map);
1317 
1318   // Perform optionalInterface map (opt_interface_set is updated within
1319   // optInterfaceResponse by map):
1320   if (opt_interface_map) {
1321     Cout << "\n----------------------------------------------------------------"
1322 	 << "--\nNestedModel Evaluation " << std::setw(4) << nestedModelEvalCntr
1323 	 << ": performing optional interface mapping\n-------------------------"
1324 	 << "-----------------------------------------\n";
1325     component_parallel_mode(INTERFACE_MODE);
1326     if (hierarchicalTagging) {
1327       String eval_tag = evalTagPrefix + '.' +
1328 	std::to_string(nestedModelEvalCntr);
1329       // don't apply a redundant interface eval id
1330       bool append_iface_tag = false;
1331       optionalInterface.eval_tag_prefix(eval_tag, append_iface_tag);
1332     }
1333 
1334     ParConfigLIter pc_iter = parallelLib.parallel_configuration_iterator();
1335     parallelLib.parallel_configuration_iterator(modelPCIter);
1336     if(interfEvaluationsDBState == EvaluationsDBState::UNINITIALIZED)
1337       interfEvaluationsDBState = evaluationsDB.interface_allocate(modelId,
1338           interface_id(), "simulation", currentVariables, optInterfaceResponse,
1339           default_interface_active_set(), optionalInterface.analysis_components());
1340 
1341     optionalInterface.map(currentVariables, opt_interface_set,
1342 			  optInterfaceResponse);
1343     if(interfEvaluationsDBState == EvaluationsDBState::ACTIVE) {
1344       evaluationsDB.store_interface_variables(modelId, interface_id(),
1345           optionalInterface.evaluation_id(), opt_interface_set, currentVariables);
1346       evaluationsDB.store_interface_response(modelId, interface_id(),
1347           optionalInterface.evaluation_id(), optInterfaceResponse);
1348     }
1349     parallelLib.parallel_configuration_iterator(pc_iter); // restore
1350 
1351     // map optInterface results into their contribution to currentResponse
1352     interface_response_overlay(optInterfaceResponse, currentResponse);
1353   }
1354 
1355   if (sub_iterator_map) {
1356     //++subIteratorJobCntr; // does not encompass blocking evals
1357 
1358     // need comm set up and master break off
1359     // (see IteratorScheduler::run_iterator())
1360     Cout << "\n-------------------------------------------------\nNestedModel "
1361 	 << "Evaluation " << std::setw(4) << nestedModelEvalCntr << ": running "
1362 	 << "sub_iterator\n-------------------------------------------------\n";
1363     component_parallel_mode(SUB_MODEL_MODE);
1364     update_sub_model(currentVariables, userDefinedConstraints);
1365     subIterator.response_results_active_set(sub_iterator_set);
1366     if (hierarchicalTagging) {
1367       String eval_tag = evalTagPrefix + '.' +
1368 	std::to_string(nestedModelEvalCntr);
1369       subIterator.eval_tag_prefix(eval_tag);
1370     }
1371 
1372     ParLevLIter pl_iter
1373       = modelPCIter->mi_parallel_level_iterator(subIteratorSched.miPLIndex);
1374     if (subIteratorSched.messagePass) {
1375       // For derived_evaluate(), subIterator scheduling would not
1376       // normally be expected, but singleton jobs could use this fn assuming
1377       // no dedicated master overload (enforced in Model::evaluate()).
1378       // Given this protection, don't schedule the job -- execute it locally.
1379       if (subIteratorSched.iteratorScheduling == PEER_SCHEDULING &&
1380 	  subIteratorSched.peerAssignJobs) {
1381 	// match 2 bcasts in IteratorScheduler::peer_static_schedule_iterators()
1382 	// needed by procs in NestedModel::serve_run()
1383 	int num_jobs = 1;
1384 	parallelLib.bcast_hs(num_jobs, *pl_iter); // over pl.hubServerIntraComm
1385 	if (subIteratorSched.iteratorCommSize > 1)
1386 	  parallelLib.bcast(num_jobs, *pl_iter);  // over pl.serverIntraComm
1387       }
1388       // run_iterator() is used since we stop subModel servers for consistency
1389       // with fall through behavior of schedule_iterators()
1390       subIteratorSched.run_iterator(subIterator, pl_iter);
1391       if (subIteratorSched.iteratorScheduling == MASTER_SCHEDULING)
1392 	subIteratorSched.stop_iterator_servers();
1393 
1394       /* This approach has 2 issues: (1) a single-processor subIterator job is
1395 	 always assigned by master to server 1 (ded master overload bypassed),
1396 	 (2) peer static init/update bookkeeping is redundant of above/below.
1397       subIteratorSched.numIteratorJobs = 1;
1398       // can use shallow copy for queue of 1 job (avoids need to copy updated
1399       // entry in subIteratorPRPQueue back to subIterator.response_results())
1400       ParamResponsePair current_pair(currentVariables, subIterator.method_id(),
1401 				     subIterator.response_results(), 1, false);
1402       subIteratorPRPQueue.insert(current_pair);
1403       subIteratorSched.schedule_iterators(*this, subIterator);
1404       */
1405     }
1406     else // run_iterator() is not used since we don't stop subModel servers
1407          // until change in component_parallel_mode
1408       subIterator.run(pl_iter);
1409 
1410     const Response& sub_iter_resp = subIterator.response_results();
1411     Cout << "\nActive response data from sub_iterator:\n"<< sub_iter_resp<<'\n';
1412     // map subIterator results into their contribution to currentResponse
1413     iterator_response_overlay(sub_iter_resp, currentResponse);
1414   }
1415 
1416   Cout << "\n---------------------------\nNestedModel Evaluation "
1417        << std::setw(4) << nestedModelEvalCntr << " results:"
1418        << "\n---------------------------\n";
1419   // if secondary variable mappings (distribution parameter insertions), the
1420   // Nested parameters are not directly reflected in the subIterator output:
1421   if (outputLevel >= VERBOSE_OUTPUT && !active2ACVarMapTargets.empty())
1422     Cout << "Nested parameters:\n" << currentVariables;
1423   Cout << "\nActive response data from nested mapping:\n" << currentResponse
1424        << '\n';
1425 }
1426 
1427 
1428 /** Asynchronous execution of subIterator on subModel and, optionally,
1429     optionalInterface. */
derived_evaluate_nowait(const ActiveSet & set)1430 void NestedModel::derived_evaluate_nowait(const ActiveSet& set)
1431 {
1432   ++nestedModelEvalCntr;
1433 
1434   // Set currentResponse asv and extract opt_interface_set/sub_iterator_set
1435   currentResponse.active_set(set);
1436   ActiveSet opt_interface_set, sub_iterator_set;
1437   bool      opt_interface_map, sub_iterator_map;
1438   set_mapping(set, opt_interface_set, opt_interface_map,
1439 	           sub_iterator_set,  sub_iterator_map);
1440 
1441   // Perform optionalInterface map (opt_interface_set is updated within
1442   // optInterfaceResponse by map):
1443   if (opt_interface_map) {
1444     Cout << "\n----------------------------------------------------------------"
1445 	 << "--\nNestedModel Evaluation " << std::setw(4)
1446 	 << nestedModelEvalCntr << ": queueing optional interface mapping\n"
1447 	 << "------------------------------------------------------------------"
1448 	 << '\n';
1449     // don't need to set component parallel mode since this only queues the job
1450     if(interfEvaluationsDBState == EvaluationsDBState::UNINITIALIZED)
1451       interfEvaluationsDBState = evaluationsDB.interface_allocate(modelId, interface_id(),
1452           "simulation", currentVariables, optInterfaceResponse, default_interface_active_set(),
1453           optionalInterface.analysis_components());
1454     optionalInterface.map(currentVariables, opt_interface_set,
1455 			  optInterfaceResponse, true);
1456     if(interfEvaluationsDBState == EvaluationsDBState::ACTIVE)
1457       evaluationsDB.store_interface_variables(modelId, interface_id(),
1458           optionalInterface.evaluation_id(), opt_interface_set, currentVariables);
1459     optInterfaceIdMap[optionalInterface.evaluation_id()] = nestedModelEvalCntr;
1460   }
1461 
1462   if (sub_iterator_map) {
1463     ++subIteratorJobCntr;
1464 
1465     // need comm set up and master break off
1466     // (see IteratorScheduler::run_iterator())
1467     Cout << "\n-------------------------------------------------\n"
1468 	 << "NestedModel Evaluation " << std::setw(4) << nestedModelEvalCntr
1469 	 << ": queueing sub_iterator"
1470 	 << "\n-------------------------------------------------\n";
1471 
1472     // load up queue of iterator jobs to be scheduled in derived synchronize:
1473     // > use the subIterator's method id as the PRP interface id
1474     // > subIterator's execNum could potentially be used as execution counter
1475     //   for id mapping but it isn't defined/incremented until run time (see
1476     //   Iterator::run())
1477     // > simplest approach of tagging with nestedModelEvalCntr is sufficient,
1478     //   since we do not need to map from a set of eval ids returned from
1479     //   lower level bookkeeping
1480     subIterator.response_results_active_set(sub_iterator_set);
1481     ParamResponsePair current_pair(currentVariables, subIterator.method_id(),
1482 				   subIterator.response_results(),
1483 				   nestedModelEvalCntr);
1484     subIteratorPRPQueue.insert(current_pair);
1485 
1486     // update bookkeeping for job_index mappings in IteratorScheduler callbacks
1487     subIteratorIdMap[subIteratorJobCntr] = nestedModelEvalCntr;
1488   }
1489 }
1490 
1491 
1492 /** Recovery of asynchronous subIterator executions and, optionally,
1493     asynchronous optionalInterface mappings. */
derived_synchronize()1494 const IntResponseMap& NestedModel::derived_synchronize()
1495 {
1496   nestedResponseMap.clear();
1497 
1498   // TO DO: optInt/subIter scheduling is currently sequential, but could be
1499   // overlapped as in HierarchSurrModel, given IteratorScheduler nowait support
1500 
1501   IntIntMIter id_it; IntRespMCIter r_cit;
1502   if (!optInterfacePointer.empty()) {
1503     component_parallel_mode(INTERFACE_MODE);
1504 
1505     ParConfigLIter pc_iter = parallelLib.parallel_configuration_iterator();
1506     parallelLib.parallel_configuration_iterator(modelPCIter);
1507     const IntResponseMap& opt_int_resp_map = optionalInterface.synchronize();
1508     parallelLib.parallel_configuration_iterator(pc_iter); // restore
1509 
1510     // overlay response sets
1511     r_cit = opt_int_resp_map.begin();
1512     while (r_cit != opt_int_resp_map.end()) {
1513       int oi_eval_id = r_cit->first;
1514       id_it = optInterfaceIdMap.find(oi_eval_id);
1515       if (id_it != optInterfaceIdMap.end()) {
1516 	interface_response_overlay(r_cit->second,
1517 				   nested_response(id_it->second));
1518 	optInterfaceIdMap.erase(id_it);
1519 	++r_cit;
1520       }
1521       else { // see also Model::rekey_synch()
1522 	++r_cit; // prior to invalidation from erase within cache_unmatched
1523 	optionalInterface.cache_unmatched_response(oi_eval_id);
1524       }
1525     }
1526   }
1527 
1528   if (!subIteratorPRPQueue.empty()) {
1529     // schedule subIteratorPRPQueue jobs
1530     component_parallel_mode(SUB_MODEL_MODE);
1531     subIteratorSched.numIteratorJobs = subIteratorPRPQueue.size();
1532     subIteratorSched.schedule_iterators(*this, subIterator);
1533     // overlay response sets (no rekey or cache necessary)
1534     for (PRPQueueIter q_it=subIteratorPRPQueue.begin();
1535 	 q_it!=subIteratorPRPQueue.end(); ++q_it)
1536       iterator_response_overlay(q_it->response(),
1537 				nested_response(q_it->eval_id()));
1538     // clear sub-iterator jobs
1539     subIteratorPRPQueue.clear();
1540     // Reset bookkeeping used in IteratorScheduler callbacks (e.g.,
1541     // {pack,unpack}_* in NestedModel.hpp); sub-iterator job counter
1542     // mirrors the passed job_index and maps to nestedModelEvalCntr
1543     // for subIteratorPRPQueue lookups.
1544     subIteratorIdMap.clear(); subIteratorJobCntr = 0;
1545   }
1546 
1547   //nestedVarsMap.clear();
1548   for (r_cit=nestedResponseMap.begin(); r_cit!=nestedResponseMap.end(); ++r_cit)
1549     Cout << "\n---------------------------\nNestedModel Evaluation "
1550 	 << std::setw(4) << r_cit->first << " total response:"
1551 	 << "\n---------------------------\n\nActive response data "
1552 	 << "from nested mapping:\n" << r_cit->second << '\n';
1553   return nestedResponseMap;
1554 }
1555 
1556 
1557 /* Asynchronous response computations are not currently supported by
1558    NestedModels.  Return a dummy to satisfy the compiler.
1559 const IntResponseMap& NestedModel::derived_synchronize_nowait()
1560 {
1561   // TO DO: will require nowait support in IteratorScheduler
1562 
1563   //nestedVarsMap.erase(eval_id);
1564   return nestedResponseMap;
1565 }
1566 */
1567 
1568 
1569 void NestedModel::
set_mapping(const ActiveSet & mapped_set,ActiveSet & opt_interface_set,bool & opt_interface_map,ActiveSet & sub_iterator_set,bool & sub_iterator_map)1570 set_mapping(const ActiveSet& mapped_set, ActiveSet& opt_interface_set,
1571 	    bool& opt_interface_map,     ActiveSet& sub_iterator_set,
1572 	    bool& sub_iterator_map)
1573 {
1574   size_t i, j, num_opt_interf_con = numOptInterfIneqCon + numOptInterfEqCon,
1575     num_mapped_primary = std::max(numOptInterfPrimary, subIterMappedPri);
1576 
1577   const ShortArray& mapped_asv = mapped_set.request_vector();
1578   const SizetArray& mapped_dvv = mapped_set.derivative_vector();
1579 
1580   if (mapped_asv.size() != num_mapped_primary + num_opt_interf_con +
1581                            subIterMappedSec) {
1582     Cerr << "\nError: mismatch is ASV lengths in NestedModel::set_mapping()."
1583 	 << "\n       expected " << mapped_asv.size() << " total, received "
1584 	 << num_mapped_primary << " primary plus " << num_opt_interf_con +
1585             subIterMappedSec << " secondary." << std::endl;
1586     abort_handler(MODEL_ERROR);
1587   }
1588 
1589   // sub_iterator_asv:
1590 
1591   sub_iterator_map = false;
1592   ShortArray sub_iterator_asv(numSubIterFns, 0);
1593   // augment sub_iterator_asv based on mapped primary asv and primaryRespCoeffs
1594   for (i=0; i<subIterMappedPri; ++i) {
1595     short mapped_asv_val = mapped_asv[i];
1596     if (mapped_asv_val) {
1597       if (identityRespMap)
1598 	sub_iterator_asv[i] |= mapped_asv_val;
1599       else
1600 	for (j=0; j<numSubIterFns; ++j)
1601 	  if (std::fabs(primaryRespCoeffs(i,j)) > DBL_MIN)
1602 	    sub_iterator_asv[j] |= mapped_asv_val;
1603       sub_iterator_map = true;
1604     }
1605   }
1606   // augment sub_iterator_asv based on mapped constr asv and secondaryRespCoeffs
1607   for (i=0; i<subIterMappedSec; ++i) {
1608     short mapped_asv_val = mapped_asv[i+num_mapped_primary+num_opt_interf_con];
1609     if (mapped_asv_val) {
1610       if (identityRespMap)
1611 	sub_iterator_asv[subIterMappedPri + i] |= mapped_asv_val;
1612       else
1613 	for (j=0; j<numSubIterFns; ++j)
1614 	  if (std::fabs(secondaryRespCoeffs(i,j)) > DBL_MIN)
1615 	    sub_iterator_asv[j] |= mapped_asv_val;
1616       sub_iterator_map = true;
1617     }
1618   }
1619   if (sub_iterator_map) {
1620     sub_iterator_set.request_vector(sub_iterator_asv);
1621 
1622     // sub_iterator_dvv:
1623 
1624     SizetMultiArrayConstView cv_ids
1625       = currentVariables.continuous_variable_ids();
1626     SizetMultiArrayConstView sm_acv_ids
1627       = subModel.all_continuous_variable_ids();
1628 
1629     /* Old Approach: subIterator must decipher/replace top-level DVV as reqd.
1630     // Note: the ordering of top-level active variables may differ from the
1631     // ordering of the inactive sub-model variables.  However, the subModel's
1632     // inactive_continuous_variable_ids() omit inserted variables, which would
1633     // result in erroneous subIterator final response array sizing.  Therefore,
1634     // use the top-level active variable ids and rely on the subIterator to
1635     // perform the mappings to the subModel variables.
1636     sub_iterator_set.derivative_vector(cv_ids);
1637     */
1638 
1639     size_t num_mapped_dvv = mapped_dvv.size();
1640     SizetArray sub_iterator_dvv;
1641     bool var_map_2 = !active2ACVarMapTargets.empty();
1642     if (!var_map_2)
1643       sub_iterator_dvv.resize(num_mapped_dvv);
1644     for (i=0; i<num_mapped_dvv; ++i) {
1645       size_t cv_index = find_index(cv_ids, mapped_dvv[i]);
1646       if (cv_index == _NPOS) {
1647 	Cerr << "\nError: NestedModel DVV component not contained within "
1648 	     << "active continuous variable ids." << std::endl;
1649 	abort_handler(MODEL_ERROR);
1650       }
1651       size_t sm_acv_index = active1ACVarMapIndices[cv_index],
1652 	     sm_acv_id    = sm_acv_ids[sm_acv_index];
1653       if (var_map_2) { // enforce uniqueness in insertion targets
1654 	if (!contains(sub_iterator_dvv, sm_acv_id))
1655 	  sub_iterator_dvv.push_back(sm_acv_id);
1656       }
1657       else
1658 	sub_iterator_dvv[i] = sm_acv_id;
1659     }
1660     sub_iterator_set.derivative_vector(sub_iterator_dvv);
1661     //Cout << "\nmapped_dvv:\n" << mapped_dvv << "\nsub_iterator_dvv:\n"
1662     //     << sub_iterator_dvv << '\n';
1663   }
1664 
1665   // opt_interface_asv:
1666 
1667   // num_mapped_primary >= numOptInterfPrimary with 1-to-1 correspondence
1668   // up to the cut off
1669   opt_interface_map = false;
1670   size_t num_opt_interf_fns = numOptInterfPrimary+num_opt_interf_con;
1671   ShortArray opt_interface_asv(num_opt_interf_fns);
1672   for (i=0; i<numOptInterfPrimary; ++i)
1673     opt_interface_asv[i] = mapped_asv[i];
1674   // num_opt_interf_con has 1-to-1 correspondence with different offsets
1675   for (i=0; i<num_opt_interf_con; ++i)
1676     opt_interface_asv[i+numOptInterfPrimary] = mapped_asv[i+num_mapped_primary];
1677   // Special case of forcing an optional interface execution that lacks fns:
1678   // this allows usage of the optional interface to generate data used only
1679   // by the sub-iterator.  Put another way, the optional interface is active
1680   // unless functions are present and all functions are inactive.
1681   if (!optInterfacePointer.empty() && num_opt_interf_fns == 0)
1682     opt_interface_map = sub_iterator_map;
1683   else // normal case of mapping optional interface fns
1684     for (i=0; i<num_opt_interf_fns; ++i)
1685       if (opt_interface_asv[i])
1686 	{ opt_interface_map = true; break; }
1687   if (opt_interface_map) {
1688     opt_interface_set.request_vector(opt_interface_asv);
1689 
1690     // opt_interface_dvv:
1691 
1692     opt_interface_set.derivative_vector(mapped_dvv);
1693   }
1694 }
1695 
1696 
1697 void NestedModel::
interface_response_overlay(const Response & opt_interface_response,Response & mapped_response)1698 interface_response_overlay(const Response& opt_interface_response,
1699 			   Response& mapped_response)
1700 {
1701   // mapped data initialization
1702   const ShortArray& mapped_asv = mapped_response.active_set_request_vector();
1703   const SizetArray& mapped_dvv = mapped_response.active_set_derivative_vector();
1704   size_t i, m_index, oi_index, num_mapped_fns = mapped_asv.size();
1705   bool deriv_flag = false;
1706   for (i=0; i<num_mapped_fns; ++i)
1707     { if (mapped_asv[i] & 6) deriv_flag = true; break; }
1708   // Sanity checks: the optional interface response must have the same DVV
1709   if ( deriv_flag &&
1710        opt_interface_response.active_set_derivative_vector() != mapped_dvv ) {
1711     Cerr << "\nError: derivative variables vector mismatch in NestedModel::"
1712          << "interface_response_overlay()." << std::endl;
1713     abort_handler(MODEL_ERROR);
1714   }
1715   check_response_map(mapped_asv);
1716 
1717   // build mapped response data:
1718 
1719   // {f}:
1720   for (i=0; i<numOptInterfPrimary; ++i) {
1721     if (mapped_asv[i] & 1)
1722       mapped_response.function_value(
1723 	opt_interface_response.function_value(i), i);
1724     if (mapped_asv[i] & 2)
1725       mapped_response.function_gradient(
1726 	opt_interface_response.function_gradient_view(i), i);
1727     if (mapped_asv[i] & 4)
1728       mapped_response.function_hessian(
1729 	opt_interface_response.function_hessian(i), i);
1730   }
1731 
1732   // {g}:
1733   size_t num_opt_interf_con = numOptInterfIneqCon + numOptInterfEqCon,
1734     num_mapped_1
1735       = std::max(numOptInterfPrimary, subIterMappedPri);
1736   for (i=0; i<num_opt_interf_con; ++i) {
1737     oi_index = numOptInterfPrimary + i;
1738     m_index  = num_mapped_1 + i; // {g_l} <= {g} <= {g_u}
1739     if (i>=numOptInterfIneqCon)  //          {g} == {g_t}
1740       m_index += numSubIterMappedIneqCon;
1741     if (mapped_asv[m_index] & 1) // mapped_vals
1742       mapped_response.function_value(
1743 	opt_interface_response.function_value(oi_index), m_index);
1744     if (mapped_asv[m_index] & 2) // mapped_grads
1745       mapped_response.function_gradient(
1746 	opt_interface_response.function_gradient_view(oi_index), m_index);
1747     if (mapped_asv[m_index] & 4) // mapped_hessians
1748       mapped_response.function_hessian(
1749 	opt_interface_response.function_hessian(oi_index), m_index);
1750   }
1751 }
1752 
1753 
1754 void NestedModel::
iterator_response_overlay(const Response & sub_iterator_response,Response & mapped_response)1755 iterator_response_overlay(const Response& sub_iterator_response,
1756 			  Response& mapped_response)
1757 {
1758   // mapped data initialization
1759   const ShortArray& mapped_asv = mapped_response.active_set_request_vector();
1760   const SizetArray& mapped_dvv = mapped_response.active_set_derivative_vector();
1761   size_t i, j, k, l, m_index, num_mapped_fns = mapped_asv.size(),
1762     num_mapped_deriv_vars = mapped_dvv.size();
1763   bool deriv_flag = false;
1764   for (i=0; i<num_mapped_fns; ++i)
1765     if (mapped_asv[i] & 6) { deriv_flag = true; break; }
1766   // Sanity checks: the derivatives in the sub-iterator response must be with
1767   // respect to the same variables; but since the numbering may be different
1768   // following insertion/augmentation, only the DVV length is verified.
1769   if ( deriv_flag && num_mapped_deriv_vars !=
1770        sub_iterator_response.active_set_derivative_vector().size() ) {
1771     Cerr << "\nError: derivative variables vector mismatch in NestedModel::"
1772          << "iterator_response_overlay()." << std::endl;
1773     abort_handler(MODEL_ERROR);
1774   }
1775   check_response_map(mapped_asv);
1776 
1777   // build mapped response data:
1778 
1779   RealVector mapped_vals = mapped_response.function_values_view();
1780   const RealVector& sub_iterator_vals = sub_iterator_response.function_values();
1781   const RealMatrix& sub_iterator_grads
1782     = sub_iterator_response.function_gradients();
1783   const RealSymMatrixArray& sub_iterator_hessians
1784     = sub_iterator_response.function_hessians();
1785   Real coeff;
1786 
1787   // [W]{S}:
1788   for (i=0; i<subIterMappedPri; ++i) {
1789     if (mapped_asv[i] & 1) { // mapped_vals
1790       if (identityRespMap)
1791 	mapped_vals[i] = sub_iterator_vals[i];
1792       else {
1793 	Real& inner_prod = mapped_vals[i];
1794 	for (j=0; j<numSubIterFns; ++j) {
1795 	  coeff = primaryRespCoeffs(i,j);
1796 	  if (coeff != 0.) // avoid propagation of nan/inf for no mapping
1797 	    inner_prod += coeff * sub_iterator_vals[j]; // [W]{S}
1798 	}
1799       }
1800     }
1801     if (mapped_asv[i] & 2) { // mapped_grads
1802       RealVector mapped_grad = mapped_response.function_gradient_view(i);
1803       for (j=0; j<num_mapped_deriv_vars; ++j) {
1804 	if (identityRespMap)
1805 	  mapped_grad[j] = sub_iterator_grads(j,i);
1806 	else {
1807 	  Real& inner_prod = mapped_grad[j]; // [W]{S}
1808 	  for (k=0; k<numSubIterFns; ++k) {
1809 	    coeff = primaryRespCoeffs(i,k);
1810 	    if (coeff != 0.) // avoid propagation of nan/inf for no mapping
1811 	      inner_prod += coeff * sub_iterator_grads(j,k);
1812 	  }
1813 	}
1814       }
1815     }
1816     if (mapped_asv[i] & 4) { // mapped_hessians
1817       RealSymMatrix mapped_hess = mapped_response.function_hessian_view(i);
1818       for (j=0; j<num_mapped_deriv_vars; ++j) {
1819 	for (k=0; k<=j; ++k) {
1820 	  if (identityRespMap)
1821 	    mapped_hess(j,k) = sub_iterator_hessians[i](j,k);
1822 	  else {
1823 	    Real& inner_prod = mapped_hess(j,k); // [W]{S}
1824 	    for (l=0; l<numSubIterFns; ++l) {
1825 	      coeff = primaryRespCoeffs(i,l);
1826 	      if (coeff != 0.) // avoid propagation of nan/inf for no mapping
1827 		inner_prod += coeff * sub_iterator_hessians[l](j,k);
1828 	    }
1829 	  }
1830 	}
1831       }
1832     }
1833   }
1834 
1835   // [A]{S}:
1836   size_t num_mapped_1 = std::max(numOptInterfPrimary, subIterMappedPri);
1837   for (i=0; i<subIterMappedSec; ++i) {
1838     m_index = num_mapped_1 + numOptInterfIneqCon + i;// {a_l} <= [A]{S} <= {a_u}
1839     if (i>=numSubIterMappedIneqCon)
1840       m_index += numOptInterfEqCon;                           // [A]{S} == {a_t}
1841     // BMA: Using m_index is only safe for the identity map case
1842     // because optional interfaces are disallowed and all the opt
1843     // interf sizes will be 0
1844     if (mapped_asv[m_index] & 1) { // mapped_vals
1845       if (identityRespMap)
1846 	mapped_vals[m_index] = sub_iterator_vals[m_index];
1847       else {
1848 	Real& inner_prod = mapped_vals[m_index]; inner_prod = 0.;
1849 	for (j=0; j<numSubIterFns; ++j) {
1850 	  coeff = secondaryRespCoeffs(i,j);
1851 	  if (coeff != 0.) // avoid propagation of nan/inf for no mapping
1852 	    inner_prod += coeff * sub_iterator_vals[j];
1853 	}
1854       }
1855     }
1856     if (mapped_asv[m_index] & 2) { // mapped_grads
1857       RealVector mapped_grad = mapped_response.function_gradient_view(m_index);
1858       for (j=0; j<num_mapped_deriv_vars; ++j) {
1859 	if (identityRespMap)
1860 	  mapped_grad[m_index] = sub_iterator_grads(j,m_index);
1861 	else {
1862 	  Real& inner_prod = mapped_grad[j]; inner_prod = 0.;
1863 	  for (k=0; k<numSubIterFns; ++k) {
1864 	    coeff = secondaryRespCoeffs(i,k);
1865 	    if (coeff != 0.) // avoid propagation of nan/inf for no mapping
1866 	      inner_prod += coeff * sub_iterator_grads(j,k);
1867 	  }
1868 	}
1869       }
1870     }
1871     if (mapped_asv[m_index] & 4) { // mapped_hessians
1872       RealSymMatrix mapped_hess
1873         = mapped_response.function_hessian_view(m_index);
1874       for (j=0; j<num_mapped_deriv_vars; ++j) {
1875         for (k=0; k<=j; ++k) {
1876 	  if (identityRespMap)
1877 	    mapped_hess(j,k) = sub_iterator_hessians[m_index](j,k);
1878 	  else {
1879 	    Real& inner_prod = mapped_hess(j,k); inner_prod = 0.;
1880 	    for (l=0; l<numSubIterFns; ++l) {
1881 	      coeff = secondaryRespCoeffs(i,l);
1882 	      if (coeff != 0.) // avoid propagation of nan/inf for no mapping
1883 		inner_prod += coeff * sub_iterator_hessians[l](j,k);
1884 	    }
1885 	  }
1886 	}
1887       }
1888     }
1889   }
1890 }
1891 
1892 
1893 void NestedModel::
iterator_error_estimation(const RealVector & sub_iterator_errors,RealVector & mapped_errors)1894 iterator_error_estimation(const RealVector& sub_iterator_errors,
1895 			  RealVector& mapped_errors)
1896 {
1897   // In the future, could be overlaid with optional interface error estimates,
1898   // but for now, assume these are zero (e.g., deterministic mappings have no
1899   // estimator variance)
1900 
1901   if (sub_iterator_errors.empty()) {
1902     Cerr << "Error: sub_iterator_errors are undefined in NestedModel::"
1903 	 << "iterator_error_estimation().\n       Check error estimation "
1904 	 << "support in sub-method." << std::endl;
1905     abort_handler(MODEL_ERROR);
1906   }
1907 
1908   size_t i, j, m_index, num_mapped_fns = currentResponse.num_functions();
1909   if (static_cast<unsigned>(mapped_errors.length()) != num_mapped_fns)
1910     mapped_errors.size(num_mapped_fns); // init to 0
1911   else
1912     mapped_errors = 0.;
1913 
1914   // Assume independent Gaussian errors: then std error from linear combination
1915   // of Gaussian errors = Sqrt[ Sum [ coeff^2 sigma_i^2 ] ]
1916   // Note: final moments may be central or standard, but error estimates are
1917   //       always standard (sqrt of estimator variance of central/std moment)
1918 
1919   // [W]{S}:
1920   Real sum, term, coeff;
1921   for (i=0; i<subIterMappedPri; ++i) {
1922     if (identityRespMap)
1923       mapped_errors[i] = sub_iterator_errors[i];
1924     else {
1925       sum = 0.;
1926       for (j=0; j<numSubIterFns; ++j) {
1927 	coeff = primaryRespCoeffs(i,j);
1928 	if (coeff != 0.) { // avoid propagation of nan/inf for no mapping
1929 	  term = coeff * sub_iterator_errors[j]; // [W]{S}
1930 	  sum += term * term;
1931 	}
1932       }
1933       mapped_errors[i] = std::sqrt(sum);
1934     }
1935   }
1936 
1937   // [A]{S}:
1938   size_t num_mapped_1 = std::max(numOptInterfPrimary, subIterMappedPri);
1939   for (i=0; i<subIterMappedSec; ++i) {
1940     m_index = num_mapped_1 + numOptInterfIneqCon + i;// {a_l} <= [A]{S} <= {a_u}
1941     if (i>=numSubIterMappedIneqCon)
1942       m_index += numOptInterfEqCon;                           // [A]{S} == {a_t}
1943     if (identityRespMap)
1944       mapped_errors[m_index] = sub_iterator_errors[m_index];
1945     else {
1946       sum = 0.;
1947       for (j=0; j<numSubIterFns; ++j) {
1948 	coeff = secondaryRespCoeffs(i,j);
1949 	if (coeff != 0.) { // avoid propagation of nan/inf for no mapping
1950 	  term = coeff * sub_iterator_errors[j]; // [W]{S}
1951 	  sum += term * term;
1952 	}
1953       }
1954       mapped_errors[m_index] = std::sqrt(sum);
1955     }
1956   }
1957 }
1958 
1959 
check_response_map(const ShortArray & mapped_asv)1960 void NestedModel::check_response_map(const ShortArray& mapped_asv)
1961 {
1962   // counter initialization & sanity checking
1963   // NOTE: numSubIterFns != subIterMappedPri + subIterMappedSec
1964   // since subIterator response is converted to sub_iter_mapped_primary/sec
1965   // through the action of [W] and [A].
1966   size_t num_opt_interf_con = numOptInterfIneqCon + numOptInterfEqCon,
1967     num_mapped_1 = std::max(numOptInterfPrimary, subIterMappedPri);
1968   if (mapped_asv.size() !=
1969       num_mapped_1 + num_opt_interf_con + subIterMappedSec ||
1970       subIterMappedSec != numSubIterMappedIneqCon + numSubIterMappedEqCon) {
1971     Cerr << "\nError: bad function counts in NestedModel::check_response_map()."
1972          << std::endl;
1973     abort_handler(MODEL_ERROR);
1974   }
1975 }
1976 
1977 
component_parallel_mode(short mode)1978 void NestedModel::component_parallel_mode(short mode)
1979 {
1980   // mode may be correct, but can't guarantee active parallel config is in sync
1981   //if (componentParallelMode == mode)
1982   //  return; // already in correct parallel mode
1983 
1984   // terminate previous serve mode (if active)
1985   if (componentParallelMode != mode) {
1986     if (componentParallelMode == INTERFACE_MODE) {
1987       size_t index = subIteratorSched.miPLIndex;
1988       if (modelPCIter->mi_parallel_level_defined(index) &&
1989 	  modelPCIter->mi_parallel_level(index).server_communicator_size() > 1){
1990 	ParConfigLIter pc_iter = parallelLib.parallel_configuration_iterator();
1991 	parallelLib.parallel_configuration_iterator(modelPCIter);
1992 	optionalInterface.stop_evaluation_servers();
1993 	parallelLib.parallel_configuration_iterator(pc_iter); // restore
1994       }
1995     }
1996     // concurrent subIterator scheduling exits on its own (see IteratorScheduler
1997     // ::schedule_iterators(), but subModel eval scheduling is terminated here.
1998     else if (componentParallelMode == SUB_MODEL_MODE &&
1999 	     !subIteratorSched.messagePass) {
2000       ParConfigLIter pc_it = subModel.parallel_configuration_iterator();
2001       size_t index = subModel.mi_parallel_level_index();
2002       if (pc_it->mi_parallel_level_defined(index) &&
2003 	  pc_it->mi_parallel_level(index).server_communicator_size() > 1)
2004 	subModel.stop_servers();
2005     }
2006   }
2007 
2008   /* Moved up a level so that config can be restored after optInterface usage
2009   // set ParallelConfiguration for new mode
2010   if (mode == INTERFACE_MODE)
2011     parallelLib.parallel_configuration_iterator(modelPCIter);
2012   else if (mode == SUB_MODEL_MODE) {
2013     // ParallelLibrary::currPCIter activation delegated to subModel
2014   }
2015   */
2016 
2017   // activate new serve mode (matches NestedModel::serve_run(pl_iter)).  This
2018   // bcast matches the outer parallel context prior to subIterator partitioning.
2019   // > INTERFACE_MODE & subModel eval scheduling only broadcasts
2020   //   for mode change
2021   // > concurrent subIterator scheduling rebroadcasts every time since this
2022   //   scheduling exits on its own (see IteratorScheduler::schedule_iterators())
2023   if ( ( componentParallelMode != mode ||
2024 	 ( mode == SUB_MODEL_MODE && subIteratorSched.messagePass ) ) &&
2025        modelPCIter->mi_parallel_level_defined(outerMIPLIndex) ) {
2026     const ParallelLevel& mi_pl = modelPCIter->mi_parallel_level(outerMIPLIndex);
2027     if (mi_pl.server_communicator_size() > 1)
2028       parallelLib.bcast(mode, mi_pl);
2029   }
2030 
2031   componentParallelMode = mode;
2032 }
2033 
2034 
serve_run(ParLevLIter pl_iter,int max_eval_concurrency)2035 void NestedModel::serve_run(ParLevLIter pl_iter, int max_eval_concurrency)
2036 {
2037   // don't recurse, as subModel.serve() will set subModel comms
2038   set_communicators(pl_iter, max_eval_concurrency, false);
2039 
2040   // manage optionalInterface and subModel servers
2041   componentParallelMode = 1;
2042   while (componentParallelMode) {
2043     // outer context: matches bcast at bottom of component_parallel_mode()
2044     parallelLib.bcast(componentParallelMode, *pl_iter);
2045     if (componentParallelMode == INTERFACE_MODE &&
2046 	!optInterfacePointer.empty()) {
2047       // store/set/restore the ParallelLibrary::currPCIter
2048       ParConfigLIter pc_iter = parallelLib.parallel_configuration_iterator();
2049       parallelLib.parallel_configuration_iterator(modelPCIter);
2050       optionalInterface.serve_evaluations();
2051       parallelLib.parallel_configuration_iterator(pc_iter); // restore
2052     }
2053     else if (componentParallelMode == SUB_MODEL_MODE) {
2054       if (subIteratorSched.messagePass) // serve concurrent subIterator execs
2055 	subIteratorSched.schedule_iterators(*this, subIterator);
2056       else { // service the subModel for a single subIterator execution
2057 	ParLevLIter si_pl_iter // inner context
2058 	  = modelPCIter->mi_parallel_level_iterator(subIteratorSched.miPLIndex);
2059 	subModel.serve_run(si_pl_iter,
2060 			   subIterator.maximum_evaluation_concurrency());
2061       }
2062     }
2063   }
2064 }
2065 
2066 
update_inactive_view(short new_view,short & view)2067 void NestedModel::update_inactive_view(short new_view, short& view)
2068 {
2069   // Perform view aggregation/separation.  Additional logic and sanity
2070   // checking occurs within Variables::inactive_view() and Variables::
2071   // check_view_compatibility().
2072 
2073   if (new_view == RELAXED_ALL || new_view == MIXED_ALL) {
2074     // outer level has an ALL view --> infer what subset of this view should
2075     // be inner inactive by computing the complement of the inner active.
2076     // Can't use inactive types/ids, since inactive is either not defined
2077     // or not up to date.
2078     const Variables& sm_vars = subModel.current_variables();
2079     size_t i, num_sm_acv = sm_vars.acv(), num_sm_cv = sm_vars.cv(),
2080       sm_cv_start = sm_vars.cv_start();
2081     UShortMultiArrayConstView sm_acv_types
2082       = sm_vars.all_continuous_variable_types();
2083     bool relaxed = (new_view == RELAXED_ALL);
2084     // TO DO: THIS IS NOT GOOD ENOUGH SINCE THIS DOES NOT DISCRIMINATE INACTIVE
2085     // VARIABLES THAT ARE ACTUALLY BEING MAPPED FROM INACTIVE VARIABLES THAT
2086     // ARE MERELY PRESENT.
2087     for (i=0; i<num_sm_acv; ++i)
2088       if (i < sm_cv_start || i >= sm_cv_start+num_sm_cv) { // inactive
2089 	unsigned short type_i = sm_acv_types[i];
2090 	if (type_i == CONTINUOUS_DESIGN ||
2091 	    ( type_i >= DISCRETE_DESIGN_RANGE &&
2092 	      type_i <= DISCRETE_DESIGN_SET_REAL ) )
2093 	  view = (relaxed) ? RELAXED_DESIGN : MIXED_DESIGN;
2094 	else if ( type_i == CONTINUOUS_STATE ||
2095 		  ( type_i >= DISCRETE_STATE_RANGE &&
2096 		    type_i <= DISCRETE_STATE_SET_REAL ) )
2097 	  view = (relaxed) ? RELAXED_STATE : MIXED_STATE;
2098 	else if (type_i >= NORMAL_UNCERTAIN &&
2099 		 type_i <= DISCRETE_UNCERTAIN_SET_REAL)
2100 	  view = (relaxed) ? RELAXED_UNCERTAIN : MIXED_UNCERTAIN;
2101       }
2102   }
2103   else if (view == EMPTY_VIEW)
2104     view = new_view;
2105   else if (view != new_view) {
2106     // there are a few acceptable view promotions
2107     if ( ( view     == MIXED_ALEATORY_UNCERTAIN &&
2108 	   new_view == MIXED_EPISTEMIC_UNCERTAIN ) ||
2109 	 ( view     == MIXED_EPISTEMIC_UNCERTAIN &&
2110 	   new_view == MIXED_ALEATORY_UNCERTAIN ) )
2111       view = MIXED_UNCERTAIN; // aggregate
2112     else if ( ( view     == RELAXED_ALEATORY_UNCERTAIN &&
2113 		new_view == RELAXED_EPISTEMIC_UNCERTAIN ) ||
2114 	      ( view     == RELAXED_EPISTEMIC_UNCERTAIN &&
2115 		new_view == RELAXED_ALEATORY_UNCERTAIN ) )
2116       view = RELAXED_UNCERTAIN; // aggregate
2117     else if ( ( view == MIXED_UNCERTAIN &&
2118 		( new_view == MIXED_ALEATORY_UNCERTAIN ||
2119 		  new_view == MIXED_EPISTEMIC_UNCERTAIN ) ) ||
2120 	      ( view == RELAXED_UNCERTAIN &&
2121 		( new_view == RELAXED_ALEATORY_UNCERTAIN ||
2122 		  new_view == RELAXED_EPISTEMIC_UNCERTAIN ) ) )
2123       { } // already aggregated
2124     else {
2125       Cerr << "\nError: inactive sub-model view discrepancy in NestedModel::"
2126 	   << "update_inactive_view()." << std::endl;
2127       abort_handler(MODEL_ERROR);
2128     }
2129   }
2130 }
2131 
2132 
update_inactive_view(unsigned short type,short & view)2133 void NestedModel::update_inactive_view(unsigned short type, short& view)
2134 {
2135   // determine RELAXED or MIXED primary view at sub-model level
2136   short new_view, active_sm_view = subModel.current_variables().view().first;
2137   bool relaxed = ( active_sm_view == RELAXED_ALL ||
2138 		   ( active_sm_view >= RELAXED_DESIGN &&
2139 		     active_sm_view <= RELAXED_STATE ) );
2140 
2141   if (type >= CONTINUOUS_DESIGN && type <= DISCRETE_DESIGN_SET_REAL) {
2142     new_view = (relaxed) ? RELAXED_DESIGN : MIXED_DESIGN;
2143     update_inactive_view(new_view, view);
2144   }
2145   else if (type >= CONTINUOUS_STATE && type <= DISCRETE_STATE_SET_REAL) {
2146     new_view = (relaxed) ? RELAXED_STATE : MIXED_STATE;
2147     update_inactive_view(new_view, view);
2148   }
2149   else if (type >= NORMAL_UNCERTAIN && type <= DISCRETE_UNCERTAIN_SET_REAL) {
2150     if (type >= CONTINUOUS_INTERVAL_UNCERTAIN)
2151       new_view = (relaxed) ? RELAXED_EPISTEMIC_UNCERTAIN :
2152 	                     MIXED_EPISTEMIC_UNCERTAIN;
2153     else
2154       new_view = (relaxed) ? RELAXED_ALEATORY_UNCERTAIN :
2155 			     MIXED_ALEATORY_UNCERTAIN;
2156     update_inactive_view(new_view, view);
2157   }
2158 }
2159 
2160 
2161 void NestedModel::
update_sub_model(const Variables & vars,const Constraints & cons)2162 update_sub_model(const Variables& vars, const Constraints& cons)
2163 {
2164   // Update subModel variables using active currentVariables through a
2165   // combination of variable insertions and augmentations.  Insertions and
2166   // augmentations are not mutually exclusive, so subModel updates must
2167   // account for both, potentially on a variable by variable basis.
2168 
2169   // Top-level active variables/bounds/labels are mapped to the subModel, but
2170   // top-level linear/nonlinear constraint coefficients/bounds/targets are not
2171   // propagated (as they are in the surrogate model classes) since they are not
2172   // relevant at the subModel level.
2173 
2174   // For defined variable mappings, insert active currentVariables into subModel
2175   // variables (1). When secondary mappings are defined, this insertion involves
2176   // the updating of sub-parameters for the subModel variables (2). When only
2177   // primary mappings are defined, the insertions update the subModel variable
2178   // values. For null primary mappings (empty strings), augmentations are
2179   // performed by updating the subModel variables of the same variable type as
2180   // the active currentVariables (3). Bounds, labels, and sub-parameters are
2181   // additionally propagated for augmentations and active complement mappings
2182   // (both based on matching consistent types), with some differences in
2183   // required updating frequency.
2184 
2185   // Distribution parameter insertion is used for second-order probability
2186   // (UQ sample inserted into distrib parameters for UQ reliability/sampling)
2187   // and for OUU with design/uncertain overlap (optimization design variables
2188   // inserted into distr. parameters for UQ reliability/sampling).  For now,
2189   // existing distribution parameters are supported for the secondary variable
2190   // mappings, and a few variables additionally support LOCATION and SCALE
2191   // (this can be expanded in the future).
2192 
2193   size_t i, curr_i, num_var_map_2 = active2ACVarMapTargets.size(),
2194     num_curr_cv  = vars.cv(),  num_curr_div = vars.div(),
2195     num_curr_dsv = vars.dsv(), num_curr_drv = vars.drv(),
2196     pacvm_index, padivm_index, padsvm_index, padrvm_index;
2197   const SharedVariablesData& svd = vars.shared_data();
2198   const SharedVariablesData& sm_svd
2199     = subModel.current_variables().shared_data();
2200   std::shared_ptr<Pecos::MarginalsCorrDistribution> sm_mvd_rep =
2201     std::static_pointer_cast<Pecos::MarginalsCorrDistribution>
2202     (subModel.multivariate_distribution().multivar_dist_rep());
2203 
2204   // Map ACTIVE CONTINUOUS VARIABLES from currentVariables
2205   if (num_curr_cv) {
2206     const RealVector& curr_c_vars   = vars.continuous_variables();
2207     const RealVector& curr_c_l_bnds = cons.continuous_lower_bounds();
2208     const RealVector& curr_c_u_bnds = cons.continuous_upper_bounds();
2209     StringMultiArrayConstView curr_c_labels
2210       = vars.continuous_variable_labels();
2211     for (i=0; i<num_curr_cv; ++i) {
2212       curr_i = svd.cv_index_to_active_index(i);
2213       pacvm_index  =  active1ACVarMapIndices[curr_i];
2214       padivm_index = active1ADIVarMapIndices[curr_i];
2215       padsvm_index = active1ADSVarMapIndices[curr_i];
2216       padrvm_index = active1ADRVarMapIndices[curr_i];
2217       if (pacvm_index != _NPOS) {
2218 	short sacvm_target
2219 	  = (num_var_map_2) ? active2ACVarMapTargets[curr_i] : Pecos::NO_TARGET;
2220 	if (sacvm_target == Pecos::NO_TARGET) {
2221 	  subModel.all_continuous_variable(curr_c_vars[i], pacvm_index);
2222 	  if (extraCVarsData[i]) { // default mapping between consistent types
2223 	    subModel.all_continuous_lower_bound(curr_c_l_bnds[i], pacvm_index);
2224 	    subModel.all_continuous_upper_bound(curr_c_u_bnds[i], pacvm_index);
2225 	    // Note: this is more general than just bounds (all dist params):
2226 	    sm_mvd_rep->pull_distribution_parameters(mvDist,
2227 	      svd.cv_index_to_all_index(i),
2228 	      sm_svd.acv_index_to_all_index(pacvm_index));
2229 	    if (firstUpdate)
2230 	      subModel.all_continuous_variable_label(curr_c_labels[i],
2231 						     pacvm_index);
2232 	  }
2233 	}
2234 	else
2235 	  real_variable_mapping(curr_c_vars[i], pacvm_index, sacvm_target);
2236       }
2237       else if (padivm_index != _NPOS) {
2238 	short sadivm_target = (num_var_map_2) ? active2ADIVarMapTargets[curr_i]
2239 	  : Pecos::NO_TARGET;
2240 	real_variable_mapping(curr_c_vars[i], padivm_index, sadivm_target);
2241       }
2242       else if (padsvm_index != _NPOS) {
2243 	short sadsvm_target = (num_var_map_2) ? active2ADSVarMapTargets[curr_i]
2244 	  : Pecos::NO_TARGET;
2245 	real_variable_mapping(curr_c_vars[i], padsvm_index, sadsvm_target);
2246       }
2247       else if (padrvm_index != _NPOS) {
2248 	short sadrvm_target = (num_var_map_2) ? active2ADRVarMapTargets[curr_i]
2249 	  : Pecos::NO_TARGET;
2250 	real_variable_mapping(curr_c_vars[i], padrvm_index, sadrvm_target);
2251       }
2252     }
2253   }
2254 
2255   // Map ACTIVE DISCRETE INTEGER VARIABLES from currentVariables
2256   if (num_curr_div) {
2257     const IntVector& curr_di_vars   = vars.discrete_int_variables();
2258     const IntVector& curr_di_l_bnds = cons.discrete_int_lower_bounds();
2259     const IntVector& curr_di_u_bnds = cons.discrete_int_upper_bounds();
2260     StringMultiArrayConstView curr_di_labels
2261       = vars.discrete_int_variable_labels();
2262     for (i=0; i<num_curr_div; ++i) {
2263       curr_i = svd.div_index_to_active_index(i);
2264       pacvm_index  =  active1ACVarMapIndices[curr_i];
2265       padivm_index = active1ADIVarMapIndices[curr_i];
2266       padsvm_index = active1ADSVarMapIndices[curr_i];
2267       padrvm_index = active1ADRVarMapIndices[curr_i];
2268       if (pacvm_index != _NPOS) {
2269 	short sacvm_target
2270 	  = (num_var_map_2) ? active2ACVarMapTargets[curr_i] : Pecos::NO_TARGET;
2271 	integer_variable_mapping(curr_di_vars[i], pacvm_index, sacvm_target);
2272       }
2273       else if (padivm_index != _NPOS) {
2274 	short sadivm_target = (num_var_map_2) ?
2275 	  active2ADIVarMapTargets[curr_i] : Pecos::NO_TARGET;
2276 	if (sadivm_target == Pecos::NO_TARGET) {
2277 	  subModel.all_discrete_int_variable(curr_di_vars[i], padivm_index);
2278 	  if (extraDIVarsData[i]) { // default mapping between consistent types
2279 	    subModel.all_discrete_int_lower_bound(curr_di_l_bnds[i],
2280 						  padivm_index);
2281 	    subModel.all_discrete_int_upper_bound(curr_di_u_bnds[i],
2282 						  padivm_index);
2283 	    // Note: this is more general than just bounds (all dist params):
2284 	    sm_mvd_rep->pull_distribution_parameters(mvDist,
2285 	      svd.div_index_to_all_index(i),
2286 	      sm_svd.adiv_index_to_all_index(padivm_index));
2287 	    if (firstUpdate)
2288 	      subModel.all_discrete_int_variable_label(curr_di_labels[i],
2289 						       padivm_index);
2290 	  }
2291 	}
2292 	else
2293 	  integer_variable_mapping(curr_di_vars[i], padivm_index,sadivm_target);
2294       }
2295       else if (padsvm_index != _NPOS) {
2296 	short sadsvm_target = (num_var_map_2) ?
2297 	  active2ADSVarMapTargets[curr_i] : Pecos::NO_TARGET;
2298 	integer_variable_mapping(curr_di_vars[i], padsvm_index, sadsvm_target);
2299       }
2300       else if (padrvm_index != _NPOS) {
2301 	short sadrvm_target = (num_var_map_2) ?
2302 	  active2ADRVarMapTargets[curr_i] : Pecos::NO_TARGET;
2303 	integer_variable_mapping(curr_di_vars[i], padrvm_index, sadrvm_target);
2304       }
2305     }
2306   }
2307 
2308   // Map ACTIVE DISCRETE STRING VARIABLES from currentVariables
2309   if (num_curr_dsv) {
2310     StringMultiArrayConstView curr_ds_vars = vars.discrete_string_variables();
2311     StringMultiArrayConstView curr_ds_labels
2312       = vars.discrete_string_variable_labels();
2313     for (i=0; i<num_curr_dsv; ++i) {
2314       curr_i = svd.dsv_index_to_active_index(i);
2315       pacvm_index  =  active1ACVarMapIndices[curr_i];
2316       padivm_index = active1ADIVarMapIndices[curr_i];
2317       padsvm_index = active1ADSVarMapIndices[curr_i];
2318       padrvm_index = active1ADRVarMapIndices[curr_i];
2319       if (pacvm_index != _NPOS) {
2320 	short sacvm_target = (num_var_map_2) ?
2321 	  active2ACVarMapTargets[curr_i] : Pecos::NO_TARGET;
2322 	string_variable_mapping(curr_ds_vars[i], pacvm_index, sacvm_target);
2323       }
2324       else if (padivm_index != _NPOS) {
2325 	short sadivm_target = (num_var_map_2) ?
2326 	  active2ADIVarMapTargets[curr_i] : Pecos::NO_TARGET;
2327 	string_variable_mapping(curr_ds_vars[i], padivm_index, sadivm_target);
2328       }
2329       else if (padsvm_index != _NPOS) {
2330 	short sadsvm_target = (num_var_map_2) ?
2331 	  active2ADSVarMapTargets[curr_i] : Pecos::NO_TARGET;
2332 	if (sadsvm_target == Pecos::NO_TARGET) {
2333 	  subModel.all_discrete_string_variable(curr_ds_vars[i], padsvm_index);
2334 	  if (extraDSVarsData[i]) { // default mapping between consistent types
2335 	    // Note: this is more general than just bounds (all dist params):
2336 	    sm_mvd_rep->pull_distribution_parameters(mvDist,
2337 	      svd.dsv_index_to_all_index(i),
2338 	      sm_svd.adsv_index_to_all_index(padsvm_index));
2339 	    if (firstUpdate)
2340 	      subModel.all_discrete_string_variable_label(curr_ds_labels[i],
2341 							  padsvm_index);
2342 	  }
2343 	}
2344 	else
2345 	  string_variable_mapping(curr_ds_vars[i], padsvm_index, sadsvm_target);
2346       }
2347       else if (padrvm_index != _NPOS) {
2348 	short sadrvm_target = (num_var_map_2) ?
2349 	  active2ADRVarMapTargets[curr_i] : Pecos::NO_TARGET;
2350 	string_variable_mapping(curr_ds_vars[i], padrvm_index, sadrvm_target);
2351       }
2352     }
2353   }
2354 
2355   // Map ACTIVE DISCRETE REAL VARIABLES from currentVariables
2356   if (num_curr_drv) {
2357     const RealVector& curr_dr_vars   = vars.discrete_real_variables();
2358     const RealVector& curr_dr_l_bnds = cons.discrete_real_lower_bounds();
2359     const RealVector& curr_dr_u_bnds = cons.discrete_real_upper_bounds();
2360     StringMultiArrayConstView curr_dr_labels
2361       = vars.discrete_real_variable_labels();
2362     for (i=0; i<num_curr_drv; ++i) {
2363       curr_i = svd.drv_index_to_active_index(i);
2364       pacvm_index  =  active1ACVarMapIndices[curr_i];
2365       padivm_index = active1ADIVarMapIndices[curr_i];
2366       padsvm_index = active1ADSVarMapIndices[curr_i];
2367       padrvm_index = active1ADRVarMapIndices[curr_i];
2368       if (pacvm_index != _NPOS) {
2369 	short sacvm_target = (num_var_map_2) ?
2370 	  active2ACVarMapTargets[curr_i] : Pecos::NO_TARGET;
2371 	real_variable_mapping(curr_dr_vars[i], pacvm_index, sacvm_target);
2372       }
2373       else if (padivm_index != _NPOS) {
2374 	short sadivm_target = (num_var_map_2) ?
2375 	  active2ADIVarMapTargets[curr_i] : Pecos::NO_TARGET;
2376 	real_variable_mapping(curr_dr_vars[i], padivm_index, sadivm_target);
2377       }
2378       else if (padsvm_index != _NPOS) {
2379 	short sadsvm_target = (num_var_map_2) ?
2380 	  active2ADSVarMapTargets[curr_i] : Pecos::NO_TARGET;
2381 	real_variable_mapping(curr_dr_vars[i], padsvm_index, sadsvm_target);
2382       }
2383       else if (padrvm_index != _NPOS) {
2384 	short sadrvm_target = (num_var_map_2) ?
2385 	  active2ADRVarMapTargets[curr_i] : Pecos::NO_TARGET;
2386 	if (sadrvm_target == Pecos::NO_TARGET) {
2387 	  subModel.all_discrete_real_variable(curr_dr_vars[i], padrvm_index);
2388 	  if (extraDRVarsData[i]) { // default mapping between consistent types
2389 	    subModel.all_discrete_real_lower_bound(curr_dr_l_bnds[i],
2390 						   padrvm_index);
2391 	    subModel.all_discrete_real_upper_bound(curr_dr_u_bnds[i],
2392 						   padrvm_index);
2393 	    // Note: this is more general than just bounds (all dist params):
2394 	    sm_mvd_rep->pull_distribution_parameters(mvDist,
2395 	      svd.drv_index_to_all_index(i),
2396 	      sm_svd.adrv_index_to_all_index(padrvm_index));
2397 	    if (firstUpdate)
2398 	      subModel.all_discrete_real_variable_label(curr_dr_labels[i],
2399 							padrvm_index);
2400 	  }
2401 	}
2402 	else
2403 	  real_variable_mapping(curr_dr_vars[i], padrvm_index, sadrvm_target);
2404       }
2405     }
2406   }
2407 
2408   // Map COMPLEMENT CONTINUOUS VARIABLES from currentVariables into
2409   // corresponding subModel type (using same logic as default active mapping)
2410   size_t num_curr_ccv = vars.acv() - num_curr_cv, c1_index;
2411   if (num_curr_ccv) {
2412     const RealVector& curr_ac_vars   = vars.all_continuous_variables();
2413     const RealVector& curr_ac_l_bnds = cons.all_continuous_lower_bounds();
2414     const RealVector& curr_ac_u_bnds = cons.all_continuous_upper_bounds();
2415     StringMultiArrayConstView curr_ac_labels
2416       = vars.all_continuous_variable_labels();
2417     for (i=0; i<num_curr_ccv; ++i) {
2418       curr_i = svd.ccv_index_to_acv_index(i);
2419       c1_index = complement1ACVarMapIndices[i];
2420       subModel.all_continuous_variable(curr_ac_vars[curr_i], c1_index);
2421       subModel.all_continuous_lower_bound(curr_ac_l_bnds[curr_i], c1_index);
2422       subModel.all_continuous_upper_bound(curr_ac_u_bnds[curr_i], c1_index);
2423       if (firstUpdate) {
2424 	subModel.all_continuous_variable_label(curr_ac_labels[curr_i],c1_index);
2425 	// Note: this is more general than just bounds (all dist params):
2426 	sm_mvd_rep->pull_distribution_parameters(mvDist,
2427 	  svd.ccv_index_to_all_index(i),
2428 	  sm_svd.acv_index_to_all_index(c1_index));
2429       }
2430     }
2431   }
2432 
2433   // Map COMPLEMENT DISCRETE INTEGER VARIABLES from currentVariables into
2434   // corresponding subModel type (using same logic as default active mapping)
2435   size_t num_curr_cdiv = vars.adiv() - num_curr_div;
2436   if (num_curr_cdiv) {
2437     const IntVector& curr_adi_vars   = vars.all_discrete_int_variables();
2438     const IntVector& curr_adi_l_bnds = cons.all_discrete_int_lower_bounds();
2439     const IntVector& curr_adi_u_bnds = cons.all_discrete_int_upper_bounds();
2440     StringMultiArrayConstView curr_adi_labels
2441       = vars.all_discrete_int_variable_labels();
2442     for (i=0; i<num_curr_cdiv; ++i) {
2443       curr_i = svd.cdiv_index_to_adiv_index(i);
2444       c1_index = complement1ADIVarMapIndices[i];
2445       subModel.all_discrete_int_variable(curr_adi_vars[curr_i], c1_index);
2446       subModel.all_discrete_int_lower_bound(curr_adi_l_bnds[curr_i], c1_index);
2447       subModel.all_discrete_int_upper_bound(curr_adi_u_bnds[curr_i], c1_index);
2448       if (firstUpdate) {
2449 	subModel.all_discrete_int_variable_label(curr_adi_labels[curr_i],
2450 						 c1_index);
2451 	// Note: this is more general than just bounds (all dist params):
2452 	sm_mvd_rep->pull_distribution_parameters(mvDist,
2453 	  svd.cdiv_index_to_all_index(i),
2454 	  sm_svd.adiv_index_to_all_index(c1_index));
2455       }
2456     }
2457   }
2458 
2459   // Map COMPLEMENT DISCRETE STRING VARIABLES from currentVariables into
2460   // corresponding subModel type (using same logic as default active mapping)
2461   size_t num_curr_cdsv = vars.adsv() - num_curr_dsv;
2462   if (num_curr_cdsv) {
2463     StringMultiArrayConstView curr_ads_vars
2464       = vars.all_discrete_string_variables();
2465     StringMultiArrayConstView curr_ads_labels
2466       = vars.all_discrete_string_variable_labels();
2467     for (i=0; i<num_curr_cdsv; ++i) {
2468       curr_i = svd.cdsv_index_to_adsv_index(i);
2469       c1_index = complement1ADSVarMapIndices[i];
2470       subModel.all_discrete_string_variable(curr_ads_vars[curr_i], c1_index);
2471       if (firstUpdate) {
2472 	subModel.all_discrete_string_variable_label(curr_ads_labels[curr_i],
2473 						    c1_index);
2474 	// Note: this is more general than just bounds (all dist params):
2475 	sm_mvd_rep->pull_distribution_parameters(mvDist,
2476 	  svd.cdsv_index_to_all_index(i),
2477 	  sm_svd.adsv_index_to_all_index(c1_index));
2478       }
2479     }
2480   }
2481 
2482   // Map COMPLEMENT DISCRETE REAL VARIABLES from currentVariables into
2483   // corresponding subModel type (using same logic as default active mapping)
2484   size_t num_curr_cdrv = vars.adrv() - num_curr_drv;
2485   if (num_curr_cdrv) {
2486     const RealVector& curr_adr_vars   = vars.all_discrete_real_variables();
2487     const RealVector& curr_adr_l_bnds = cons.all_discrete_real_lower_bounds();
2488     const RealVector& curr_adr_u_bnds = cons.all_discrete_real_upper_bounds();
2489     StringMultiArrayConstView curr_adr_labels
2490       = vars.all_discrete_real_variable_labels();
2491     for (i=0; i<num_curr_cdrv; ++i) {
2492       curr_i = svd.cdrv_index_to_adrv_index(i);
2493       c1_index = complement1ADRVarMapIndices[i];
2494       subModel.all_discrete_real_variable(curr_adr_vars[curr_i], c1_index);
2495       subModel.all_discrete_real_lower_bound(curr_adr_l_bnds[curr_i], c1_index);
2496       subModel.all_discrete_real_upper_bound(curr_adr_u_bnds[curr_i], c1_index);
2497       if (firstUpdate) {
2498 	subModel.all_discrete_real_variable_label(curr_adr_labels[curr_i],
2499 						  c1_index);
2500 	// Note: this is more general than just bounds (all dist params):
2501 	sm_mvd_rep->pull_distribution_parameters(mvDist,
2502 	  svd.cdrv_index_to_all_index(i),
2503 	  sm_svd.adrv_index_to_all_index(c1_index));
2504       }
2505     }
2506   }
2507 
2508   firstUpdate = false;
2509 }
2510 
2511 
2512 void NestedModel::
real_variable_mapping(Real r_var,size_t av_index,short svm_target)2513 real_variable_mapping(Real r_var, size_t av_index, short svm_target)
2514 {
2515   Pecos::MultivariateDistribution& sm_mvd
2516     = subModel.multivariate_distribution();
2517   std::shared_ptr<Pecos::MarginalsCorrDistribution> sm_mvd_rep =
2518     std::static_pointer_cast<Pecos::MarginalsCorrDistribution>
2519     (sm_mvd.multivar_dist_rep());
2520 
2521   const SharedVariablesData& sm_svd
2522     = subModel.current_variables().shared_data();
2523 
2524   switch (svm_target) {
2525   case Pecos::CR_LWR_BND:  case Pecos::N_LWR_BND:  case Pecos::LN_LWR_BND:
2526   case Pecos::U_LWR_BND:   case Pecos::LU_LWR_BND: case Pecos::T_LWR_BND:
2527   case Pecos::BE_LWR_BND:
2528     sm_mvd_rep->push_parameter(sm_svd.acv_index_to_all_index(av_index),
2529 			       svm_target, r_var);
2530     subModel.all_continuous_lower_bound(r_var, av_index);
2531     break;
2532   case Pecos::CR_UPR_BND:  case Pecos::N_UPR_BND:  case Pecos::LN_UPR_BND:
2533   case Pecos::U_UPR_BND:   case Pecos::LU_UPR_BND: case Pecos::T_UPR_BND:
2534   case Pecos::BE_UPR_BND:
2535     sm_mvd_rep->push_parameter(sm_svd.acv_index_to_all_index(av_index),
2536 			       svm_target, r_var);
2537     subModel.all_continuous_upper_bound(r_var, av_index);
2538     break;
2539   case Pecos::N_MEAN:      case Pecos::N_STD_DEV:  case Pecos::LN_MEAN:
2540   case Pecos::LN_STD_DEV:  case Pecos::LN_LAMBDA:  case Pecos::LN_ZETA:
2541   case Pecos::LN_ERR_FACT: case Pecos::T_MODE:     case Pecos::E_BETA:
2542   case Pecos::BE_ALPHA:    case Pecos::BE_BETA:    case Pecos::GA_ALPHA:
2543   case Pecos::GA_BETA:     case Pecos::GU_ALPHA:   case Pecos::GU_BETA:
2544   case Pecos::F_ALPHA:     case Pecos::F_BETA:     case Pecos::W_ALPHA:
2545   case Pecos::W_BETA:
2546     sm_mvd_rep->push_parameter(sm_svd.acv_index_to_all_index(av_index),
2547 			       svm_target, r_var);
2548     break;
2549   case Pecos::P_LAMBDA:        case Pecos::BI_P_PER_TRIAL:
2550   case Pecos::NBI_P_PER_TRIAL: case Pecos::GE_P_PER_TRIAL:
2551     sm_mvd_rep->push_parameter(sm_svd.adiv_index_to_all_index(av_index),
2552 			       svm_target, r_var);
2553     break;
2554   // N_{MEAN,STD_DEV,LWR_BND,UPR_BND} change individual dist parameters only.
2555   // N_{LOCATION,SCALE} change multiple parameters to accomplish a translation
2556   // or scaling.  They are mapped to a convenient definition believed to be the
2557   // most natural for a user, where the definition changes from distribution to
2558   // distribution (a consistent meaning of mu/sigma would be more awkward for a
2559   // user to convert).  For Normal, location & scale are mean & std deviation.
2560   case Pecos::N_LOCATION: { // a translation with no change in shape/scale
2561     size_t rv_index = sm_svd.acv_index_to_all_index(av_index);
2562     Real mean, l_bnd, u_bnd;
2563     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::N_MEAN, mean);
2564     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::N_LWR_BND, l_bnd);
2565     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::N_UPR_BND, u_bnd);
2566     Real delta = r_var - mean;
2567     // translate: change bounds by same amount as mean
2568     sm_mvd_rep->push_parameter(rv_index, Pecos::N_MEAN, r_var);
2569     Real dbl_inf = std::numeric_limits<Real>::infinity();
2570     if (l_bnd > -dbl_inf) {
2571       Real new_l_bnd = l_bnd + delta;
2572       sm_mvd_rep->push_parameter(rv_index, Pecos::N_LWR_BND, new_l_bnd);
2573       subModel.all_continuous_lower_bound(new_l_bnd, av_index);
2574     }
2575     if (u_bnd <  dbl_inf) {
2576       Real new_u_bnd = u_bnd + delta;
2577       sm_mvd_rep->push_parameter(rv_index, Pecos::N_UPR_BND, new_u_bnd);
2578       subModel.all_continuous_upper_bound(new_u_bnd, av_index);
2579     }
2580     break;
2581   }
2582   case Pecos::N_SCALE: { // change in shape/scale without translation
2583     size_t rv_index = sm_svd.acv_index_to_all_index(av_index);
2584     Real mean, stdev, l_bnd, u_bnd;
2585     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::N_MEAN, mean);
2586     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::N_STD_DEV, stdev);
2587     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::N_LWR_BND, l_bnd);
2588     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::N_UPR_BND, u_bnd);
2589     // scale: preserve number of std deviations where l,u bound occurs
2590     sm_mvd_rep->push_parameter(rv_index, Pecos::N_STD_DEV, r_var);
2591     Real dbl_inf = std::numeric_limits<Real>::infinity();
2592     if (l_bnd > -dbl_inf) {
2593       Real num_sig_l = (mean - l_bnd) / stdev,
2594 	   new_l_bnd = mean - num_sig_l * r_var;
2595       sm_mvd_rep->push_parameter(rv_index, Pecos::N_LWR_BND, new_l_bnd);
2596       subModel.all_continuous_lower_bound(new_l_bnd, av_index);
2597     }
2598     if (u_bnd <  dbl_inf) {
2599       Real num_sig_u = (u_bnd - mean) / stdev,
2600 	   new_u_bnd = mean + num_sig_u * r_var;
2601       sm_mvd_rep->push_parameter(rv_index, Pecos::N_UPR_BND, new_u_bnd);
2602       subModel.all_continuous_upper_bound(new_u_bnd, av_index);
2603     }
2604     break;
2605   }
2606   // U_{LWR_BND,UPR_BND} change individual dist parameters only.
2607   // U_{LOCATION,SCALE} change multiple parameters to accomplish a translation
2608   // or scaling.  They are mapped to a convenient definition believed to be the
2609   // most natural for a user, where the definition changes from distribution to
2610   // distribution (a consistent meaning of mu/sigma would be more awkward for a
2611   // user to convert).  For Uniform, location & scale are center & range.
2612   case Pecos::U_LOCATION: {
2613     size_t rv_index = sm_svd.acv_index_to_all_index(av_index);
2614     // translate: change both bounds by same amount
2615     Real l_bnd, u_bnd;
2616     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::U_LWR_BND, l_bnd);
2617     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::U_UPR_BND, u_bnd);
2618     Real center = (u_bnd + l_bnd) / 2., delta = r_var - center,
2619       new_l_bnd = l_bnd + delta, new_u_bnd = u_bnd + delta;
2620     sm_mvd_rep->push_parameter(rv_index, Pecos::U_LWR_BND, new_l_bnd);
2621     sm_mvd_rep->push_parameter(rv_index, Pecos::U_UPR_BND, new_u_bnd);
2622     subModel.all_continuous_lower_bound(new_l_bnd, av_index);
2623     subModel.all_continuous_upper_bound(new_u_bnd, av_index);
2624     break;
2625   }
2626   case Pecos::U_SCALE: {
2627     size_t rv_index = sm_svd.acv_index_to_all_index(av_index);
2628     // scale: move bounds in/out by same amount about consistent center
2629     Real l_bnd, u_bnd;
2630     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::U_LWR_BND, l_bnd);
2631     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::U_UPR_BND, u_bnd);
2632     Real center = (u_bnd + l_bnd) / 2., half_range = r_var / 2.,
2633       new_l_bnd = center-half_range, new_u_bnd = center+half_range;
2634     sm_mvd_rep->push_parameter(rv_index, Pecos::U_LWR_BND, new_l_bnd);
2635     sm_mvd_rep->push_parameter(rv_index, Pecos::U_UPR_BND, new_u_bnd);
2636     subModel.all_continuous_lower_bound(new_l_bnd, av_index);
2637     subModel.all_continuous_upper_bound(new_u_bnd, av_index);
2638     break;
2639   }
2640   // T_{MODE,LWR_BND,UPR_BND} change individual dist parameters only.
2641   // T_{LOCATION,SCALE} change multiple parameters to accomplish a translation
2642   // or scaling.  They are mapped to a convenient definition believed to be the
2643   // most natural for a user, where the definition changes from distribution to
2644   // distribution (a consistent meaning of mu/sigma would be more awkward for a
2645   // user to convert).  For Triangular, location & scale are mode & range.
2646   case Pecos::T_LOCATION: {
2647     size_t rv_index = sm_svd.acv_index_to_all_index(av_index);
2648     // translate: change mode and both bounds by same amount
2649     Real mode, l_bnd, u_bnd;
2650     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::T_MODE, mode);
2651     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::T_LWR_BND, l_bnd);
2652     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::T_UPR_BND, u_bnd);
2653     Real  delta = r_var - mode, new_l_bnd = l_bnd + delta,
2654       new_u_bnd = u_bnd + delta;
2655     sm_mvd_rep->push_parameter(rv_index, Pecos::T_MODE,    r_var);
2656     sm_mvd_rep->push_parameter(rv_index, Pecos::T_LWR_BND, new_l_bnd);
2657     sm_mvd_rep->push_parameter(rv_index, Pecos::T_UPR_BND, new_u_bnd);
2658     subModel.all_continuous_lower_bound(new_l_bnd, av_index);
2659     subModel.all_continuous_upper_bound(new_u_bnd, av_index);
2660     break;
2661   }
2662   case Pecos::T_SCALE: {
2663     size_t rv_index = sm_svd.acv_index_to_all_index(av_index);
2664     // scale: preserve L/M/U proportions while scaling range
2665     Real mode, l_bnd, u_bnd;
2666     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::T_MODE, mode);
2667     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::T_LWR_BND, l_bnd);
2668     sm_mvd_rep->pull_parameter<Real>(rv_index, Pecos::T_UPR_BND, u_bnd);
2669     Real  range = u_bnd - l_bnd, perc_l = (mode - l_bnd) / range,
2670          perc_u = (u_bnd - mode) / range, new_l_bnd = mode - perc_l * r_var,
2671       new_u_bnd = mode + perc_u * r_var;
2672     sm_mvd_rep->push_parameter(rv_index, Pecos::T_LWR_BND, new_l_bnd);
2673     sm_mvd_rep->push_parameter(rv_index, Pecos::T_UPR_BND, new_u_bnd);
2674     subModel.all_continuous_lower_bound(new_l_bnd, av_index);
2675     subModel.all_continuous_upper_bound(new_u_bnd, av_index);
2676     break;
2677   }
2678   case Pecos::NO_TARGET: default:
2679     Cerr << "\nError: secondary mapping target unmatched for real value "
2680 	 << "insertion in NestedModel::real_variable_mapping()." << std::endl;
2681     abort_handler(MODEL_ERROR);
2682   }
2683 }
2684 
2685 
2686 void NestedModel::
integer_variable_mapping(int i_var,size_t av_index,short svm_target)2687 integer_variable_mapping(int i_var, size_t av_index, short svm_target)
2688 {
2689   Pecos::MultivariateDistribution& sm_mvd
2690     = subModel.multivariate_distribution();
2691   std::shared_ptr<Pecos::MarginalsCorrDistribution> sm_mvd_rep =
2692     std::static_pointer_cast<Pecos::MarginalsCorrDistribution>
2693     (sm_mvd.multivar_dist_rep());
2694 
2695   const SharedVariablesData& sm_svd
2696     = subModel.current_variables().shared_data();
2697 
2698   switch (svm_target) {
2699   case Pecos::DR_LWR_BND:
2700     sm_mvd_rep->push_parameter(sm_svd.adiv_index_to_all_index(av_index),
2701 			       svm_target, i_var);
2702     subModel.all_discrete_int_lower_bound(i_var, av_index);
2703     break;
2704   case Pecos::DR_UPR_BND:
2705     sm_mvd_rep->push_parameter(sm_svd.adiv_index_to_all_index(av_index),
2706 			       svm_target, i_var);
2707     subModel.all_discrete_int_upper_bound(i_var, av_index);
2708     break;
2709   case Pecos::BI_TRIALS:    case Pecos::NBI_TRIALS:
2710   case Pecos::HGE_TOT_POP:  case Pecos::HGE_SEL_POP:  case Pecos::HGE_DRAWN: {
2711     unsigned int ui_var = (unsigned int)i_var;
2712     sm_mvd_rep->push_parameter(sm_svd.adiv_index_to_all_index(av_index),
2713 			       svm_target, ui_var);
2714     break;
2715   }
2716   case Pecos::NO_TARGET: default:
2717     Cerr << "\nError: secondary mapping target unmatched for integer value "
2718 	 << "insertion in NestedModel::integer_variable_mapping()" << std::endl;
2719     abort_handler(MODEL_ERROR);
2720   }
2721 }
2722 
2723 
2724 void NestedModel::
string_variable_mapping(const String & s_var,size_t av_index,short svm_target)2725 string_variable_mapping(const String& s_var, size_t av_index,
2726 			short svm_target)
2727 {
2728   Pecos::MultivariateDistribution& sm_mvd
2729     = subModel.multivariate_distribution();
2730   std::shared_ptr<Pecos::MarginalsCorrDistribution> sm_mvd_rep =
2731     std::static_pointer_cast<Pecos::MarginalsCorrDistribution>
2732     (sm_mvd.multivar_dist_rep());
2733 
2734   const SharedVariablesData& sm_svd
2735     = subModel.current_variables().shared_data();
2736 
2737   switch (svm_target) {
2738   case Pecos::NO_TARGET: default:
2739     Cerr << "\nError: secondary mapping target unmatched for string value "
2740 	 << "insertion in NestedModel::string_variable_mapping()" << std::endl;
2741     abort_handler(MODEL_ERROR);
2742   }
2743 }
2744 
2745 
default_interface_active_set()2746 ActiveSet NestedModel::default_interface_active_set()
2747 {
2748   size_t num_fun = numOptInterfPrimary + numOptInterfIneqCon + numOptInterfEqCon;
2749   ActiveSet set;
2750   set.derivative_vector(currentVariables.all_continuous_variable_ids());
2751   bool has_deriv_vars = set.derivative_vector().size() != 0;
2752   ShortArray asv(num_fun, 1);
2753   if(has_deriv_vars) {
2754     if(optInterfGradientType == "analytic") {
2755       for(auto &a : asv)
2756         a |=  2;
2757     } else if(optInterfGradientType == "mixed") {
2758       for(const auto &gi : optInterfGradIdAnalytic)
2759         asv[gi-1] |= 2;
2760     }
2761 
2762     if(optInterfHessianType == "analytic") {
2763       for(auto &a : asv)
2764         a |=  4;
2765     } else if(optInterfHessianType == "mixed") {
2766       for(const auto &hi : optInterfHessIdAnalytic)
2767         asv[hi-1] |= 4;
2768     }
2769   }
2770   set.request_vector(asv);
2771   return set;
2772 }
2773 
2774 } // namespace Dakota
2775