44 #ifndef ROL_STDVECTOR_H
45 #define ROL_STDVECTOR_H
59 template <
class Real,
class Element=Real>
73 std_vec_ = makePtr<std::vector<Element>>(dim,val);
79 std::invalid_argument,
80 "Error: Vectors must have the same dimension." );
83 const std::vector<Element>& xval = *ex.
getVector();
84 std::copy(xval.begin(),xval.end(),
std_vec_->begin());
90 std::invalid_argument,
91 "Error: Vectors must have the same dimension." );
94 const std::vector<Element>& xval = *ex.
getVector();
96 for (
uint i=0; i<dim; i++) {
97 (*std_vec_)[i] += xval[i];
104 std::invalid_argument,
105 "Error: Vectors must have the same dimension." );
108 const std::vector<Element>& xval = *ex.
getVector();
110 for (
uint i=0; i<dim; i++) {
111 (*std_vec_)[i] += alpha*xval[i];
117 for (
uint i=0; i<dim; i++) {
118 (*std_vec_)[i] *= alpha;
125 std::invalid_argument,
126 "Error: Vectors must have the same dimension." );
129 const std::vector<Element>& xval = *ex.
getVector();
132 for (
uint i=0; i<dim; i++) {
133 val += (*std_vec_)[i]*xval[i];
140 val = std::sqrt(
dot(*
this) );
144 virtual Ptr<Vector<Real> >
clone()
const {
145 return makePtr<StdVector>( makePtr<std::vector<Element>>(
std_vec_->size(),
static_cast<Element
>(0)));
156 Ptr<Vector<Real> >
basis(
const int i )
const {
158 ROL_TEST_FOR_EXCEPTION( i >=
dimension() || i<0,
159 std::invalid_argument,
160 "Error: Basis index must be between 0 and vector dimension." );
162 Ptr<Vector<Real>> e =
clone();
163 (*staticPtrCast<StdVector>(e)->
getVector())[i] = 1.0;
168 return static_cast<int>(
std_vec_->size());
171 void applyUnary(
const Elementwise::UnaryFunction<Real> &f ) {
173 for(
uint i=0; i<dim; ++i) {
174 (*std_vec_)[i] = f.apply((*
std_vec_)[i]);
181 std::invalid_argument,
182 "Error: Vectors must have the same dimension." );
185 const std::vector<Element>& xval = *ex.
getVector();
187 for (
uint i=0; i<dim; i++) {
188 (*std_vec_)[i] = f.apply((*
std_vec_)[i],xval[i]);
193 Real
reduce(
const Elementwise::ReductionOp<Real> &r )
const {
194 Real result = r.initialValue();
196 for(
uint i=0; i<dim; ++i) {
207 void randomize(
const Real l = 0.0,
const Real u = 1.0 ) {
212 for (
uint i=0; i<dim; ++i) {
213 x =
static_cast<Real
>(rand())/static_cast<Real>(RAND_MAX);
214 (*std_vec_)[i] = a*x + b;
218 virtual void print( std::ostream &outStream )
const {
220 for(
uint i=0; i<dim; ++i) {
221 outStream << (*std_vec_)[i] <<
" ";
223 outStream << std::endl;
typename PV< Real >::size_type size_type
void axpy(const Real alpha, const Vector< Real > &x)
Compute where .
void scale(const Real alpha)
Compute where .
StdVector(const int dim, const Element val=0.0)
virtual int dimension() const
Return dimension of the vector space.
Ptr< Vector< Real > > basis(const int i) const
Return i-th basis vector.
void randomize(const Real l=0.0, const Real u=1.0)
Set vector to be uniform random between [l,u].
Ptr< const std::vector< Element > > getVector() const
StdVector(const Ptr< std::vector< Element > > &std_vec)
void applyBinary(const Elementwise::BinaryFunction< Real > &f, const Vector< Real > &x)
Defines the linear algebra or vector space interface.
virtual void print(std::ostream &outStream) const
Real norm() const
Returns where .
Provides the ROL::Vector interface for scalar values, to be used, for example, with scalar constraint...
void plus(const Vector< Real > &x)
Compute , where .
virtual Ptr< Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
Real reduce(const Elementwise::ReductionOp< Real > &r) const
int dimension() const
Return dimension of the vector space.
std::vector< Real >::size_type uint
void set(const Vector< Real > &x)
Set where .
void setScalar(const Real C)
Set where .
Ptr< std::vector< Element > > std_vec_
void applyUnary(const Elementwise::UnaryFunction< Real > &f)
virtual Real dot(const Vector< Real > &x) const
Compute where .
Ptr< std::vector< Element > > getVector()