Thyra
Version of the Day
|
Simple parallel response-only ModelEvaluator. More...
#include <Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_decl.hpp>
Constructors/Initializers/Accessors. | |
DiagonalQuadraticResponseOnlyModelEvaluator (const int localDim, const RCP< const Teuchos::Comm< Ordinal > > &comm=Teuchos::null) | |
void | setSolutionVector (const RCP< const VectorBase< Scalar > > &ps) |
Set the solution vector ps . More... | |
const RCP< const VectorBase < Scalar > > | getSolutionVector () const |
Get the solution vector ps . More... | |
void | setDiagonalVector (const RCP< const VectorBase< Scalar > > &diag) |
Set the diagonal vector diag. More... | |
void | setDiagonalBarVector (const RCP< const VectorBase< Scalar > > &diag_bar) |
Set the diagonal vector diag_bar. More... | |
void | setNonlinearTermFactor (const Scalar &nonlinearTermFactor) |
Set nonlinear term factory. More... | |
void | setScalarOffset (const Scalar &g_offset) |
Set offset scalar g_offset . More... | |
Public functions overridden from ModelEvaulator. | |
int | Np () const |
int | Ng () const |
RCP< const VectorSpaceBase < Scalar > > | get_p_space (int l) const |
RCP< const VectorSpaceBase < Scalar > > | get_g_space (int j) const |
ModelEvaluatorBase::InArgs < Scalar > | createInArgs () const |
Simple parallel response-only ModelEvaluator.
The representation of the model in coefficient form is:
g(p) = 0.5 * sum( diag[i] * (p[i] - ps[i])^2, i=0...n-1 ) + 0.5 * nonlinearTermFactor * sum( (p[i] - ps[i])^3, i=0...n-1 ) + g_offset
In vector coefficient form it becomes:
g(p) = 0.5 * (p-ps)^T * D * (p-ps) + cubicTerm(p) + g_offset
where D = diag(diag)
and cubitTerm(p)
is given above.
This test model implementation also supports the definition of diagonal scalar product implementation. The impact of this is the definition of new abstract Thyra Euclidean vectors of the form:
p_bar = E_p * p g_grad_bar = E_p * g_grad
In this notation, we say that the Thyra vectors themselves p_bar
are vectors in a Euclidean space (from the standpoint of the ANA clients) and that p
are the coefficient vectors in the new non-Euclidean space S_p
with scalar product:
y_bar^T * x_bar = y^T * E_p^T * E_p * x = y^T * S_p * x = <y, x>_p
The ideas is that by providing this new vector space basis E_p
and inner product <.,.>_p
, we change the function that the ANA sees.
This testing class allows a client to specify the desirable diagonal D_bar
matrix by passing in the vector diag_bar
. From diag_bar
and diag
, the diagonal s_diag
for S_p = diag(s_diag)
as:
s_diag[i] = diag_bar[i] / diag[i]
The inner product is therefore:
scalarProd(y_bar, x_bar) = sum( s_diag[i] * y[i] * x[i], i=0...n-1 )
Definition at line 126 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_decl.hpp.
Thyra::DiagonalQuadraticResponseOnlyModelEvaluator< Scalar >::DiagonalQuadraticResponseOnlyModelEvaluator | ( | const int | localDim, |
const RCP< const Teuchos::Comm< Ordinal > > & | comm = Teuchos::null |
||
) |
Definition at line 70 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
void Thyra::DiagonalQuadraticResponseOnlyModelEvaluator< Scalar >::setSolutionVector | ( | const RCP< const VectorBase< Scalar > > & | ps | ) |
Set the solution vector ps .
Definition at line 117 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
const RCP< const Thyra::VectorBase< Scalar > > Thyra::DiagonalQuadraticResponseOnlyModelEvaluator< Scalar >::getSolutionVector | ( | ) | const |
Get the solution vector ps .
Definition at line 126 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
void Thyra::DiagonalQuadraticResponseOnlyModelEvaluator< Scalar >::setDiagonalVector | ( | const RCP< const VectorBase< Scalar > > & | diag | ) |
Set the diagonal vector diag.
Definition at line 133 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
void Thyra::DiagonalQuadraticResponseOnlyModelEvaluator< Scalar >::setDiagonalBarVector | ( | const RCP< const VectorBase< Scalar > > & | diag_bar | ) |
Set the diagonal vector diag_bar.
NOTE: You must call setDiagonalVector(diag) first in order to set the objective diagonal.
Definition at line 142 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
void Thyra::DiagonalQuadraticResponseOnlyModelEvaluator< Scalar >::setNonlinearTermFactor | ( | const Scalar & | nonlinearTermFactor | ) |
Set nonlinear term factory.
Definition at line 174 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
void Thyra::DiagonalQuadraticResponseOnlyModelEvaluator< Scalar >::setScalarOffset | ( | const Scalar & | g_offset | ) |
Set offset scalar g_offset .
Definition at line 182 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
|
virtual |
Reimplemented from Thyra::ModelEvaluatorDefaultBase< Scalar >.
Definition at line 193 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
|
virtual |
Reimplemented from Thyra::ModelEvaluatorDefaultBase< Scalar >.
Definition at line 200 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
|
virtual |
Implements Thyra::ModelEvaluator< Scalar >.
Definition at line 208 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
|
virtual |
Implements Thyra::ModelEvaluator< Scalar >.
Definition at line 221 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.
|
virtual |
Implements Thyra::ModelEvaluator< Scalar >.
Definition at line 234 of file Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_def.hpp.