20 #include "Teuchos_LAPACK.hpp"
48 void update(std::vector<Real> &u,
const std::vector<Real> &s,
const Real alpha=1.0)
const {
49 for (
unsigned i=0; i<u.size(); i++) {
54 void axpy(std::vector<Real> &out,
const Real a,
const std::vector<Real> &x,
const std::vector<Real> &y)
const {
55 for (
unsigned i=0; i < x.size(); i++) {
56 out[i] = a*x[i] + y[i];
60 void scale(std::vector<Real> &u,
const Real alpha=0.0)
const {
61 for (
unsigned i=0; i<u.size(); i++) {
66 void linear_solve(std::vector<Real> &u, std::vector<Real> &dl, std::vector<Real> &d, std::vector<Real> &du,
67 const std::vector<Real> &r,
const bool transpose =
false)
const {
68 if ( r.size() == 1 ) {
69 u.resize(1,r[0]/d[0]);
71 else if ( r.size() == 2 ) {
73 Real det = d[0]*d[1] - dl[0]*du[0];
74 u[0] = (d[1]*r[0] - du[0]*r[1])/det;
75 u[1] = (d[0]*r[1] - dl[0]*r[0])/det;
78 u.assign(r.begin(),r.end());
80 Teuchos::LAPACK<int,Real> lp;
81 std::vector<Real> du2(r.size()-2,0.0);
82 std::vector<int> ipiv(r.size(),0);
87 lp.GTTRF(dim,&dl[0],&d[0],&du[0],&du2[0],&ipiv[0],&info);
92 lp.GTTRS(trans,dim,nhrs,&dl[0],&d[0],&du[0],&du2[0],&ipiv[0],&u[0],ldb,&info);
97 BurgersFEM(
int nx = 128, Real nl = 1.0, Real cH1 = 1.0, Real cL2 = 1.0)
106 nu_ = std::pow(10.0,nu-2.0);
124 Real
compute_L2_dot(
const std::vector<Real> &x,
const std::vector<Real> &y)
const {
126 Real c = (((int)x.size()==
nx_) ? 4.0 : 2.0);
127 for (
unsigned i=0; i<x.size(); i++) {
129 ip +=
dx_/6.0*(c*x[i] + x[i+1])*y[i];
131 else if ( i == x.size()-1 ) {
132 ip +=
dx_/6.0*(x[i-1] + c*x[i])*y[i];
135 ip +=
dx_/6.0*(x[i-1] + 4.0*x[i] + x[i+1])*y[i];
147 void apply_mass(std::vector<Real> &Mu,
const std::vector<Real> &u )
const {
148 Mu.resize(u.size(),0.0);
149 Real c = (((int)u.size()==
nx_) ? 4.0 : 2.0);
150 for (
unsigned i=0; i<u.size(); i++) {
152 Mu[i] =
dx_/6.0*(c*u[i] + u[i+1]);
154 else if ( i == u.size()-1 ) {
155 Mu[i] =
dx_/6.0*(u[i-1] + c*u[i]);
158 Mu[i] =
dx_/6.0*(u[i-1] + 4.0*u[i] + u[i+1]);
165 unsigned nx = u.size();
167 std::vector<Real> dl(nx-1,
dx_/6.0);
168 std::vector<Real> d(nx,2.0*
dx_/3.0);
169 std::vector<Real> du(nx-1,
dx_/6.0);
170 if ( (
int)nx !=
nx_ ) {
179 std::vector<Real> u(
nx_,0.0), Mu(
nx_,0.0), iMMu(
nx_,0.0), diff(
nx_,0.0);
180 for (
int i = 0; i <
nx_; i++) {
181 u[i] = 2.0*(Real)rand()/(Real)RAND_MAX - 1.0;
185 axpy(diff,-1.0,iMMu,u);
188 outStream <<
"\nTest Inverse State Mass Matrix\n";
189 outStream <<
" ||u - inv(M)Mu|| = " << error <<
"\n";
190 outStream <<
" ||u|| = " << normu <<
"\n";
191 outStream <<
" Relative Error = " << error/normu <<
"\n";
194 u.resize(nx_+2,0.0); Mu.resize(nx_+2,0.0); iMMu.resize(nx_+2,0.0); diff.resize(nx_+2,0.0);
195 for (
int i = 0; i < nx_+2; i++) {
196 u[i] = 2.0*(Real)rand()/(Real)RAND_MAX - 1.0;
200 axpy(diff,-1.0,iMMu,u);
203 outStream <<
"\nTest Inverse Control Mass Matrix\n";
204 outStream <<
" ||z - inv(M)Mz|| = " << error <<
"\n";
205 outStream <<
" ||z|| = " << normu <<
"\n";
206 outStream <<
" Relative Error = " << error/normu <<
"\n";
214 Real
compute_H1_dot(
const std::vector<Real> &x,
const std::vector<Real> &y)
const {
216 for (
int i=0; i<
nx_; i++) {
218 ip +=
cL2_*
dx_/6.0*(4.0*x[i] + x[i+1])*y[i];
219 ip +=
cH1_*(2.0*x[i] - x[i+1])/
dx_*y[i];
221 else if ( i == nx_-1 ) {
222 ip +=
cL2_*
dx_/6.0*(x[i-1] + 4.0*x[i])*y[i];
223 ip +=
cH1_*(2.0*x[i] - x[i-1])/
dx_*y[i];
226 ip +=
cL2_*
dx_/6.0*(x[i-1] + 4.0*x[i] + x[i+1])*y[i];
227 ip +=
cH1_*(2.0*x[i] - x[i-1] - x[i+1])/
dx_*y[i];
239 void apply_H1(std::vector<Real> &Mu,
const std::vector<Real> &u )
const {
241 for (
int i=0; i<
nx_; i++) {
243 Mu[i] =
cL2_*
dx_/6.0*(4.0*u[i] + u[i+1])
244 +
cH1_*(2.0*u[i] - u[i+1])/
dx_;
246 else if ( i == nx_-1 ) {
247 Mu[i] =
cL2_*
dx_/6.0*(u[i-1] + 4.0*u[i])
248 +
cH1_*(2.0*u[i] - u[i-1])/
dx_;
251 Mu[i] =
cL2_*
dx_/6.0*(u[i-1] + 4.0*u[i] + u[i+1])
252 +
cH1_*(2.0*u[i] - u[i-1] - u[i+1])/
dx_;
267 std::vector<Real> u(
nx_,0.0), Mu(
nx_,0.0), iMMu(
nx_,0.0), diff(
nx_,0.0);
268 for (
int i = 0; i <
nx_; i++) {
269 u[i] = 2.0*(Real)rand()/(Real)RAND_MAX - 1.0;
273 axpy(diff,-1.0,iMMu,u);
276 outStream <<
"\nTest Inverse State H1 Matrix\n";
277 outStream <<
" ||u - inv(M)Mu|| = " << error <<
"\n";
278 outStream <<
" ||u|| = " << normu <<
"\n";
279 outStream <<
" Relative Error = " << error/normu <<
"\n";
288 const std::vector<Real> &z)
const {
291 for (
int i=0; i<
nx_; i++) {
294 r[i] =
nu_/
dx_*(2.0*u[i]-u[i+1]);
297 r[i] =
nu_/
dx_*(2.0*u[i]-u[i-1]);
300 r[i] =
nu_/
dx_*(2.0*u[i]-u[i-1]-u[i+1]);
304 r[i] +=
nl_*u[i+1]*(u[i]+u[i+1])/6.0;
307 r[i] -=
nl_*u[i-1]*(u[i-1]+u[i])/6.0;
310 r[i] -=
dx_/6.0*(z[i]+4.0*z[i+1]+z[i+2]);
324 std::vector<Real> &d,
325 std::vector<Real> &du,
326 const std::vector<Real> &u)
const {
335 for (
int i=0; i<
nx_; i++) {
337 dl[i] +=
nl_*(-2.0*u[i]-u[i+1])/6.0;
338 d[i] +=
nl_*u[i+1]/6.0;
341 d[i] -=
nl_*u[i-1]/6.0;
342 du[i-1] +=
nl_*(u[i-1]+2.0*u[i])/6.0;
352 const std::vector<Real> &v,
353 const std::vector<Real> &u,
354 const std::vector<Real> &z)
const {
356 for (
int i = 0; i <
nx_; i++) {
359 jv[i] += -
nu_/
dx_*v[i-1]-
nl_*(u[i-1]/6.0*v[i]+(u[i]+2.0*u[i-1])/6.0*v[i-1]);
362 jv[i] += -
nu_/
dx_*v[i+1]+
nl_*(u[i+1]/6.0*v[i]+(u[i]+2.0*u[i+1])/6.0*v[i+1]);
365 jv[ 0] -=
nl_*
u0_/6.0*v[0];
366 jv[nx_-1] +=
nl_*
u1_/6.0*v[nx_-1];
371 const std::vector<Real> &v,
372 const std::vector<Real> &u,
373 const std::vector<Real> &z)
const {
375 std::vector<Real> d(
nx_,0.0);
376 std::vector<Real> dl(
nx_-1,0.0);
377 std::vector<Real> du(
nx_-1,0.0);
385 const std::vector<Real> &v,
386 const std::vector<Real> &u,
387 const std::vector<Real> &z)
const {
389 for (
int i = 0; i <
nx_; i++) {
390 ajv[i] =
nu_/
dx_*2.0*v[i];
392 ajv[i] += -
nu_/
dx_*v[i-1]-
nl_*(u[i-1]/6.0*v[i]
393 -(u[i-1]+2.0*u[i])/6.0*v[i-1]);
396 ajv[i] += -
nu_/
dx_*v[i+1]+
nl_*(u[i+1]/6.0*v[i]
397 -(u[i+1]+2.0*u[i])/6.0*v[i+1]);
400 ajv[ 0] -=
nl_*
u0_/6.0*v[0];
401 ajv[nx_-1] +=
nl_*
u1_/6.0*v[nx_-1];
406 const std::vector<Real> &v,
407 const std::vector<Real> &u,
408 const std::vector<Real> &z)
const {
410 std::vector<Real> d(
nx_,0.0);
411 std::vector<Real> du(
nx_-1,0.0);
412 std::vector<Real> dl(
nx_-1,0.0);
423 const std::vector<Real> &v,
424 const std::vector<Real> &u,
425 const std::vector<Real> &z)
const {
426 for (
int i=0; i<
nx_; i++) {
428 jv[i] = -
dx_/6.0*(v[i]+4.0*v[i+1]+v[i+2]);
434 const std::vector<Real> &v,
435 const std::vector<Real> &u,
436 const std::vector<Real> &z)
const {
437 for (
int i=0; i<
nx_+2; i++) {
439 jv[i] = -
dx_/6.0*v[i];
442 jv[i] = -
dx_/6.0*(4.0*v[i-1]+v[i]);
444 else if ( i == nx_ ) {
445 jv[i] = -
dx_/6.0*(4.0*v[i-1]+v[i-2]);
447 else if ( i == nx_+1 ) {
448 jv[i] = -
dx_/6.0*v[i-2];
451 jv[i] = -
dx_/6.0*(v[i-2]+4.0*v[i-1]+v[i]);
460 const std::vector<Real> &w,
461 const std::vector<Real> &v,
462 const std::vector<Real> &u,
463 const std::vector<Real> &z)
const {
464 for (
int i=0; i<
nx_; i++) {
468 ahwv[i] += (w[i]*v[i+1] - w[i+1]*(2.0*v[i]+v[i+1]))/6.0;
471 ahwv[i] += (w[i-1]*(v[i-1]+2.0*v[i]) - w[i]*v[i-1])/6.0;
477 const std::vector<Real> &w,
478 const std::vector<Real> &v,
479 const std::vector<Real> &u,
480 const std::vector<Real> &z) {
481 ahwv.assign(u.size(),0.0);
484 const std::vector<Real> &w,
485 const std::vector<Real> &v,
486 const std::vector<Real> &u,
487 const std::vector<Real> &z) {
488 ahwv.assign(z.size(),0.0);
491 const std::vector<Real> &w,
492 const std::vector<Real> &v,
493 const std::vector<Real> &u,
494 const std::vector<Real> &z) {
495 ahwv.assign(z.size(),0.0);
502 ROL::Ptr<BurgersFEM<Real>>
fem_;
514 const std::vector<Real>& xval = *ex.
getVector();
516 return fem_->compute_L2_dot(xval,yval);
519 ROL::Ptr<ROL::Vector<Real>>
clone()
const override {
521 return ROL::makePtr<L2VectorPrimal>(
528 dual_vec_ = ROL::makePtr<L2VectorDual<Real>>(
532 fem_->apply_mass(*ROL::constPtrCast<std::vector<Real>>(
dual_vec_->getVector()),
541 ROL::Ptr<BurgersFEM<Real>>
fem_;
553 const std::vector<Real>& xval = *ex.
getVector();
556 std::vector<Real> Mx(dimension,Real(0));
557 fem_->apply_inverse_mass(Mx,xval);
560 val += Mx[i]*yval[i];
565 ROL::Ptr<ROL::Vector<Real>>
clone()
const override {
567 return ROL::makePtr<L2VectorDual>(
574 dual_vec_ = ROL::makePtr<L2VectorPrimal<Real>>(
578 fem_->apply_inverse_mass(*ROL::constPtrCast<std::vector<Real>>(
dual_vec_->getVector()),
587 ROL::Ptr<BurgersFEM<Real>>
fem_;
599 const std::vector<Real>& xval = *ex.
getVector();
601 return fem_->compute_H1_dot(xval,yval);
604 ROL::Ptr<ROL::Vector<Real>>
clone()
const override {
606 return ROL::makePtr<H1VectorPrimal>(
613 dual_vec_ = ROL::makePtr<H1VectorDual<Real>>(
617 fem_->apply_H1(*ROL::constPtrCast<std::vector<Real>>(
dual_vec_->getVector()),
626 ROL::Ptr<BurgersFEM<Real>>
fem_;
638 const std::vector<Real>& xval = *ex.
getVector();
641 std::vector<Real> Mx(dimension,0.0);
642 fem_->apply_inverse_H1(Mx,xval);
645 val += Mx[i]*yval[i];
650 ROL::Ptr<ROL::Vector<Real>>
clone()
const override {
652 return ROL::makePtr<H1VectorDual>(
659 dual_vec_ = ROL::makePtr<H1VectorPrimal<Real>>(
663 fem_->apply_inverse_H1(*ROL::constPtrCast<std::vector<Real>>(
dual_vec_->getVector()),
682 ROL::Ptr<BurgersFEM<Real> >
fem_;
687 const bool useHessian =
true)
692 ROL::Ptr<std::vector<Real> > cp =
694 ROL::Ptr<const std::vector<Real> > up =
696 ROL::Ptr<const std::vector<Real> > zp =
699 fem_->compute_residual(*cp,*up,*zp);
704 ROL::Ptr<std::vector<Real> > jvp =
706 ROL::Ptr<const std::vector<Real> > vp =
708 ROL::Ptr<const std::vector<Real> > up =
710 ROL::Ptr<const std::vector<Real> > zp =
713 fem_->apply_pde_jacobian(*jvp,*vp,*up,*zp);
718 ROL::Ptr<std::vector<Real> > jvp =
720 ROL::Ptr<const std::vector<Real> > vp =
722 ROL::Ptr<const std::vector<Real> > up =
724 ROL::Ptr<const std::vector<Real> > zp =
727 fem_->apply_control_jacobian(*jvp,*vp,*up,*zp);
732 ROL::Ptr<std::vector<Real> > ijvp =
734 ROL::Ptr<const std::vector<Real> > vp =
736 ROL::Ptr<const std::vector<Real> > up =
738 ROL::Ptr<const std::vector<Real> > zp =
741 fem_->apply_inverse_pde_jacobian(*ijvp,*vp,*up,*zp);
746 ROL::Ptr<std::vector<Real> > jvp =
748 ROL::Ptr<const std::vector<Real> > vp =
750 ROL::Ptr<const std::vector<Real> > up =
752 ROL::Ptr<const std::vector<Real> > zp =
755 fem_->apply_adjoint_pde_jacobian(*jvp,*vp,*up,*zp);
760 ROL::Ptr<std::vector<Real> > jvp =
762 ROL::Ptr<const std::vector<Real> > vp =
764 ROL::Ptr<const std::vector<Real> > up =
766 ROL::Ptr<const std::vector<Real> > zp =
769 fem_->apply_adjoint_control_jacobian(*jvp,*vp,*up,*zp);
774 ROL::Ptr<std::vector<Real> > iajvp =
776 ROL::Ptr<const std::vector<Real> > vp =
778 ROL::Ptr<const std::vector<Real> > up =
780 ROL::Ptr<const std::vector<Real> > zp =
783 fem_->apply_inverse_adjoint_pde_jacobian(*iajvp,*vp,*up,*zp);
789 ROL::Ptr<std::vector<Real> > ahwvp =
791 ROL::Ptr<const std::vector<Real> > wp =
793 ROL::Ptr<const std::vector<Real> > vp =
795 ROL::Ptr<const std::vector<Real> > up =
797 ROL::Ptr<const std::vector<Real> > zp =
800 fem_->apply_adjoint_pde_hessian(*ahwvp,*wp,*vp,*up,*zp);
810 ROL::Ptr<std::vector<Real> > ahwvp =
812 ROL::Ptr<const std::vector<Real> > wp =
814 ROL::Ptr<const std::vector<Real> > vp =
816 ROL::Ptr<const std::vector<Real> > up =
818 ROL::Ptr<const std::vector<Real> > zp =
821 fem_->apply_adjoint_control_pde_hessian(*ahwvp,*wp,*vp,*up,*zp);
830 ROL::Ptr<std::vector<Real> > ahwvp =
832 ROL::Ptr<const std::vector<Real> > wp =
834 ROL::Ptr<const std::vector<Real> > vp =
836 ROL::Ptr<const std::vector<Real> > up =
838 ROL::Ptr<const std::vector<Real> > zp =
841 fem_->apply_adjoint_pde_control_hessian(*ahwvp,*wp,*vp,*up,*zp);
850 ROL::Ptr<std::vector<Real> > ahwvp =
852 ROL::Ptr<const std::vector<Real> > wp =
854 ROL::Ptr<const std::vector<Real> > vp =
856 ROL::Ptr<const std::vector<Real> > up =
858 ROL::Ptr<const std::vector<Real> > zp =
861 fem_->apply_adjoint_control_hessian(*ahwvp,*wp,*vp,*up,*zp);
H1VectorPrimal< Real > DualConstraintVector
BurgersFEM(int nx=128, Real nl=1.0, Real cH1=1.0, Real cL2=1.0)
Real compute_H1_dot(const std::vector< Real > &x, const std::vector< Real > &y) const
void applyAdjointJacobian_1(ROL::Vector< Real > &ajv, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the adjoint of the partial constraint Jacobian at , , to the vector . This is the primary inter...
void apply_mass(std::vector< Real > &Mu, const std::vector< Real > &u) const
void axpy(std::vector< Real > &out, const Real a, const std::vector< Real > &x, const std::vector< Real > &y) const
ROL::Ptr< ROL::Vector< Real > > clone() const override
Clone to make a new (uninitialized) vector.
Real dot(const ROL::Vector< Real > &x) const override
Compute where .
void apply_adjoint_pde_control_hessian(std::vector< Real > &ahwv, const std::vector< Real > &w, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z)
const ROL::Vector< Real > & dual() const override
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
int dimension() const
Return dimension of the vector space.
ROL::Ptr< L2VectorDual< Real > > dual_vec_
StdVector(const Ptr< std::vector< Real >> &std_vec)
void applyAdjointHessian_12(ROL::Vector< Real > &ahwv, const ROL::Vector< Real > &w, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the optimization-space derivative of the adjoint of the constraint simulation-space Jacobian at...
Constraint_BurgersControl(const ROL::Ptr< BurgersFEM< Real > > &fem, const bool useHessian=true)
void apply_adjoint_pde_hessian(std::vector< Real > &ahwv, const std::vector< Real > &w, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z) const
ROL::Ptr< BurgersFEM< Real > > fem_
Contains definitions of custom data types in ROL.
Ptr< const std::vector< Element > > getVector() const
void applyAdjointHessian_11(ROL::Vector< Real > &ahwv, const ROL::Vector< Real > &w, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the simulation-space derivative of the adjoint of the constraint simulation-space Jacobian at ...
ROL::Ptr< const std::vector< Real > > getVector() const
ROL::Ptr< BurgersFEM< Real > > fem_
Real compute_H1_norm(const std::vector< Real > &r) const
void applyInverseAdjointJacobian_1(ROL::Vector< Real > &iajv, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the inverse of the adjoint of the partial constraint Jacobian at , , to the vector ...
const ROL::Vector< Real > & dual() const override
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
ROL::Ptr< ROL::Vector< Real > > clone() const override
Clone to make a new (uninitialized) vector.
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
H1VectorDual< Real > PrimalConstraintVector
L2VectorPrimal(const ROL::Ptr< std::vector< Real >> &vec, const ROL::Ptr< BurgersFEM< Real >> &fem)
Real dot(const ROL::Vector< Real > &x) const override
Compute where .
ROL::Ptr< ROL::Vector< Real > > clone() const override
Clone to make a new (uninitialized) vector.
H1VectorPrimal< Real > PrimalStateVector
Real mesh_spacing(void) const
int dimension() const
Return dimension of the vector space.
void test_inverse_mass(std::ostream &outStream=std::cout)
void apply_pde_jacobian(std::vector< Real > &jv, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z) const
H1VectorDual(const ROL::Ptr< std::vector< Real >> &vec, const ROL::Ptr< BurgersFEM< Real >> &fem)
ROL::Ptr< BurgersFEM< Real > > fem_
void applyAdjointHessian_21(ROL::Vector< Real > &ahwv, const ROL::Vector< Real > &w, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the simulation-space derivative of the adjoint of the constraint optimization-space Jacobian at...
Provides the ROL::Vector interface for scalar values, to be used, for example, with scalar constraint...
ROL::Ptr< const std::vector< Real > > getVector() const
ROL::Ptr< const std::vector< Real > > getVector() const
void compute_pde_jacobian(std::vector< Real > &dl, std::vector< Real > &d, std::vector< Real > &du, const std::vector< Real > &u) const
H1VectorDual< Real > DualStateVector
void applyJacobian_1(ROL::Vector< Real > &jv, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the partial constraint Jacobian at , , to the vector .
Real compute_L2_dot(const std::vector< Real > &x, const std::vector< Real > &y) const
void linear_solve(std::vector< Real > &u, std::vector< Real > &dl, std::vector< Real > &d, std::vector< Real > &du, const std::vector< Real > &r, const bool transpose=false) const
L2VectorPrimal< Real > PrimalControlVector
void apply_inverse_mass(std::vector< Real > &Mu, const std::vector< Real > &u) const
L2VectorDual(const ROL::Ptr< std::vector< Real >> &vec, const ROL::Ptr< BurgersFEM< Real >> &fem)
ROL::Ptr< BurgersFEM< Real > > fem_
ROL::Ptr< const std::vector< Real > > getVector() const
ROL::Ptr< BurgersFEM< Real > > fem_
void update(std::vector< Real > &u, const std::vector< Real > &s, const Real alpha=1.0) const
const ROL::Vector< Real > & dual() const override
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
void apply_adjoint_control_hessian(std::vector< Real > &ahwv, const std::vector< Real > &w, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z)
void applyInverseJacobian_1(ROL::Vector< Real > &ijv, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the inverse partial constraint Jacobian at , , to the vector .
void compute_residual(std::vector< Real > &r, const std::vector< Real > &u, const std::vector< Real > &z) const
void set_problem_data(const Real nu, const Real f, const Real u0, const Real u1)
void applyJacobian_2(ROL::Vector< Real > &jv, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the partial constraint Jacobian at , , to the vector .
void apply_control_jacobian(std::vector< Real > &jv, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z) const
Real dot(const ROL::Vector< Real > &x) const override
Compute where .
void test_inverse_H1(std::ostream &outStream=std::cout)
void apply_adjoint_control_jacobian(std::vector< Real > &jv, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z) const
H1VectorPrimal(const ROL::Ptr< std::vector< Real >> &vec, const ROL::Ptr< BurgersFEM< Real >> &fem)
void apply_inverse_pde_jacobian(std::vector< Real > &ijv, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z) const
L2VectorDual< Real > DualControlVector
void applyAdjointHessian_22(ROL::Vector< Real > &ahwv, const ROL::Vector< Real > &w, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the optimization-space derivative of the adjoint of the constraint optimization-space Jacobian ...
void apply_H1(std::vector< Real > &Mu, const std::vector< Real > &u) const
void apply_inverse_H1(std::vector< Real > &Mu, const std::vector< Real > &u) const
void apply_adjoint_pde_jacobian(std::vector< Real > &ajv, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z) const
ROL::Ptr< ROL::Vector< Real > > clone() const override
Clone to make a new (uninitialized) vector.
void apply_inverse_adjoint_pde_jacobian(std::vector< Real > &iajv, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z) const
void apply_adjoint_control_pde_hessian(std::vector< Real > &ahwv, const std::vector< Real > &w, const std::vector< Real > &v, const std::vector< Real > &u, const std::vector< Real > &z)
void applyAdjointJacobian_2(ROL::Vector< Real > &jv, const ROL::Vector< Real > &v, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Apply the adjoint of the partial constraint Jacobian at , , to vector . This is the primary interface...
int dimension() const
Return dimension of the vector space.
Defines the constraint operator interface for simulation-based optimization.
ROL::Ptr< H1VectorPrimal< Real > > dual_vec_
void value(ROL::Vector< Real > &c, const ROL::Vector< Real > &u, const ROL::Vector< Real > &z, Real &tol)
Evaluate the constraint operator at .
void scale(std::vector< Real > &u, const Real alpha=0.0) const
const ROL::Vector< Real > & dual() const override
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
ROL::Ptr< H1VectorDual< Real > > dual_vec_
Real dot(const ROL::Vector< Real > &x) const override
Compute where .
int dimension() const
Return dimension of the vector space.
Real compute_L2_norm(const std::vector< Real > &r) const
ROL::Ptr< L2VectorPrimal< Real > > dual_vec_