Thyra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Thyra_EpetraModelEvaluator.hpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Thyra: Interfaces and Support for Abstract Numerical Algorithms
5 // Copyright (2004) 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 Roscoe A. Bartlett (bartlettra@ornl.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #ifndef THYRA_EPETRA_MODEL_EVALUATOR_HPP
43 #define THYRA_EPETRA_MODEL_EVALUATOR_HPP
44 
45 #include "Thyra_ModelEvaluatorDefaultBase.hpp"
47 #include "Thyra_LinearOpWithSolveFactoryBase.hpp"
48 #include "EpetraExt_ModelEvaluator.h"
49 #include "Epetra_Map.h"
50 #include "Teuchos_Array.hpp"
51 
52 
53 namespace Thyra {
54 
55 
176  : public ModelEvaluatorDefaultBase<double>,
177  virtual public Teuchos::ParameterListAcceptor
178 {
179 public:
180 
183 
186 
189  const RCP<const EpetraExt::ModelEvaluator> &epetraModel,
190  const RCP<LinearOpWithSolveFactoryBase<double> > &W_factory
191  );
192 
194  void initialize(
195  const RCP<const EpetraExt::ModelEvaluator> &epetraModel,
196  const RCP<LinearOpWithSolveFactoryBase<double> > &W_factory
197  );
198 
201 
207  void setNominalValues( const ModelEvaluatorBase::InArgs<double>& nominalValues );
208 
217  const RCP<const Epetra_Vector> &stateVariableScalingVec
218  );
219 
224 
229 
233  const RCP<const Epetra_Vector> &stateFunctionScalingVec
234  );
235 
240 
242  void uninitialize(
243  RCP<const EpetraExt::ModelEvaluator> *epetraModel = NULL,
244  RCP<LinearOpWithSolveFactoryBase<double> > *W_factory = NULL
245  );
246 
249 
251  bool finalPointWasSolved() const;
252 
254 
257 
259  std::string description() const;
260 
262 
265 
267  void setParameterList(RCP<Teuchos::ParameterList> const& paramList);
276 
278 
281 
283  int Np() const;
285  int Ng() const;
313  void reportFinalPoint(
314  const ModelEvaluatorBase::InArgs<double> &finalPoint
315  ,const bool wasSolved
316  );
317 
319 
320  // Made public to simplify implementation but this is harmless to be public.
321  // Clients should not deal with this type.
322  enum EStateFunctionScaling { STATE_FUNC_SCALING_NONE, STATE_FUNC_SCALING_ROW_SUM };
323 
324 private:
325 
328 
330  RCP<LinearOpBase<double> > create_DfDp_op_impl(int l) const;
332  RCP<LinearOpBase<double> > create_DgDx_dot_op_impl(int j) const;
334  RCP<LinearOpBase<double> > create_DgDx_op_impl(int j) const;
336  RCP<LinearOpBase<double> > create_DgDp_op_impl(int j, int l) const;
338  ModelEvaluatorBase::OutArgs<double> createOutArgsImpl() const;
340  void evalModelImpl(
343  ) const;
344 
346 
347 private:
348 
349  // ////////////////////
350  // Private types
351 
354  typedef std::vector<bool> p_map_is_local_t;
355  typedef std::vector<bool> g_map_is_local_t;
356 
358  p_space_t;
360  g_space_t;
361 
362  // /////////////////////
363  // Private data members
364 
366 
367  RCP<Teuchos::ParameterList> paramList_;
368 
370 
371  RCP<const Epetra_Map> x_map_;
372  p_map_t p_map_;
373  g_map_t g_map_;
374  p_map_is_local_t p_map_is_local_;
375  p_map_is_local_t g_map_is_local_;
376  RCP<const Epetra_Map> f_map_;
377 
379  p_space_t p_space_;
381  g_space_t g_space_;
382 
383  mutable ModelEvaluatorBase::InArgs<double> nominalValues_;
384  mutable ModelEvaluatorBase::InArgs<double> lowerBounds_;
385  mutable ModelEvaluatorBase::InArgs<double> upperBounds_;
386  mutable bool nominalValuesAndBoundsAreUpdated_;
387 
389 
390  EStateFunctionScaling stateFunctionScaling_;
391  mutable RCP<const Epetra_Vector> stateFunctionScalingVec_;
392 
393  RCP<const Epetra_Vector> stateVariableScalingVec_; // S_x
394  mutable RCP<const Epetra_Vector> invStateVariableScalingVec_; // inv(S_x)
395  mutable EpetraExt::ModelEvaluator::InArgs epetraInArgsScaling_;
396  mutable EpetraExt::ModelEvaluator::OutArgs epetraOutArgsScaling_;
397 
398  mutable RCP<Epetra_Vector> x_unscaled_;
399  mutable RCP<Epetra_Vector> x_dot_unscaled_;
400 
401  mutable ModelEvaluatorBase::InArgs<double> prototypeInArgs_;
402  mutable ModelEvaluatorBase::OutArgs<double> prototypeOutArgs_;
403  mutable bool currentInArgsOutArgs_;
404 
405  bool finalPointWasSolved_;
406 
407  // //////////////////////////
408  // Private member functions
409 
411  void convertInArgsFromEpetraToThyra(
412  const EpetraExt::ModelEvaluator::InArgs &epetraInArgs,
414  ) const;
415 
417  void convertInArgsFromThyraToEpetra(
419  EpetraExt::ModelEvaluator::InArgs *epetraInArgs
420  ) const;
421 
423  void convertOutArgsFromThyraToEpetra(
424  // Thyra form of the outArgs
426  // Epetra form of the unscaled output arguments
427  EpetraExt::ModelEvaluator::OutArgs *epetraUnscaledOutArgs,
428  // The passed-in form of W
429  RCP<LinearOpBase<double> > *W_op,
430  RCP<EpetraLinearOp> *efwdW,
431  // The actual Epetra object passed to the underylying EpetraExt::ModelEvaluator
433  ) const;
434 
436  void preEvalScalingSetup(
437  EpetraExt::ModelEvaluator::InArgs *epetraInArgs,
438  EpetraExt::ModelEvaluator::OutArgs *epetraUnscaledOutArgs,
439  const RCP<Teuchos::FancyOStream> &out,
440  const Teuchos::EVerbosityLevel verbLevel
441  ) const;
442 
444  void postEvalScalingSetup(
445  const EpetraExt::ModelEvaluator::OutArgs &epetraUnscaledOutArgs,
446  const RCP<Teuchos::FancyOStream> &out,
447  const Teuchos::EVerbosityLevel verbLevel
448  ) const;
449 
451  void finishConvertingOutArgsFromEpetraToThyra(
452  const EpetraExt::ModelEvaluator::OutArgs &epetraOutArgs,
453  RCP<LinearOpBase<double> > &W_op,
454  RCP<EpetraLinearOp> &efwdW,
456  const ModelEvaluatorBase::OutArgs<double> &outArgs // Output!
457  ) const;
458  // 2007/08/03: rabartl: Above, I pass many of the RCP objects by non-const
459  // reference since I don't want the compiler to perform any implicit
460  // conversions on this RCP objects.
461 
463  void updateNominalValuesAndBounds() const;
464 
466  void updateInArgsOutArgs() const;
467 
469  RCP<EpetraLinearOp> create_epetra_W_op() const;
470 
471 };
472 
473 
474 //
475 // Utility functions
476 //
477 
478 
483 epetraModelEvaluator(
484  const RCP<const EpetraExt::ModelEvaluator> &epetraModel,
485  const RCP<LinearOpWithSolveFactoryBase<double> > &W_factory
486  );
487 
488 
493 convert( const EpetraExt::ModelEvaluator::EDerivativeMultiVectorOrientation &mvOrientation );
494 
495 
499 EpetraExt::ModelEvaluator::EDerivativeMultiVectorOrientation
500 convert( const ModelEvaluatorBase::EDerivativeMultiVectorOrientation &mvOrientation );
501 
502 
507 convert( const EpetraExt::ModelEvaluator::DerivativeProperties &derivativeProperties );
508 
509 
514 convert( const EpetraExt::ModelEvaluator::DerivativeSupport &derivativeSupport );
515 
516 
520 EpetraExt::ModelEvaluator::Derivative
521 convert(
522  const ModelEvaluatorBase::Derivative<double> &derivative,
523  const RCP<const Epetra_Map> &fnc_map,
524  const RCP<const Epetra_Map> &var_map
525  );
526 
527 EpetraExt::ModelEvaluator::MPDerivative
528 convert(
529  const ModelEvaluatorBase::MPDerivative &derivative,
530  const RCP<const Epetra_Map> &fnc_map,
531  const RCP<const Epetra_Map> &var_map
532  );
533 
534 } // namespace Thyra
535 
536 
537 #endif // THYRA_EPETRA_MODEL_EVALUATOR_HPP
const ModelEvaluatorBase::InArgs< double > & getFinalPoint() const
ModelEvaluatorBase::InArgs< double > getLowerBounds() const
void setNominalValues(const ModelEvaluatorBase::InArgs< double > &nominalValues)
Set the nominal values.
RCP< const EpetraExt::ModelEvaluator > getEpetraModel() const
Default base class for concrete model evaluators.
void reportFinalPoint(const ModelEvaluatorBase::InArgs< double > &finalPoint, const bool wasSolved)
RCP< const Teuchos::ParameterList > getValidParameters() const
RCP< LinearOpBase< double > > create_W_op() const
ModelEvaluatorBase::InArgs< double > getNominalValues() const
RCP< Teuchos::ParameterList > getNonconstParameterList()
RCP< PreconditionerBase< double > > create_W_prec() const
Returns null currently.
Simple aggregate class that stores a derivative object as a general linear operator or as a multi-vec...
void setStateVariableScalingVec(const RCP< const Epetra_Vector > &stateVariableScalingVec)
Set the state variable scaling vector s_x (see above).
ModelEvaluatorBase::InArgs< double > createInArgs() const
RCP< const Epetra_Vector > getStateVariableScalingVec() const
Get the inverse state variable scaling vector inv_s_x (see above).
RCP< const Teuchos::ParameterList > getParameterList() const
RCP< const LinearOpWithSolveFactoryBase< double > > get_W_factory() const
RCP< const Epetra_Vector > getStateVariableInvScalingVec() const
Get the state variable scaling vector s_x (see above).
ModelEvaluatorBase::InArgs< double > getUpperBounds() const
void setStateFunctionScalingVec(const RCP< const Epetra_Vector > &stateFunctionScalingVec)
Set the state function scaling vector s_f (see above).
RCP< const Epetra_Vector > getStateFunctionScalingVec() const
Get the state function scaling vector s_f (see above).
void initialize(const RCP< const EpetraExt::ModelEvaluator > &epetraModel, const RCP< LinearOpWithSolveFactoryBase< double > > &W_factory)
RCP< const VectorSpaceBase< double > > get_f_space() const
RCP< const VectorSpaceBase< double > > get_x_space() const
RCP< const VectorSpaceBase< double > > get_g_space(int j) const
Simple public strict containing properties of a derivative object.
RCP< const VectorSpaceBase< double > > get_p_space(int l) const
Concrete Adapter subclass that takes an EpetraExt::ModelEvaluator object and wraps it as a Thyra::Mod...
RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
void setParameterList(RCP< Teuchos::ParameterList > const &paramList)
void uninitialize(RCP< const EpetraExt::ModelEvaluator > *epetraModel=NULL, RCP< LinearOpWithSolveFactoryBase< double > > *W_factory=NULL)
Teuchos::ArrayView< const std::string > get_g_names(int j) const
Determines the forms of a general derivative that are supported.
Simple aggregate class that stores a derivative object as a general linear operator or as a multi-vec...
RCP< Teuchos::ParameterList > unsetParameterList()