EpetraExt Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
EpetraModelEval2DSim.cpp
Go to the documentation of this file.
1/*
2//@HEADER
3// ***********************************************************************
4//
5// EpetraExt: Epetra Extended - Linear Algebra Services Package
6// Copyright (2011) Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39//
40// ***********************************************************************
41//@HEADER
42*/
43
45#include "Teuchos_ScalarTraits.hpp"
46#include "Epetra_SerialComm.h"
47#include "Epetra_CrsMatrix.h"
48
50 const double d
51 ,const double p0
52 ,const double p1
53 ,const double x00
54 ,const double x01
55 ,const bool showGetInvalidArg
56 )
57 :d_(d),showGetInvalidArg_(showGetInvalidArg)
58{
59 using Teuchos::rcp;
60
62
63 const int nx = 2;
64
65 map_x_ = rcp(new Epetra_Map(nx,0,*epetra_comm_));
66 x0_ = rcp(new Epetra_Vector(*map_x_)); (*x0_)[0] = x00; (*x0_)[1] = x01;
67 p_ = rcp(new Epetra_Vector(*map_x_)); (*p_)[0] = p0; (*p_)[1] = p1;
68
69 // Initialize the graph for W CrsMatrix object
70 W_graph_ = rcp(new Epetra_CrsGraph(::Copy,*map_x_,nx));
71 {
72 int indices[nx] = { 0, 1 };
73 for( int i = 0; i < nx; ++i )
74 W_graph_->InsertGlobalIndices(i,nx,indices);
75 }
76 W_graph_->FillComplete();
77
78 isInitialized_ = true;
79
80}
81
82// Overridden from EpetraExt::ModelEvaluator
83
84Teuchos::RCP<const Epetra_Map>
86{
87 return map_x_;
88}
89
90Teuchos::RCP<const Epetra_Map>
92{
93 return map_x_;
94}
95
96Teuchos::RCP<const Epetra_Vector>
98{
99 return x0_;
100}
101
102Teuchos::RCP<Epetra_Operator>
104{
105 return Teuchos::rcp(new Epetra_CrsMatrix(::Copy,*W_graph_));
106}
107
110{
111 InArgsSetup inArgs;
112 inArgs.setModelEvalDescription(this->description());
113 inArgs.setSupports(IN_ARG_x,true);
114 return inArgs;
115}
116
119{
120 OutArgsSetup outArgs;
121 outArgs.setModelEvalDescription(this->description());
122 outArgs.setSupports(OUT_ARG_f,true);
123 outArgs.setSupports(OUT_ARG_W,true);
124 outArgs.set_W_properties(
128 ,true // supportsAdjoint
129 )
130 );
131 return outArgs;
132}
133
134void EpetraModelEval2DSim::evalModel( const InArgs& inArgs, const OutArgs& outArgs ) const
135{
136 using Teuchos::dyn_cast;
137 using Teuchos::rcp_dynamic_cast;
138 //
139 // Get the input arguments
140 //
141 const Epetra_Vector &x = *inArgs.get_x();
142 //
143 // Get the output arguments
144 //
145 Epetra_Vector *f_out = outArgs.get_f().get();
146 Epetra_Operator *W_out = outArgs.get_W().get();
148 Epetra_Vector *g_out = outArgs.get_g(0).get();
149 }
150 //
151 // Compute the functions
152 //
153 const Epetra_Vector &p = *p_;
154 if(f_out) {
155 Epetra_Vector &f = *f_out;
156 f[0] = x[0] + x[1]*x[1] - p[0];
157 f[1] = d_ * ( x[0]*x[0] - x[1] - p[1] );
158 }
159 if(W_out) {
160 Epetra_CrsMatrix &DfDx = dyn_cast<Epetra_CrsMatrix>(*W_out);
161 DfDx.PutScalar(0.0);
162 //
163 // Fill W = DfDx
164 //
165 // W = DfDx = [ 1.0, 2*x[1] ]
166 // [ 2*d*x[0], -d ]
167 //
168 double values[2];
169 int indexes[2];
170 // Row [0]
171 values[0] = 1.0; indexes[0] = 0;
172 values[1] = 2.0*x[1]; indexes[1] = 1;
173 DfDx.SumIntoGlobalValues( 0, 2, values, indexes );
174 // Row [1]
175 values[0] = 2.0*d_*x[0]; indexes[0] = 0;
176 values[1] = -d_; indexes[1] = 1;
177 DfDx.SumIntoGlobalValues( 1, 2, values, indexes );
178 }
179}
void setSupports(EInArgsMembers arg, bool supports=true)
void setModelEvalDescription(const std::string &modelEvalDescription)
Teuchos::RCP< const Epetra_Vector > get_x() const
Set solution vector Taylor polynomial.
void setModelEvalDescription(const std::string &modelEvalDescription)
void set_W_properties(const DerivativeProperties &properties)
void setSupports(EOutArgsMembers arg, bool supports=true)
Teuchos::RCP< Epetra_Operator > get_W() const
Evaluation< Epetra_Vector > get_g(int j) const
Get g(j) where 0 <= j && j < this->Ng().
Evaluation< Epetra_Vector > get_f() const
Teuchos::RCP< Epetra_Vector > x0_
EpetraModelEval2DSim(const double d=10.0, const double p0=2.0, const double p1=0.0, const double x00=1.0, const double x01=1.0, const bool showGetInvalidArg=false)
Teuchos::RCP< Epetra_Vector > p_
Teuchos::RCP< const Epetra_Comm > epetra_comm_
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Teuchos::RCP< const Epetra_Map > get_x_map() const
Teuchos::RCP< Epetra_Operator > create_W() const
Teuchos::RCP< const Epetra_Map > get_f_map() const
Teuchos::RCP< const Epetra_Map > map_x_
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Teuchos::RCP< Epetra_CrsGraph > W_graph_
int PutScalar(double ScalarConstant)
virtual int SumIntoGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)