44 #ifndef OPTIPACK_DIAGONAL_QUADRATIC_RESPONSE_ONLY_MODEL_EVALUATOR_DEF_HPP 45 #define OPTIPACK_DIAGONAL_QUADRATIC_RESPONSE_ONLY_MODEL_EVALUATOR_DEF_HPP 48 #include "Thyra_DiagonalQuadraticResponseOnlyModelEvaluator_decl.hpp" 49 #include "Thyra_DiagonalScalarProd.hpp" 50 #include "Thyra_VectorStdOps.hpp" 51 #include "Thyra_DefaultSpmdVectorSpace.hpp" 52 #include "Thyra_DetachedSpmdVectorView.hpp" 53 #include "Teuchos_DefaultComm.hpp" 54 #include "Teuchos_CommHelpers.hpp" 55 #include "Teuchos_Assert.hpp" 69 template<
class Scalar>
72 const RCP<
const Teuchos::Comm<Thyra::Ordinal> > &comm
74 :Np_(1), Ng_(1), comm_(comm), localDim_(localDim),
75 nonlinearTermFactor_(0.0), g_offset_(0.0)
78 typedef ScalarTraits<Scalar> ST;
79 using Thyra::createMember;
81 TEUCHOS_ASSERT( localDim > 0 );
85 comm_ = Teuchos::DefaultComm<Thyra::Ordinal>::getComm();
89 g_space_ = Thyra::locallyReplicatedDefaultSpmdVectorSpace<Scalar>(comm_, 1);
92 p_space_ = Thyra::defaultSpmdVectorSpace<Scalar>(comm_, localDim, -1);
95 const RCP<Thyra::VectorBase<Scalar> > ps = createMember<Scalar>(p_space_);
96 V_S(ps.ptr(), ST::zero());
100 const RCP<Thyra::VectorBase<Scalar> > diag = createMember<Scalar>(p_space_);
101 V_S(diag.ptr(), ST::one());
106 const RCP<Thyra::VectorBase<Scalar> > s_bar = createMember<Scalar>(p_space_);
107 V_S(s_bar.ptr(), ST::one());
111 g_offset_ = ST::zero();
116 template<
class Scalar>
120 ps_ = ps.assert_not_null();
124 template<
class Scalar>
125 const RCP<const Thyra::VectorBase<Scalar> >
132 template<
class Scalar>
141 template<
class Scalar>
146 typedef Teuchos::ScalarTraits<Scalar> ST;
147 using Teuchos::rcp_dynamic_cast;
148 using Thyra::createMember;
149 using Thyra::ele_wise_divide;
152 diag_bar_ = diag_bar.assert_not_null();
156 RCP<Thyra::VectorBase<Scalar> > s_bar = createMember<Scalar>(p_space_->clone());
162 V_S( s_bar.ptr(), ST::zero() );
163 ele_wise_divide( ST::one(), *diag_, *diag_bar_, s_bar.ptr() );
166 const RCP<Thyra::ScalarProdVectorSpaceBase<Scalar> > sp_p_space =
168 sp_p_space->
setScalarProd(diagonalScalarProd<Scalar>(s_bar_));
173 template<
class Scalar>
175 const Scalar &nonlinearTermFactor)
177 nonlinearTermFactor_ = nonlinearTermFactor;
181 template<
class Scalar>
183 const Scalar &g_offset)
185 g_offset_ = g_offset;
192 template<
class Scalar>
199 template<
class Scalar>
206 template<
class Scalar>
207 RCP<const Thyra::VectorSpaceBase<Scalar> >
211 TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE( l, 0, Np_ );
217 template<
class Scalar>
218 RCP<const Thyra::VectorSpaceBase<Scalar> >
222 TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE( j, 0, Ng_ );
228 template<
class Scalar>
233 MEB::InArgsSetup<Scalar> inArgs;
234 inArgs.setModelEvalDescription(this->description());
243 template<
class Scalar>
248 MEB::OutArgsSetup<Scalar> outArgs;
249 outArgs.setModelEvalDescription(this->description());
250 outArgs.set_Np_Ng(Np_,Ng_);
251 outArgs.setSupports(MEB::OUT_ARG_DgDp, 0 ,0, MEB::DERIV_TRANS_MV_BY_ROW);
256 template<
class Scalar>
257 void DiagonalQuadraticResponseOnlyModelEvaluator<Scalar>::evalModelImpl(
264 using Teuchos::outArg;
265 typedef Teuchos::ScalarTraits<Scalar> ST;
271 const ConstDetachedSpmdVectorView<Scalar> p(inArgs.
get_p(0));
272 const ConstDetachedSpmdVectorView<Scalar> ps(ps_);
273 const ConstDetachedSpmdVectorView<Scalar> diag(diag_);
274 const ConstDetachedSpmdVectorView<Scalar> s_bar(s_bar_);
277 if (!is_null(outArgs.
get_g(0))) {
278 Scalar g_val = ST::zero();
279 for (
Ordinal i = 0; i < p.subDim(); ++i) {
280 const Scalar p_ps = p[i] - ps[i];
281 g_val += diag[i] * p_ps*p_ps;
282 if (nonlinearTermFactor_ != ST::zero()) {
283 g_val += nonlinearTermFactor_ * p_ps * p_ps * p_ps;
287 Teuchos::reduceAll<Ordinal, Scalar>(*comm_, Teuchos::REDUCE_SUM,
288 g_val, outArg(global_g_val) );
289 DetachedSpmdVectorView<Scalar>(outArgs.
get_g(0))[0] =
290 as<Scalar>(0.5) * global_g_val + g_offset_;
294 if (!outArgs.
get_DgDp(0,0).isEmpty()) {
295 const RCP<Thyra::MultiVectorBase<Scalar> > DgDp_trans_mv =
296 get_mv<Scalar>(outArgs.
get_DgDp(0,0),
"DgDp^T", MEB::DERIV_TRANS_MV_BY_ROW);
297 const DetachedSpmdVectorView<Scalar> DgDp_grad(DgDp_trans_mv->col(0));
299 const Scalar p_ps = p[i] - ps[i];
300 Scalar DgDp_grad_i = diag[i] * p_ps;
301 if (nonlinearTermFactor_ != ST::zero()) {
302 DgDp_grad_i += as<Scalar>(1.5) * nonlinearTermFactor_ * p_ps * p_ps;
304 DgDp_grad[i] = DgDp_grad_i / s_bar[i];
315 #endif // OPTIPACK_DIAGONAL_QUADRATIC_RESPONSE_ONLY_MODEL_EVALUATOR_DEF_HPP Create an explicit detached mutable (non-const) view of all of the local elements on this process of ...
Concrete aggregate class for all output arguments computable by a ModelEvaluator subclass object...
const RCP< const VectorBase< Scalar > > getSolutionVector() const
Get the solution vector ps .
Derivative< Scalar > get_DgDp(int j, int l) const
Precondition: supports(OUT_ARG_DgDp,j,l)==true.
void setSolutionVector(const RCP< const VectorBase< Scalar > > &ps)
Set the solution vector ps .
RCP< const VectorSpaceBase< Scalar > > get_p_space(int l) const
RCP< const VectorSpaceBase< Scalar > > get_g_space(int j) const
DiagonalQuadraticResponseOnlyModelEvaluator(const int localDim, const RCP< const Teuchos::Comm< Ordinal > > &comm=Teuchos::null)
void setScalarOffset(const Scalar &g_offset)
Set offset scalar g_offset .
RCP< const VectorBase< Scalar > > get_p(int l) const
Get p(l) where 0 <= l && l < this->Np().
void setNonlinearTermFactor(const Scalar &nonlinearTermFactor)
Set nonlinear term factory.
Teuchos::Ordinal Ordinal
Type for the dimension of a vector space. `*.
Abstract interface for finite-dimensional dense vectors.
Evaluation< VectorBase< Scalar > > get_g(int j) const
Precondition: supports(OUT_ARG_g)==true..
void setDiagonalVector(const RCP< const VectorBase< Scalar > > &diag)
Set the diagonal vector diag.
Base subclass for ModelEvaluator that defines some basic types.
ModelEvaluatorBase::InArgs< Scalar > createInArgs() const
virtual void setScalarProd(const RCP< const ScalarProdBase< Scalar > > &scalarProd)
Set a different scalar product.
void setDiagonalBarVector(const RCP< const VectorBase< Scalar > > &diag_bar)
Set the diagonal vector diag_bar.
Create an explicit detached non-mutable (const) view of all of the local elements on this process of ...
Simple parallel response-only ModelEvaluator.
Concrete aggregate class for all output arguments computable by a ModelEvaluator subclass object...