Thyra Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Simple2DTpetraModelEvaluator_def.hpp
Go to the documentation of this file.
1/*
2// @HEADER
3// ***********************************************************************
4//
5// Thyra: Interfaces and Support for Abstract Numerical Algorithms
6// Copyright (2004) Sandia Corporation
7//
8// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9// license for use of this work by or on behalf of the U.S. Government.
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 Roscoe A. Bartlett (bartlettra@ornl.gov)
39//
40// ***********************************************************************
41// @HEADER
42*/
43
44
45#ifndef SIMPLE_2D_TPETRA_MODEL_EVALUATOR_HPP
46#define SIMPLE_2D_TPETRA_MODEL_EVALUATOR_HPP
47
48
51#include "Tpetra_CrsMatrix.hpp"
53
54
55// Constructors/Initializers/Accessors
56
57
58template<class Scalar>
60 : d_(0.0)
61{
62
63 using Teuchos::RCP;
64 using Teuchos::rcp;
65 using Teuchos::tuple;
66 using Thyra::VectorBase;
67 using Thyra::createMember;
68 typedef Thyra::ModelEvaluatorBase MEB;
70
71 const int dim = 2;
72
73 //
74 // A) Create the structure for the problem
75 //
76
77 MEB::InArgsSetup<Scalar> inArgs;
78 inArgs.setModelEvalDescription(this->description());
79 inArgs.setSupports(MEB::IN_ARG_x);
80 prototypeInArgs_ = inArgs;
81
82 MEB::OutArgsSetup<Scalar> outArgs;
83 outArgs.setModelEvalDescription(this->description());
84 outArgs.setSupports(MEB::OUT_ARG_f);
85 outArgs.setSupports(MEB::OUT_ARG_W_op);
86 prototypeOutArgs_ = outArgs;
87
88 //
89 // B) Create the Tpetra objects
90 //
91
92 const RCP<const Teuchos::Comm<int> > comm =
94
95 const RCP<const Tpetra::Map<> > map = rcp(new Tpetra::Map<> (dim, 0, comm));
96
97 typedef Tpetra::Map<>::global_ordinal_type GO;
98
99 W_op_graph_ = rcp(new Tpetra::CrsGraph<> (map, dim));
100 W_op_graph_->insertGlobalIndices(0, tuple<GO>(0, 1)());
101 W_op_graph_->insertGlobalIndices(1, tuple<GO>(0, 1)());
102 W_op_graph_->fillComplete();
103
104 p_.resize(dim, ST::zero());
105
106 x0_ = rcp(new Tpetra::Vector<Scalar>(map));
107 x0_->putScalar(ST::zero());
108
109 //
110 // C) Create the Thyra wrapped objects
111 //
112
115
116 nominalValues_ = inArgs;
118
119 //
120 // D) Set initial values through interface functions
121 //
122
123 set_d(10.0);
124 set_p(Teuchos::tuple<Scalar>(2.0, 0.0)());
125 set_x0(Teuchos::tuple<Scalar>(1.0, 1.0)());
126
127}
128
129
130template<class Scalar>
132{
133 d_ = d;
134}
135
136
137template<class Scalar>
139{
140#ifdef TEUCHOS_DEBUG
141 TEUCHOS_ASSERT_EQUALITY(p_.size(), p.size());
142#endif
143 p_().assign(p);
144}
145
146
147template<class Scalar>
149{
150#ifdef TEUCHOS_DEBUG
151 TEUCHOS_ASSERT_EQUALITY(x_space_->dim(), x0_in.size());
152#endif
153 x0_->get1dViewNonConst()().assign(x0_in);
154}
155
156
157// Public functions overridden from ModelEvaulator
158
159
160template<class Scalar>
163{
164 return x_space_;
165}
166
167
168template<class Scalar>
171{
172 return f_space_;
173}
174
175
176template<class Scalar>
177Thyra::ModelEvaluatorBase::InArgs<Scalar>
179{
180 return nominalValues_;
181}
182
183
184template<class Scalar>
187{
189 Teuchos::RCP<Tpetra::Operator<Scalar,int> >(
190 Teuchos::rcp(new Tpetra::CrsMatrix<Scalar,int>(W_op_graph_))
191 )
192 );
193}
194
195
196template<class Scalar>
197Thyra::ModelEvaluatorBase::InArgs<Scalar>
199{
200 return prototypeInArgs_;
201}
202
203
204// Private functions overridden from ModelEvaulatorDefaultBase
205
206
207template<class Scalar>
208Thyra::ModelEvaluatorBase::OutArgs<Scalar>
210{
211 return prototypeOutArgs_;
212}
213
214
215template<class Scalar>
217 const Thyra::ModelEvaluatorBase::InArgs<Scalar> &inArgs,
218 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs
219 ) const
220{
221 using Teuchos::RCP;
222 using Teuchos::ArrayRCP;
223 using Teuchos::Array;
224 using Teuchos::tuple;
225 using Teuchos::rcp_dynamic_cast;
227 typedef Tpetra::Map<>::global_ordinal_type GO;
229
230 const RCP<const Tpetra::Vector<Scalar, int> > x_vec =
231 ConverterT::getConstTpetraVector(inArgs.get_x());
232 const ArrayRCP<const Scalar> x = x_vec->get1dView();
233
234 const RCP<Tpetra::Vector<Scalar, int> > f_vec =
235 ConverterT::getTpetraVector(outArgs.get_f());
236
237 const RCP<Tpetra::CrsMatrix<Scalar, int> > W =
238 rcp_dynamic_cast<Tpetra::CrsMatrix<Scalar,int> >(
239 ConverterT::getTpetraOperator(outArgs.get_W_op()),
240 true
241 );
242
243 if (nonnull(f_vec)) {
244 const ArrayRCP<Scalar> f = f_vec->get1dViewNonConst();
245 f[0] = x[0] + x[1]*x[1] - p_[0];
246 f[1] = d_ * (x[0]*x[0] -x[1] - p_[1]);
247 }
248
249 if (nonnull(W)) {
250 W->setAllToScalar(ST::zero());
251 W->sumIntoGlobalValues(0, tuple<GO>(0, 1), tuple<Scalar>(1.0, 2.0*x[1]));
252 W->sumIntoGlobalValues(1, tuple<GO>(0, 1), tuple<Scalar>(2.0*d_*x[0], -d_));
253 }
254
255}
256
257
258#endif // SIMPLE_2D_TPETRA_MODEL_EVALUATOR_HPP
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > x_space_
Thyra::ModelEvaluatorBase::InArgs< Scalar > getNominalValues() const
Thyra::ModelEvaluatorBase::OutArgs< Scalar > prototypeOutArgs_
void set_p(const Teuchos::ArrayView< const Scalar > &p)
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_f_space() const
Thyra::ModelEvaluatorBase::InArgs< Scalar > prototypeInArgs_
Teuchos::RCP< Tpetra::CrsGraph<> > W_op_graph_
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const
void evalModelImpl(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_W_op() const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > f_space_
Thyra::ModelEvaluatorBase::OutArgs< Scalar > createOutArgsImpl() const
Teuchos::RCP< Tpetra::Vector< Scalar > > x0_
Thyra::ModelEvaluatorBase::InArgs< Scalar > nominalValues_
Thyra::ModelEvaluatorBase::InArgs< Scalar > createInArgs() const
void set_x0(const Teuchos::ArrayView< const Scalar > &x0)
void resize(size_type new_size, const value_type &x=value_type())
static Teuchos::RCP< const Comm< OrdinalType > > getComm()
void f()
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
RCP< VectorBase< Scalar > > createVector(const RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVector, const RCP< const VectorSpaceBase< Scalar > > space=Teuchos::null)
RCP< LinearOpBase< Scalar > > createLinearOp(const RCP< Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraOperator, const RCP< const VectorSpaceBase< Scalar > > rangeSpace=Teuchos::null, const RCP< const VectorSpaceBase< Scalar > > domainSpace=Teuchos::null)