PDE-constrained inverse problem based on a 2D discretization of a diffusion/reaction system. More...

#include <GLpApp_AdvDiffReactOptModel.hpp>

Inheritance diagram for GLpApp::AdvDiffReactOptModel:

Public Member Functions
	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.
void	set_q (Teuchos::RCP< const Epetra_Vector > const &q)
Teuchos::RCP< GLpApp::GLpYUEpetraDataPool >	getDataPool ()
Teuchos::RCP< const Epetra_MultiVector >	get_B_bar () const
virtual	~ModelEvaluator ()
virtual Teuchos::RCP< const Teuchos::Array< std::string > >	get_p_names (int l) const
	Get the names of the parameters associated with parameter subvector l if available.
virtual Teuchos::ArrayView< const std::string >	get_g_names (int j) const
	Get the names of the response functions associated with response subvector j if available.
virtual Teuchos::RCP< const Epetra_Vector >	get_x_dot_init () const
virtual Teuchos::RCP< const Epetra_Vector >	get_x_dotdot_init () const
virtual double	get_t_init () const
virtual double	getInfBound () const
	Return the value of an infinite bound.
virtual double	get_t_lower_bound () const
virtual double	get_t_upper_bound () const
virtual Teuchos::RCP< EpetraExt::ModelEvaluator::Preconditioner >	create_WPrec () const
virtual Teuchos::RCP< Epetra_Operator >	create_DgDx_dot_op (int j) const
virtual Teuchos::RCP< Epetra_Operator >	create_DgDx_dotdot_op (int j) const
virtual Teuchos::RCP< Epetra_Operator >	create_DgDx_op (int j) const
virtual Teuchos::RCP< Epetra_Operator >	create_DgDp_op (int j, int l) const

Overridden from EpetraExt::ModelEvaluator .
Teuchos::RCP< const Epetra_Map >	get_x_map () const
Teuchos::RCP< const Epetra_Map >	get_f_map () const
Teuchos::RCP< const Epetra_Map >	get_p_map (int l) const
	\breif .
Teuchos::RCP< const Epetra_Map >	get_g_map (int j) const
	\breif .
Teuchos::RCP< const Epetra_Vector >	get_x_init () const
Teuchos::RCP< const Epetra_Vector >	get_p_init (int l) const
Teuchos::RCP< const Epetra_Vector >	get_x_lower_bounds () const
Teuchos::RCP< const Epetra_Vector >	get_x_upper_bounds () const
Teuchos::RCP< const Epetra_Vector >	get_p_lower_bounds (int l) const
Teuchos::RCP< const Epetra_Vector >	get_p_upper_bounds (int l) const
Teuchos::RCP< Epetra_Operator >	create_W () const
Teuchos::RCP< Epetra_Operator >	create_DfDp_op (int l) const
InArgs	createInArgs () const
OutArgs	createOutArgs () const
void	evalModel (const InArgs &inArgs, const OutArgs &outArgs) const

Additional Inherited Members
Public Types inherited from EpetraExt::ModelEvaluator
enum	EInArgsMembers { IN_ARG_x_dot , IN_ARG_x , IN_ARG_x_dot_poly , IN_ARG_x_poly , IN_ARG_x_dot_sg , IN_ARG_x_sg , IN_ARG_x_dot_mp , IN_ARG_x_mp , IN_ARG_t , IN_ARG_alpha , IN_ARG_beta , IN_ARG_step_size , IN_ARG_stage_number , IN_ARG_x_dotdot , IN_ARG_x_dotdot_poly , IN_ARG_x_dotdot_sg , IN_ARG_x_dotdot_mp , IN_ARG_omega , IN_ARG_sg_basis , IN_ARG_sg_quadrature , IN_ARG_sg_expansion }
enum	EInArgs_p_sg { IN_ARG_p_sg }
enum	EInArgs_p_mp { IN_ARG_p_mp }
enum	EEvalType { EVAL_TYPE_EXACT , EVAL_TYPE_APPROX_DERIV , EVAL_TYPE_VERY_APPROX_DERIV }
enum	EDerivativeMultiVectorOrientation { DERIV_MV_BY_COL , DERIV_TRANS_MV_BY_ROW }
enum	EDerivativeLinearOp { DERIV_LINEAR_OP }
enum	EDerivativeLinearity { DERIV_LINEARITY_UNKNOWN , DERIV_LINEARITY_CONST , DERIV_LINEARITY_NONCONST }
enum	ERankStatus { DERIV_RANK_UNKNOWN , DERIV_RANK_FULL , DERIV_RANK_DEFICIENT }
enum	EOutArgsMembers { OUT_ARG_f , OUT_ARG_W , OUT_ARG_f_poly , OUT_ARG_f_sg , OUT_ARG_W_sg , OUT_ARG_f_mp , OUT_ARG_W_mp , OUT_ARG_WPrec }
enum	EOutArgsDfDp { OUT_ARG_DfDp }
enum	EOutArgsDgDx_dot { OUT_ARG_DgDx_dot }
enum	EOutArgsDgDx_dotdot { OUT_ARG_DgDx_dotdot }
enum	EOutArgsDgDx { OUT_ARG_DgDx }
enum	EOutArgsDgDp { OUT_ARG_DgDp }
enum	EOutArgsDfDp_sg { OUT_ARG_DfDp_sg }
enum	EOutArgs_g_sg { OUT_ARG_g_sg }
enum	EOutArgsDgDx_dot_sg { OUT_ARG_DgDx_dot_sg }
enum	EOutArgsDgDx_dotdot_sg { OUT_ARG_DgDx_dotdot_sg }
enum	EOutArgsDgDx_sg { OUT_ARG_DgDx_sg }
enum	EOutArgsDgDp_sg { OUT_ARG_DgDp_sg }
enum	EOutArgsDfDp_mp { OUT_ARG_DfDp_mp }
enum	EOutArgs_g_mp { OUT_ARG_g_mp }
enum	EOutArgsDgDx_dot_mp { OUT_ARG_DgDx_dot_mp }
enum	EOutArgsDgDx_dotdot_mp { OUT_ARG_DgDx_dotdot_mp }
enum	EOutArgsDgDx_mp { OUT_ARG_DgDx_mp }
enum	EOutArgsDgDp_mp { OUT_ARG_DgDp_mp }
typedef Teuchos::RCP< const Stokhos::ProductEpetraVector >	mp_const_vector_t
typedef Teuchos::RCP< const Stokhos::ProductEpetraMultiVector >	mp_const_multivector_t
typedef Teuchos::RCP< const Stokhos::ProductEpetraOperator >	mp_const_operator_t
typedef Teuchos::RCP< Stokhos::ProductEpetraVector >	mp_vector_t
typedef Teuchos::RCP< Stokhos::ProductEpetraMultiVector >	mp_multivector_t
typedef Teuchos::RCP< Stokhos::ProductEpetraOperator >	mp_operator_t
static const int	NUM_E_IN_ARGS_MEMBERS =21
static const int	NUM_E_OUT_ARGS_MEMBERS =9

Detailed Description

PDE-constrained inverse problem based on a 2D discretization of a diffusion/reaction system.

The model evaluator subclass is used to represent the simulation-constrained optimization problem:

 min   g(x,p)
 s.t.  f(x,p) = 0;

where:

x is the vector of discretized concentations of the species in the 2D domain.
p is the global vector of coefficients of a sine series basis (see B_bar below).
f(x,p) = A*x + reationRate*Ny(x) + B*(B_bar*p) is the discretized 2D diffusion/reaction PDE.
g(x,p) = 0.5 * (x-q)'H(x-q) + 0.5*regBeta*(B_bar*p)'R(B_bar*p) is the least squares objective function.
A is the discretized Laplacian operator for the diffusion part of the PDE state equation. This matrix is constant, square and singular.
B is the sensitivity of the flux boundary conditions. This is a constant rectangular matrix.
B_bar are the sine series coefficients with a column dimension of np.
Ny(x) is the nonlinear terms for the discretized reaction over the 2D domain.
reactionRate is the relative reaction rate which must take on a non-zero value to form a solvable problem.
H is the symmetric positive definite mass matrix for the problem (i.e. the discretization of the inner product operator over the 2D domain).
q is a matching or target vector for the state x over the 2D domain of the problem.
R is the symmetric positive definite discretization of the inner product of the flux function over the boundary of the 2D domain.
regBeta is a regularization parameter that must be greater than zero.

The nuts and bolts of the implementation for this problem are contained in the C++ class GLpApp::GLpYUEpetraDataPool that was originally implemented by Denis Ridzal while a student at Rice University. The class GLpApp::GLpYUEpetraDataPool implements the basic operators and nonlinear functions but this class puts them together to form a valid the model in terms of a model evaluator interface.

This example problem demonstrates a few different aspects of the EpetraExt::ModelEvaluator interface:

How to manage parallel vector data. The state variables in x are managed as fully distributed parallel data while the flux sine-series parameter coefficients p are managed as locally replicated data.

Demonstrates shared compuation between the objective function g(x,p) and the simulation equality constraints f(x,p) and their derivatives. The intermediate vector B_bar*p is computed only once and is shared with the computation of g and f. The intermediate vector R*(B_bar*p) is computed once and shared between the computation of g and DgDp.

The functions AdvDiffReactOptModel() createInArgs(), createOutArgs() and evalModel() are fairly cleanly written and are appropriate to be studied in order to show how to implement other parallel simulation-constrained problems based on Epetra objects.

The mesh for the 2D domain can either be read in as a mesh data file give the files name or can be generated automatically on a square 2D domain.

The program triangle can be used to generate meshes for arbitary 2D geometries and then metis can be used to partition the mesh to multiple domains. Instructions for how to use triangle and metis to generate meshes is described ???here???.

Instead of reading in a mesh file, a square 2D mesh can be automatically generated given just the length in the x and y directions and the number of local elements in each direction. Currently, the square mesh is only partitioned in the x direction and therefore will not demonstrate great parallel scalability for large numbers of processors due to excessive amounts of shared boundary between processes.

ToDo: Finish Documentation!

Definition at line 162 of file GLpApp_AdvDiffReactOptModel.hpp.

Constructor & Destructor Documentation

◆ AdvDiffReactOptModel()

GLpApp::AdvDiffReactOptModel::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.

Definition at line 95 of file GLpApp_AdvDiffReactOptModel.cpp.

Member Function Documentation

◆ set_q()

void GLpApp::AdvDiffReactOptModel::set_q ( Teuchos::RCP< const Epetra_Vector > const & q )

Definition at line 246 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ getDataPool()

Teuchos::RCP< GLpApp::GLpYUEpetraDataPool > GLpApp::AdvDiffReactOptModel::getDataPool ( )

Definition at line 259 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_B_bar()

Teuchos::RCP< const Epetra_MultiVector > GLpApp::AdvDiffReactOptModel::get_B_bar ( ) const

Definition at line 265 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_x_map()

Teuchos::RCP< const Epetra_Map > GLpApp::AdvDiffReactOptModel::get_x_map ( ) const

virtual

Implements EpetraExt::ModelEvaluator.

Definition at line 273 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_f_map()

Teuchos::RCP< const Epetra_Map > GLpApp::AdvDiffReactOptModel::get_f_map ( ) const

virtual

Implements EpetraExt::ModelEvaluator.

Definition at line 279 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_p_map()

Teuchos::RCP< const Epetra_Map > GLpApp::AdvDiffReactOptModel::get_p_map ( int l ) const

virtual

\breif .

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 285 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_g_map()

Teuchos::RCP< const Epetra_Map > GLpApp::AdvDiffReactOptModel::get_g_map ( int j ) const

virtual

\breif .

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 292 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_x_init()

Teuchos::RCP< const Epetra_Vector > GLpApp::AdvDiffReactOptModel::get_x_init ( ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 299 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_p_init()

Teuchos::RCP< const Epetra_Vector > GLpApp::AdvDiffReactOptModel::get_p_init ( int l ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 305 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_x_lower_bounds()

Teuchos::RCP< const Epetra_Vector > GLpApp::AdvDiffReactOptModel::get_x_lower_bounds ( ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 312 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_x_upper_bounds()

Teuchos::RCP< const Epetra_Vector > GLpApp::AdvDiffReactOptModel::get_x_upper_bounds ( ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 318 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_p_lower_bounds()

Teuchos::RCP< const Epetra_Vector > GLpApp::AdvDiffReactOptModel::get_p_lower_bounds ( int l ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 324 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ get_p_upper_bounds()

Teuchos::RCP< const Epetra_Vector > GLpApp::AdvDiffReactOptModel::get_p_upper_bounds ( int l ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 331 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ create_W()

Teuchos::RCP< Epetra_Operator > GLpApp::AdvDiffReactOptModel::create_W ( ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 338 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ create_DfDp_op()

Teuchos::RCP< Epetra_Operator > GLpApp::AdvDiffReactOptModel::create_DfDp_op ( int l ) const

virtual

Reimplemented from EpetraExt::ModelEvaluator.

Definition at line 344 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ createInArgs()

EpetraExt::ModelEvaluator::InArgs GLpApp::AdvDiffReactOptModel::createInArgs ( ) const

virtual

Implements EpetraExt::ModelEvaluator.

Definition at line 352 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ createOutArgs()

EpetraExt::ModelEvaluator::OutArgs GLpApp::AdvDiffReactOptModel::createOutArgs ( ) const

virtual

Implements EpetraExt::ModelEvaluator.

Definition at line 362 of file GLpApp_AdvDiffReactOptModel.cpp.

◆ evalModel()

void GLpApp::AdvDiffReactOptModel::evalModel	(	const InArgs &	inArgs,
		const OutArgs &	outArgs ) const

virtual

Implements EpetraExt::ModelEvaluator.

Definition at line 411 of file GLpApp_AdvDiffReactOptModel.cpp.

The documentation for this class was generated from the following files:

Public Member Functions

Overridden from EpetraExt::ModelEvaluator .

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ AdvDiffReactOptModel()

Member Function Documentation

◆ set_q()

◆ getDataPool()

◆ get_B_bar()

◆ get_x_map()

◆ get_f_map()

◆ get_p_map()

◆ get_g_map()

◆ get_x_init()

◆ get_p_init()

◆ get_x_lower_bounds()

◆ get_x_upper_bounds()

◆ get_p_lower_bounds()

◆ get_p_upper_bounds()

◆ create_W()

◆ create_DfDp_op()

◆ createInArgs()

◆ createOutArgs()

◆ evalModel()