53 #ifndef ROL_LEASTSQUARES_HPP
54 #define ROL_LEASTSQUARES_HPP
70 Teuchos::RCP<const std::vector<Real> > xp = ex.
getVector();
73 Real h = 1.0/((Real)n+1.0);
76 for (
int i=0; i<n; i++) {
78 res = 2.0*h*(5.0/6.0) + 1.0/h*((*xp)[i+1]-2.0*(*xp)[i]);
80 else if ( i == n-1 ) {
81 res = 2.0*h*(5.0/6.0) + 1.0/h*((*xp)[i-1]-2.0*(*xp)[i]);
84 res = 2.0*h + 1.0/h*((*xp)[i-1]-2.0*(*xp)[i]+(*xp)[i+1]);
93 Teuchos::RCP<const std::vector<Real> > xp =
95 Teuchos::RCP<std::vector<Real> > gp =
96 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(g)).getVector());
99 Real h = 1.0/((Real)n+1.0);
100 std::vector<Real> res(n,0.0);
101 for (
int i=0; i<n; i++) {
103 res[i] = 2.0*h*(5.0/6.0) + 1.0/h*((*xp)[i+1]-2.0*(*xp)[i]);
105 else if ( i == n-1 ) {
106 res[i] = 2.0*h*(5.0/6.0) + 1.0/h*((*xp)[i-1]-2.0*(*xp)[i]);
109 res[i] = 2.0*h + 1.0/h*((*xp)[i-1]-2.0*(*xp)[i]+(*xp)[i+1]);
113 for (
int i=0; i<n; i++) {
115 (*gp)[i] = 1.0/h*(res[i+1]-2.0*res[i]);
117 else if ( i == n-1 ) {
118 (*gp)[i] = 1.0/h*(res[i-1]-2.0*res[i]);
121 (*gp)[i] = 1.0/h*(res[i-1]-2.0*res[i]+res[i+1]);
127 Teuchos::RCP<const std::vector<Real> > xp =
129 Teuchos::RCP<const std::vector<Real> > vp =
131 Teuchos::RCP<std::vector<Real> > hvp =
132 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(hv)).getVector());
135 Real h = 1.0/((Real)n+1.0);
136 std::vector<Real> res(n,0.0);
137 for (
int i=0; i<n; i++) {
139 res[i] = 1.0/h*((*vp)[i+1]-2.0*(*vp)[i]);
141 else if ( i == n-1 ) {
142 res[i] = 1.0/h*((*vp)[i-1]-2.0*(*vp)[i]);
145 res[i] = 1.0/h*((*vp)[i-1]-2.0*(*vp)[i]+(*vp)[i+1]);
149 for (
int i=0; i<n; i++) {
151 (*hvp)[i] = 1.0/h*(res[i+1]-2.0*res[i]);
153 else if ( i == n-1 ) {
154 (*hvp)[i] = 1.0/h*(res[i-1]-2.0*res[i]);
157 (*hvp)[i] = 1.0/h*(res[i-1]-2.0*res[i]+res[i+1]);
167 Teuchos::RCP<std::vector<Real> > x0p =
168 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(x0)).getVector());
169 Teuchos::RCP<std::vector<Real> > xp =
170 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(x)).getVector());
179 for (
int i=0; i<n; i++) {
183 Real h = 1.0/((Real)n+1.0);
185 for(
int i=0; i<n; i++ ) {
187 (*xp)[i] = pt*(1.0-pt);
Provides the interface to evaluate objective functions.
Teuchos::RCP< const std::vector< Element > > getVector() const
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Defines the linear algebra or vector space interface.
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Provides the std::vector implementation of the ROL::Vector interface.
void getLeastSquares(Teuchos::RCP< Objective< Real > > &obj, Vector< Real > &x0, Vector< Real > &x)