Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_MPModelEvaluatorAdapter.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Stokhos Package
5 // Copyright (2009) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
38 //
39 // ***********************************************************************
40 // @HEADER
41 
46 #include "Teuchos_Assert.hpp"
47 #include "Epetra_Map.h"
48 
52  me(me_)
53 {
54 }
55 
56 // Overridden from EpetraExt::ModelEvaluator
57 
60 get_x_map() const
61 {
62  return me->get_x_map();
63 }
64 
67 get_f_map() const
68 {
69  return me->get_f_map();
70 }
71 
74 get_p_map(int l) const
75 {
76  return me->get_p_map(l);
77 }
78 
81 get_g_map(int l) const
82 {
83  return me->get_g_map(l);
84 }
85 
88 get_p_names(int l) const
89 {
90  return me->get_p_names(l);
91 }
92 
95 get_x_init() const
96 {
97  return me->get_x_init();
98 }
99 
102 get_p_init(int l) const
103 {
104  return me->get_p_init(l);
105 }
106 
109 create_W() const
110 {
111  return me->create_W();
112 }
113 
114 EpetraExt::ModelEvaluator::InArgs
117 {
118  InArgsSetup inArgs;
119  InArgs me_inargs = me->createInArgs();
120 
121  inArgs.setModelEvalDescription(this->description());
122  inArgs.set_Np(me_inargs.Np());
123  inArgs.setSupports(IN_ARG_x_dot, me_inargs.supports(IN_ARG_x_dot));
124  inArgs.setSupports(IN_ARG_x, me_inargs.supports(IN_ARG_x));
125  inArgs.setSupports(IN_ARG_t, me_inargs.supports(IN_ARG_t));
126  inArgs.setSupports(IN_ARG_alpha, me_inargs.supports(IN_ARG_alpha));
127  inArgs.setSupports(IN_ARG_beta, me_inargs.supports(IN_ARG_beta));
128 
129  for (int i=0; i<me_inargs.Np(); i++)
130  inArgs.setSupports(IN_ARG_p_mp, i, true);
131  inArgs.setSupports(IN_ARG_x_mp, me_inargs.supports(IN_ARG_x));
132  inArgs.setSupports(IN_ARG_x_dot_mp, me_inargs.supports(IN_ARG_x_dot));
133 
134  return inArgs;
135 }
136 
137 EpetraExt::ModelEvaluator::OutArgs
140 {
141  OutArgsSetup outArgs;
142  OutArgs me_outargs = me->createOutArgs();
143 
144  outArgs.setModelEvalDescription(this->description());
145  outArgs.set_Np_Ng(me_outargs.Np(), me_outargs.Ng());
146  outArgs.setSupports(OUT_ARG_f, me_outargs.supports(OUT_ARG_f));
147  outArgs.setSupports(OUT_ARG_W, me_outargs.supports(OUT_ARG_W));
148  for (int j=0; j<me_outargs.Np(); j++)
149  outArgs.setSupports(OUT_ARG_DfDp, j,
150  me_outargs.supports(OUT_ARG_DfDp, j));
151  for (int i=0; i<me_outargs.Ng(); i++) {
152  outArgs.setSupports(OUT_ARG_DgDx, i,
153  me_outargs.supports(OUT_ARG_DgDx, i));
154  outArgs.setSupports(OUT_ARG_DgDx_dot, i,
155  me_outargs.supports(OUT_ARG_DgDx_dot, i));
156  for (int j=0; j<me_outargs.Np(); j++)
157  outArgs.setSupports(OUT_ARG_DgDp, i, j,
158  me_outargs.supports(OUT_ARG_DgDp, i, j));
159  }
160 
161  outArgs.setSupports(OUT_ARG_f_mp, me_outargs.supports(OUT_ARG_f));
162  if (me_outargs.supports(OUT_ARG_W)) {
163  outArgs.set_W_properties(me_outargs.get_W_properties());
164  outArgs.setSupports(OUT_ARG_W_mp, true);
165  }
166  for (int j=0; j<me_outargs.Np(); j++)
167  outArgs.setSupports(OUT_ARG_DfDp_mp, j,
168  me_outargs.supports(OUT_ARG_DfDp, j));
169  for (int i=0; i<me_outargs.Ng(); i++) {
170  outArgs.setSupports(OUT_ARG_g_mp, i, true);
171  outArgs.setSupports(OUT_ARG_DgDx_mp, i,
172  me_outargs.supports(OUT_ARG_DgDx, i));
173  outArgs.setSupports(OUT_ARG_DgDx_dot_mp, i,
174  me_outargs.supports(OUT_ARG_DgDx_dot, i));
175  for (int j=0; j<me_outargs.Np(); j++)
176  outArgs.setSupports(OUT_ARG_DgDp_mp, i, j,
177  me_outargs.supports(OUT_ARG_DgDp, i, j));
178  }
179 
180  return outArgs;
181 }
182 
183 void
185 evalModel(const InArgs& inArgs, const OutArgs& outArgs) const
186 {
187  // Create underlying inargs
188  InArgs me_inargs = me->createInArgs();
189  if (me_inargs.supports(IN_ARG_x))
190  me_inargs.set_x(inArgs.get_x());
191  if (me_inargs.supports(IN_ARG_x_dot))
192  me_inargs.set_x_dot(inArgs.get_x_dot());
193  if (me_inargs.supports(IN_ARG_alpha))
194  me_inargs.set_alpha(inArgs.get_alpha());
195  if (me_inargs.supports(IN_ARG_beta))
196  me_inargs.set_beta(inArgs.get_beta());
197  if (me_inargs.supports(IN_ARG_t))
198  me_inargs.set_t(inArgs.get_t());
199  for (int i=0; i<inArgs.Np(); i++)
200  me_inargs.set_p(i, inArgs.get_p(i));
201 
202  // Create underlying outargs
203  OutArgs me_outargs = me->createOutArgs();
204  if (me_outargs.supports(OUT_ARG_f))
205  me_outargs.set_f(outArgs.get_f());
206  if (me_outargs.supports(OUT_ARG_W))
207  me_outargs.set_W(outArgs.get_W());
208  for (int j=0; j<outArgs.Np(); j++)
209  if (!outArgs.supports(OUT_ARG_DfDp, j).none())
210  me_outargs.set_DfDp(j, outArgs.get_DfDp(j));
211  for (int i=0; i<outArgs.Ng(); i++) {
212  me_outargs.set_g(i, outArgs.get_g(i));
213  if (!outArgs.supports(OUT_ARG_DgDx, i).none())
214  me_outargs.set_DgDx(i, outArgs.get_DgDx(i));
215  if (!outArgs.supports(OUT_ARG_DgDx_dot, i).none())
216  me_outargs.set_DgDx(i, outArgs.get_DgDx_dot(i));
217  for (int j=0; j<outArgs.Np(); j++)
218  if (!outArgs.supports(OUT_ARG_DgDp, i, j).none())
219  me_outargs.set_DgDp(i, j, outArgs.get_DgDp(i,j));
220  }
221 
222  mp_const_vector_t x_mp;
223  mp_const_vector_t x_dot_mp;
224  Teuchos::Array<mp_const_vector_t> p_mp(inArgs.Np());
225  mp_vector_t f_mp;
226  mp_operator_t W_mp;
227  Teuchos::Array<MPDerivative> dfdp_mp(outArgs.Np());
228  Teuchos::Array<mp_vector_t> g_mp(outArgs.Ng());
229  Teuchos::Array<MPDerivative> dgdx_mp(outArgs.Ng());
230  Teuchos::Array<MPDerivative> dgdx_dot_mp(outArgs.Ng());
231  Teuchos::Array< Teuchos::Array<MPDerivative> > dgdp_mp(outArgs.Ng());
232  int num_mp = 0;
233 
234  if (inArgs.supports(IN_ARG_x_mp)) {
235  x_mp = inArgs.get_x_mp();
236  if (x_mp != Teuchos::null) {
237  num_mp = x_mp->size();
238  }
239  }
240  if (inArgs.supports(IN_ARG_x_dot_mp)) {
241  x_dot_mp = inArgs.get_x_dot_mp();
242  if (x_dot_mp != Teuchos::null) {
243  num_mp = x_dot_mp->size();
244  }
245  }
246  for (int i=0; i<inArgs.Np(); i++) {
247  p_mp[i] = inArgs.get_p_mp(i);
248  if (p_mp[i] != Teuchos::null) {
249  num_mp = p_mp[i]->size();
250  }
251  }
252  if (outArgs.supports(OUT_ARG_f_mp)) {
253  f_mp = outArgs.get_f_mp();
254  if (f_mp != Teuchos::null)
255  f_mp->init(0.0);
256  }
257  if (outArgs.supports(OUT_ARG_W_mp)) {
258  W_mp = outArgs.get_W_mp();
259  if (W_mp != Teuchos::null)
260  W_mp->init(0.0);
261  }
262  for (int i=0; i<inArgs.Np(); i++) {
263  if (!outArgs.supports(OUT_ARG_DfDp_mp, i).none()) {
264  dfdp_mp[i] = outArgs.get_DfDp_mp(i);
265  if (dfdp_mp[i].getMultiVector() != Teuchos::null)
266  dfdp_mp[i].getMultiVector()->init(0.0);
267  else if (dfdp_mp[i].getLinearOp() != Teuchos::null)
268  dfdp_mp[i].getLinearOp()->init(0.0);
269  }
270  }
271 
272  for (int i=0; i<outArgs.Ng(); i++) {
273  g_mp[i] = outArgs.get_g_mp(i);
274  if (g_mp[i] != Teuchos::null)
275  g_mp[i]->init(0.0);
276 
277  if (!outArgs.supports(OUT_ARG_DgDx_mp, i).none()) {
278  dgdx_mp[i] = outArgs.get_DgDx_mp(i);
279  if (dgdx_mp[i].getMultiVector() != Teuchos::null)
280  dgdx_mp[i].getMultiVector()->init(0.0);
281  else if (dgdx_mp[i].getLinearOp() != Teuchos::null)
282  dgdx_mp[i].getLinearOp()->init(0.0);
283  }
284 
285  if (!outArgs.supports(OUT_ARG_DgDx_dot_mp, i).none()) {
286  dgdx_dot_mp[i] = outArgs.get_DgDx_dot_mp(i);
287  if (dgdx_dot_mp[i].getMultiVector() != Teuchos::null)
288  dgdx_dot_mp[i].getMultiVector()->init(0.0);
289  else if (dgdx_dot_mp[i].getLinearOp() != Teuchos::null)
290  dgdx_dot_mp[i].getLinearOp()->init(0.0);
291  }
292 
293  dgdp_mp[i].resize(outArgs.Np());
294  for (int j=0; j<outArgs.Np(); j++) {
295  if (!outArgs.supports(OUT_ARG_DgDp_mp, i, j).none()) {
296  dgdp_mp[i][j] = outArgs.get_DgDp_mp(i,j);
297  if (dgdp_mp[i][j].getMultiVector() != Teuchos::null)
298  dgdp_mp[i][j].getMultiVector()->init(0.0);
299  else if (dgdp_mp[i][j].getLinearOp() != Teuchos::null)
300  dgdp_mp[i][j].getLinearOp()->init(0.0);
301  }
302  }
303  }
304 
305  for (int qp=0; qp<num_mp; qp++) {
306 
307  // Set InArgs
308  if (x_mp != Teuchos::null)
309  me_inargs.set_x(x_mp->getCoeffPtr(qp));
310  if (x_dot_mp != Teuchos::null)
311  me_inargs.set_x_dot(x_dot_mp->getCoeffPtr(qp));
312  for (int i=0; i<inArgs.Np(); i++)
313  if (p_mp[i] != Teuchos::null)
314  me_inargs.set_p(i, p_mp[i]->getCoeffPtr(qp));
315 
316  // Set OutArgs
317  if (f_mp != Teuchos::null)
318  me_outargs.set_f(f_mp->getCoeffPtr(qp));
319  if (W_mp != Teuchos::null)
320  me_outargs.set_W(W_mp->getCoeffPtr(qp));
321 
322  for (int i=0; i<inArgs.Np(); i++) {
323  if (!dfdp_mp[i].isEmpty()) {
324  if (dfdp_mp[i].getMultiVector() != Teuchos::null) {
325  Derivative deriv(dfdp_mp[i].getMultiVector()->getCoeffPtr(qp),
326  dfdp_mp[i].getMultiVectorOrientation());
327  me_outargs.set_DfDp(i, deriv);
328  }
329  else if (dfdp_mp[i].getLinearOp() != Teuchos::null) {
330  Derivative deriv(dfdp_mp[i].getLinearOp()->getCoeffPtr(qp));
331  me_outargs.set_DfDp(i, deriv);
332  }
333  }
334  }
335 
336  for (int i=0; i<outArgs.Ng(); i++) {
337  if (g_mp[i] != Teuchos::null)
338  me_outargs.set_g(i, g_mp[i]->getCoeffPtr(qp));
339  if (!dgdx_dot_mp[i].isEmpty()) {
340  if (dgdx_dot_mp[i].getMultiVector() != Teuchos::null) {
341  Derivative deriv(dgdx_dot_mp[i].getMultiVector()->getCoeffPtr(qp),
342  dgdx_dot_mp[i].getMultiVectorOrientation());
343  me_outargs.set_DgDx_dot(i, deriv);
344  }
345  else if (dgdx_dot_mp[i].getLinearOp() != Teuchos::null) {
346  Derivative deriv(dgdx_dot_mp[i].getLinearOp()->getCoeffPtr(qp));
347  me_outargs.set_DgDx_dot(i, deriv);
348  }
349  }
350  if (!dgdx_mp[i].isEmpty()) {
351  if (dgdx_mp[i].getMultiVector() != Teuchos::null) {
352  Derivative deriv(dgdx_mp[i].getMultiVector()->getCoeffPtr(qp),
353  dgdx_mp[i].getMultiVectorOrientation());
354  me_outargs.set_DgDx(i, deriv);
355  }
356  else if (dgdx_mp[i].getLinearOp() != Teuchos::null) {
357  Derivative deriv(dgdx_mp[i].getLinearOp()->getCoeffPtr(qp));
358  me_outargs.set_DgDx(i, deriv);
359  }
360  }
361  for (int j=0; j<outArgs.Np(); j++) {
362  if (!dgdp_mp[i][j].isEmpty()) {
363  if (dgdp_mp[i][j].getMultiVector() != Teuchos::null) {
364  Derivative deriv(dgdp_mp[i][j].getMultiVector()->getCoeffPtr(qp),
365  dgdp_mp[i][j].getMultiVectorOrientation());
366  me_outargs.set_DgDp(i, j, deriv);
367  }
368  else if (dgdp_mp[i][j].getLinearOp() != Teuchos::null) {
369  Derivative deriv(dgdp_mp[i][j].getLinearOp()->getCoeffPtr(qp));
370  me_outargs.set_DgDp(i, j, deriv);
371  }
372  }
373  }
374  }
375 
376  // Evaluate model
377  me->evalModel(me_inargs, me_outargs);
378 
379  }
380 
381  // Evaluate single-point model
382  if (num_mp == 0)
383  me->evalModel(me_inargs, me_outargs);
384 }
Teuchos::RCP< const Epetra_Map > get_g_map(int l) const
Return observation vector map.
MPModelEvaluatorAdapter(const Teuchos::RCP< EpetraExt::ModelEvaluator > &me)
Teuchos::RCP< const Epetra_Map > get_f_map() const
Return residual vector map.
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Evaluate model on InArgs.
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Return initial solution.
Teuchos::RCP< const Epetra_Map > get_x_map() const
Return solution vector map.
Teuchos::RCP< const Epetra_Map > get_p_map(int l) const
Return parameter vector map.
void resize(size_type new_size, const value_type &x=value_type())
Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Return initial parameters.
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Return array of parameter names.
size_type size() const
Teuchos::RCP< Epetra_Operator > create_W() const
Create W = alpha*M + beta*J matrix.
OutArgs createOutArgs() const
Create OutArgs.