EpetraExt Development
EpetraModelEval2DSim.cpp
Go to the documentation of this file.
00001 /*
00002 //@HEADER
00003 // ***********************************************************************
00004 //
00005 //     EpetraExt: Epetra Extended - Linear Algebra Services Package
00006 //                 Copyright (2011) Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00039 //
00040 // ***********************************************************************
00041 //@HEADER
00042 */
00043 
00044 #include "EpetraModelEval2DSim.hpp"
00045 #include "Teuchos_ScalarTraits.hpp"
00046 #include "Epetra_SerialComm.h"
00047 #include "Epetra_CrsMatrix.h"
00048 
00049 EpetraModelEval2DSim::EpetraModelEval2DSim(
00050   const double         d
00051   ,const double        p0
00052   ,const double        p1
00053   ,const double        x00
00054   ,const double        x01
00055   ,const bool          showGetInvalidArg
00056   )
00057   :d_(d),showGetInvalidArg_(showGetInvalidArg)
00058 {
00059   using Teuchos::rcp;
00060 
00061   epetra_comm_ = rcp(new Epetra_SerialComm());
00062 
00063   const int nx = 2;
00064 
00065   map_x_ = rcp(new Epetra_Map(nx,0,*epetra_comm_));
00066   x0_ = rcp(new Epetra_Vector(*map_x_));  (*x0_)[0] = x00; (*x0_)[1] = x01;
00067   p_  = rcp(new Epetra_Vector(*map_x_));  (*p_)[0] = p0; (*p_)[1] = p1;
00068 
00069   // Initialize the graph for W CrsMatrix object
00070   W_graph_ = rcp(new Epetra_CrsGraph(::Copy,*map_x_,nx));
00071   {
00072     int indices[nx] = { 0, 1 };
00073     for( int i = 0; i < nx; ++i )
00074       W_graph_->InsertGlobalIndices(i,nx,indices);
00075   }
00076   W_graph_->FillComplete();
00077 
00078   isInitialized_ = true;
00079 
00080 }
00081 
00082 // Overridden from EpetraExt::ModelEvaluator
00083 
00084 Teuchos::RefCountPtr<const Epetra_Map>
00085 EpetraModelEval2DSim::get_x_map() const
00086 {
00087   return map_x_;
00088 }
00089 
00090 Teuchos::RefCountPtr<const Epetra_Map>
00091 EpetraModelEval2DSim::get_f_map() const
00092 {
00093   return map_x_;
00094 }
00095 
00096 Teuchos::RefCountPtr<const Epetra_Vector>
00097 EpetraModelEval2DSim::get_x_init() const
00098 {
00099   return x0_;
00100 }
00101 
00102 Teuchos::RefCountPtr<Epetra_Operator>
00103 EpetraModelEval2DSim::create_W() const
00104 {
00105   return Teuchos::rcp(new Epetra_CrsMatrix(::Copy,*W_graph_));
00106 }
00107 
00108 EpetraExt::ModelEvaluator::InArgs
00109 EpetraModelEval2DSim::createInArgs() const
00110 {
00111   InArgsSetup inArgs;
00112   inArgs.setModelEvalDescription(this->description());
00113   inArgs.setSupports(IN_ARG_x,true);
00114   return inArgs;
00115 }
00116 
00117 EpetraExt::ModelEvaluator::OutArgs
00118 EpetraModelEval2DSim::createOutArgs() const
00119 {
00120   OutArgsSetup outArgs;
00121   outArgs.setModelEvalDescription(this->description());
00122   outArgs.setSupports(OUT_ARG_f,true);
00123   outArgs.setSupports(OUT_ARG_W,true);
00124   outArgs.set_W_properties(
00125     DerivativeProperties(
00126       DERIV_LINEARITY_NONCONST
00127       ,DERIV_RANK_FULL
00128       ,true // supportsAdjoint
00129       )
00130     );
00131   return outArgs;
00132 }
00133 
00134 void EpetraModelEval2DSim::evalModel( const InArgs& inArgs, const OutArgs& outArgs ) const
00135 {
00136   using Teuchos::dyn_cast;
00137   using Teuchos::rcp_dynamic_cast;
00138   //
00139   // Get the input arguments
00140   //
00141   const Epetra_Vector &x = *inArgs.get_x();
00142   //
00143   // Get the output arguments
00144   //
00145   Epetra_Vector       *f_out = outArgs.get_f().get();
00146   Epetra_Operator     *W_out = outArgs.get_W().get();
00147   if(showGetInvalidArg_) {
00148     Epetra_Vector *g_out = outArgs.get_g(0).get();
00149   }
00150   //
00151   // Compute the functions
00152   //
00153   const Epetra_Vector &p = *p_;
00154   if(f_out) {
00155     Epetra_Vector &f = *f_out;
00156     f[0] =        x[0]      + x[1]*x[1] - p[0];
00157     f[1] = d_ * ( x[0]*x[0] - x[1]      - p[1] );
00158   }
00159   if(W_out) {
00160     Epetra_CrsMatrix &DfDx = dyn_cast<Epetra_CrsMatrix>(*W_out);
00161     DfDx.PutScalar(0.0);
00162     //
00163     // Fill W = DfDx
00164     //
00165     // W = DfDx = [      1.0,  2*x[1] ]
00166     //            [ 2*d*x[0],     -d  ]
00167     //
00168     double values[2];
00169     int indexes[2];
00170     // Row [0]
00171     values[0] = 1.0;           indexes[0] = 0;
00172     values[1] = 2.0*x[1];      indexes[1] = 1;
00173     DfDx.SumIntoGlobalValues( 0, 2, values, indexes );
00174     // Row [1]
00175     values[0] = 2.0*d_*x[0];   indexes[0] = 0;
00176     values[1] = -d_;           indexes[1] = 1;
00177     DfDx.SumIntoGlobalValues( 1, 2, values, indexes );
00178   }
00179 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines