EpetraExt Development
Loading...
Searching...
No Matches
GLpApp_AdvDiffReactOptModel.hpp
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
44#ifndef GLP_APP_ADV_DIFF_REACT_OPT_MODEL_HPP
45#define GLP_APP_ADV_DIFF_REACT_OPT_MODEL_HPP
46
49#include "Epetra_Map.h"
50#include "Epetra_Vector.h"
51#include "Epetra_Comm.h"
52#include "Epetra_CrsGraph.h"
53#include "Teuchos_VerboseObject.hpp"
54#include "Teuchos_Array.hpp"
55
56namespace GLpApp {
57
164 , public Teuchos::VerboseObject<AdvDiffReactOptModel>
165{
166public:
167
170 const Teuchos::RCP<const Epetra_Comm> &comm
171 ,const double beta
172 ,const double len_x // Ignored if meshFile is *not* empty
173 ,const double len_y // Ignored if meshFile is *not* empty
174 ,const int local_nx // Ignored if meshFile is *not* empty
175 ,const int local_ny // Ignored if meshFile is *not* empty
176 ,const char meshFile[]
177 ,const int np
178 ,const double x0
179 ,const double p0
180 ,const double reactionRate
181 ,const bool normalizeBasis
182 ,const bool supportDerivatives
183 );
184
186 void set_q( Teuchos::RCP<const Epetra_Vector> const& q );
187
189 Teuchos::RCP<GLpApp::GLpYUEpetraDataPool> getDataPool();
190
192 Teuchos::RCP<const Epetra_MultiVector> get_B_bar() const;
193
196
198 Teuchos::RCP<const Epetra_Map> get_x_map() const;
200 Teuchos::RCP<const Epetra_Map> get_f_map() const;
202 Teuchos::RCP<const Epetra_Map> get_p_map(int l) const;
204 Teuchos::RCP<const Epetra_Map> get_g_map(int j) const;
206 Teuchos::RCP<const Epetra_Vector> get_x_init() const;
208 Teuchos::RCP<const Epetra_Vector> get_p_init(int l) const;
210 Teuchos::RCP<const Epetra_Vector> get_x_lower_bounds() const;
212 Teuchos::RCP<const Epetra_Vector> get_x_upper_bounds() const;
214 Teuchos::RCP<const Epetra_Vector> get_p_lower_bounds(int l) const;
216 Teuchos::RCP<const Epetra_Vector> get_p_upper_bounds(int l) const;
218 Teuchos::RCP<Epetra_Operator> create_W() const;
220 Teuchos::RCP<Epetra_Operator> create_DfDp_op(int l) const;
222 InArgs createInArgs() const;
224 OutArgs createOutArgs() const;
226 void evalModel( const InArgs& inArgs, const OutArgs& outArgs ) const;
227
229
230private:
231
232 // /////////////////////////////////////
233 // Private types
234
235 typedef Teuchos::Array<Teuchos::RCP<const Epetra_Map> > RCP_Eptra_Map_Array_t;
236 typedef Teuchos::Array<Teuchos::RCP<Epetra_Vector> > RCP_Eptra_Vector_Array_t;
237
238 // /////////////////////////////////////
239 // Private member data
240
241 static const int Np_ = 2; // Number of axiliary parameters
242 static const int p_bndy_idx = 0; // index for boundary flux parameters
243 static const int p_rx_idx = 1; // index for reaction rate parameter
244
245 bool supportDerivatives_;
246
247 bool isInitialized_;
248
249 Teuchos::RCP<GLpApp::GLpYUEpetraDataPool> dat_;
250 int np_;
251 Teuchos::RCP<const Epetra_Vector> q_;
252
253 Teuchos::RCP<const Epetra_Map> map_p_bar_;
254 Teuchos::RCP<Epetra_MultiVector> B_bar_;
255
256 Teuchos::RCP<const Epetra_Comm> epetra_comm_;
257 Teuchos::RCP<const Epetra_Map> map_x_;
258 RCP_Eptra_Map_Array_t map_p_;
259 Teuchos::RCP<const Epetra_Map> map_f_;
260 Teuchos::RCP<const Epetra_Map> map_g_;
261
262 Teuchos::RCP<Epetra_Vector> x0_;
263 Teuchos::RCP<Epetra_Vector> xL_;
264 Teuchos::RCP<Epetra_Vector> xU_;
265 RCP_Eptra_Vector_Array_t p0_;
266 RCP_Eptra_Vector_Array_t pL_;
267 RCP_Eptra_Vector_Array_t pU_;
268 Teuchos::RCP<Epetra_Vector> gL_;
269 Teuchos::RCP<Epetra_Vector> gU_;
270
271 Teuchos::RCP<Epetra_CrsGraph> W_graph_;
272
273};
274
275} // namespace GLpApp
276
277#endif // GLP_APP_ADV_DIFF_REACT_OPT_MODEL_HPP
Base interface for evaluating a stateless "model".
PDE-constrained inverse problem based on a 2D discretization of a diffusion/reaction system.
Teuchos::RCP< Epetra_Operator > create_W() const
Teuchos::RCP< const Epetra_Map > get_p_map(int l) const
\breif .
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Teuchos::RCP< const Epetra_Vector > get_p_upper_bounds(int l) const
Teuchos::RCP< const Epetra_Map > get_g_map(int j) const
\breif .
Teuchos::RCP< GLpApp::GLpYUEpetraDataPool > getDataPool()
Teuchos::RCP< const Epetra_Vector > get_x_upper_bounds() const
Teuchos::RCP< const Epetra_Vector > get_p_lower_bounds(int l) const
Teuchos::RCP< Epetra_Operator > create_DfDp_op(int l) const
Teuchos::RCP< const Epetra_MultiVector > get_B_bar() const
Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Teuchos::RCP< const Epetra_Map > get_f_map() const
void set_q(Teuchos::RCP< const Epetra_Vector > const &q)
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
AdvDiffReactOptModel(const Teuchos::RCP< const Epetra_Comm > &comm, const double beta, const double len_x, const double len_y, const int local_nx, const int local_ny, const char meshFile[], const int np, const double x0, const double p0, const double reactionRate, const bool normalizeBasis, const bool supportDerivatives)
Constructor.
Teuchos::RCP< const Epetra_Vector > get_x_lower_bounds() const
Teuchos::RCP< const Epetra_Map > get_x_map() const