diff options
Diffstat (limited to 'thirdparty/linux/include/coin/IpVector.hpp')
| -rw-r--r-- | thirdparty/linux/include/coin/IpVector.hpp | 774 |
1 files changed, 774 insertions, 0 deletions
diff --git a/thirdparty/linux/include/coin/IpVector.hpp b/thirdparty/linux/include/coin/IpVector.hpp new file mode 100644 index 0000000..a903558 --- /dev/null +++ b/thirdparty/linux/include/coin/IpVector.hpp @@ -0,0 +1,774 @@ +// Copyright (C) 2004, 2008 International Business Machines and others. +// All Rights Reserved. +// This code is published under the Eclipse Public License. +// +// $Id: IpVector.hpp 2472 2014-04-05 17:47:20Z stefan $ +// +// Authors: Carl Laird, Andreas Waechter IBM 2004-08-13 + +#ifndef __IPVECTOR_HPP__ +#define __IPVECTOR_HPP__ + +#include "IpTypes.hpp" +#include "IpTaggedObject.hpp" +#include "IpCachedResults.hpp" +#include "IpSmartPtr.hpp" +#include "IpJournalist.hpp" +#include "IpException.hpp" + +#include <vector> + +namespace Ipopt +{ + /** Exception that can be used to flag unimplemented linear algebra + * methods */ + DECLARE_STD_EXCEPTION(UNIMPLEMENTED_LINALG_METHOD_CALLED); + + /* forward declarations */ + class VectorSpace; + + /** Vector Base Class. + * This is the base class for all derived vector types. Those vectors + * are meant to store entities like iterates, Lagrangian multipliers, + * constraint values etc. The implementation of a vector type depends + * on the computational environment (e.g. just a double array on a shared + * memory machine, or distributed double arrays for a distributed + * memory machine.) + * + * Deriving from Vector: This class inherits from tagged object to + * implement an advanced caching scheme. Because of this, the + * TaggedObject method ObjectChanged() must be called each time the + * Vector changes. If you overload the XXXX_Impl protected methods, + * this taken care of (along with caching if possible) for you. If + * you have additional methods in your derived class that change the + * underlying data (vector values), you MUST remember to call + * ObjectChanged() AFTER making the change! + */ + class Vector : public TaggedObject + { + public: + /** @name Constructor/Destructor */ + //@{ + /** Constructor. It has to be given a pointer to the + * corresponding VectorSpace. + */ + inline + Vector(const VectorSpace* owner_space); + + /** Destructor */ + inline + virtual ~Vector(); + //@} + + /** Create new Vector of the same type with uninitialized data */ + inline + Vector* MakeNew() const; + + /** Create new Vector of the same type and copy the data over */ + inline + Vector* MakeNewCopy() const; + + /**@name Standard BLAS-1 Operations + * (derived classes do NOT overload these + * methods, instead, overload the + * protected versions of these methods). */ + //@{ + /** Copy the data of the vector x into this vector (DCOPY). */ + inline + void Copy(const Vector& x); + + /** Scales the vector by scalar alpha (DSCAL) */ + void Scal(Number alpha); + + /** Add the multiple alpha of vector x to this vector (DAXPY) */ + inline + void Axpy(Number alpha, const Vector &x); + + /** Computes inner product of vector x with this (DDOT) */ + inline + Number Dot(const Vector &x) const; + + /** Computes the 2-norm of this vector (DNRM2) */ + inline + Number Nrm2() const; + + /** Computes the 1-norm of this vector (DASUM) */ + inline + Number Asum() const; + + /** Computes the max-norm of this vector (based on IDAMAX) */ + inline + Number Amax() const; + //@} + + /** @name Additional (Non-BLAS) Vector Methods + * (derived classes do NOT overload these + * methods, instead, overload the + * protected versions of these methods). */ + //@{ + /** Set each element in the vector to the scalar alpha. */ + inline + void Set(Number alpha); + + /** Element-wise division \f$y_i \gets y_i/x_i\f$*/ + inline + void ElementWiseDivide(const Vector& x); + + /** Element-wise multiplication \f$y_i \gets y_i*x_i\f$ */ + inline + void ElementWiseMultiply(const Vector& x); + + /** Element-wise max against entries in x */ + inline + void ElementWiseMax(const Vector& x); + + /** Element-wise min against entries in x */ + inline + void ElementWiseMin(const Vector& x); + + /** Reciprocates the entries in the vector */ + inline + void ElementWiseReciprocal(); + + /** Absolute values of the entries in the vector */ + inline + void ElementWiseAbs(); + + /** Element-wise square root of the entries in the vector */ + inline + void ElementWiseSqrt(); + + /** Replaces the vector values with their sgn values + ( -1 if x_i < 0, 0 if x_i == 0, and 1 if x_i > 0) + */ + inline + void ElementWiseSgn(); + + /** Add scalar to every vector component */ + inline + void AddScalar(Number scalar); + + /** Returns the maximum value in the vector */ + inline + Number Max() const; + + /** Returns the minimum value in the vector */ + inline + Number Min() const; + + /** Returns the sum of the vector entries */ + inline + Number Sum() const; + + /** Returns the sum of the logs of each vector entry */ + inline + Number SumLogs() const; + //@} + + /** @name Methods for specialized operations. A prototype + * implementation is provided, but for efficient implementation + * those should be specially implemented. + */ + //@{ + /** Add one vector, y = a * v1 + c * y. This is automatically + * reduced to call AddTwoVectors. */ + inline + void AddOneVector(Number a, const Vector& v1, Number c); + + /** Add two vectors, y = a * v1 + b * v2 + c * y. Here, this + * vector is y */ + inline void AddTwoVectors(Number a, const Vector& v1, + Number b, const Vector& v2, Number c); + /** Fraction to the boundary parameter. Computes \f$\alpha = + * \max\{\bar\alpha\in(0,1] : x + \bar\alpha \Delta \geq (1-\tau)x\}\f$ + */ + inline + Number FracToBound(const Vector& delta, Number tau) const; + /** Add the quotient of two vectors, y = a * z/s + c * y. */ + inline + void AddVectorQuotient(Number a, const Vector& z, const Vector& s, + Number c); + //@} + + /** Method for determining if all stored numbers are valid (i.e., + * no Inf or Nan). */ + inline + bool HasValidNumbers() const; + + /** @name Accessor methods */ + //@{ + /** Dimension of the Vector */ + inline + Index Dim() const; + + /** Return the owner VectorSpace*/ + inline + SmartPtr<const VectorSpace> OwnerSpace() const; + //@} + + /** @name Output methods + * (derived classes do NOT overload these + * methods, instead, overload the + * protected versions of these methods). */ + //@{ + /** Print the entire vector */ + void Print(SmartPtr<const Journalist> jnlst, + EJournalLevel level, + EJournalCategory category, + const std::string& name, + Index indent=0, + const std::string& prefix="") const; + void Print(const Journalist& jnlst, + EJournalLevel level, + EJournalCategory category, + const std::string& name, + Index indent=0, + const std::string& prefix="") const; + //@} + + protected: + /** @name implementation methods (derived classes MUST + * overload these pure virtual protected methods.) + */ + //@{ + /** Copy the data of the vector x into this vector (DCOPY). */ + virtual void CopyImpl(const Vector& x)=0; + + /** Scales the vector by scalar alpha (DSCAL) */ + virtual void ScalImpl(Number alpha)=0; + + /** Add the multiple alpha of vector x to this vector (DAXPY) */ + virtual void AxpyImpl(Number alpha, const Vector &x)=0; + + /** Computes inner product of vector x with this (DDOT) */ + virtual Number DotImpl(const Vector &x) const =0; + + /** Computes the 2-norm of this vector (DNRM2) */ + virtual Number Nrm2Impl() const =0; + + /** Computes the 1-norm of this vector (DASUM) */ + virtual Number AsumImpl() const =0; + + /** Computes the max-norm of this vector (based on IDAMAX) */ + virtual Number AmaxImpl() const =0; + + /** Set each element in the vector to the scalar alpha. */ + virtual void SetImpl(Number alpha)=0; + + /** Element-wise division \f$y_i \gets y_i/x_i\f$*/ + virtual void ElementWiseDivideImpl(const Vector& x)=0; + + /** Element-wise multiplication \f$y_i \gets y_i*x_i\f$ */ + virtual void ElementWiseMultiplyImpl(const Vector& x)=0; + + /** Element-wise max against entries in x */ + virtual void ElementWiseMaxImpl(const Vector& x)=0; + + /** Element-wise min against entries in x */ + virtual void ElementWiseMinImpl(const Vector& x)=0; + + /** Reciprocates the elements of the vector */ + virtual void ElementWiseReciprocalImpl()=0; + + /** Take elementwise absolute values of the elements of the vector */ + virtual void ElementWiseAbsImpl()=0; + + /** Take elementwise square-root of the elements of the vector */ + virtual void ElementWiseSqrtImpl()=0; + + /** Replaces entries with sgn of the entry */ + virtual void ElementWiseSgnImpl()=0; + + /** Add scalar to every component of vector */ + virtual void AddScalarImpl(Number scalar)=0; + + /** Max value in the vector */ + virtual Number MaxImpl() const=0; + + /** Min number in the vector */ + virtual Number MinImpl() const=0; + + /** Sum of entries in the vector */ + virtual Number SumImpl() const=0; + + /** Sum of logs of entries in the vector */ + virtual Number SumLogsImpl() const=0; + + /** Add two vectors (a * v1 + b * v2). Result is stored in this + vector. */ + virtual void AddTwoVectorsImpl(Number a, const Vector& v1, + Number b, const Vector& v2, Number c); + + /** Fraction to boundary parameter. */ + virtual Number FracToBoundImpl(const Vector& delta, Number tau) const; + + /** Add the quotient of two vectors */ + virtual void AddVectorQuotientImpl(Number a, const Vector& z, + const Vector& s, Number c); + + /** Method for determining if all stored numbers are valid (i.e., + * no Inf or Nan). A default implementation using Asum is + * provided. */ + virtual bool HasValidNumbersImpl() const; + + /** Print the entire vector */ + virtual void PrintImpl(const Journalist& jnlst, + EJournalLevel level, + EJournalCategory category, + const std::string& name, + Index indent, + const std::string& prefix) const =0; + //@} + + private: + /**@name Default Compiler Generated Methods + * (Hidden to avoid implicit creation/calling). + * These methods are not implemented and + * we do not want the compiler to implement + * them for us, so we declare them private + * and do not define them. This ensures that + * they will not be implicitly created/called. */ + //@{ + /** Default constructor */ + Vector(); + + /** Copy constructor */ + Vector(const Vector&); + + /** Overloaded Equals Operator */ + Vector& operator=(const Vector&); + //@} + + /** Vector Space */ + const SmartPtr<const VectorSpace> owner_space_; + + /**@name CachedResults data members */ + //@{ + /** Cache for dot products */ + mutable CachedResults<Number> dot_cache_; + + mutable TaggedObject::Tag nrm2_cache_tag_; + mutable Number cached_nrm2_; + + mutable TaggedObject::Tag asum_cache_tag_; + mutable Number cached_asum_; + + mutable TaggedObject::Tag amax_cache_tag_; + mutable Number cached_amax_; + + mutable TaggedObject::Tag max_cache_tag_; + mutable Number cached_max_; + + mutable TaggedObject::Tag min_cache_tag_; + mutable Number cached_min_; + + mutable TaggedObject::Tag sum_cache_tag_; + mutable Number cached_sum_; + + mutable TaggedObject::Tag sumlogs_cache_tag_; + mutable Number cached_sumlogs_; + + mutable TaggedObject::Tag valid_cache_tag_; + mutable bool cached_valid_; + + // AW: I removed this cache since it gets in the way for the + // quality function search + // /** Cache for FracToBound */ + // mutable CachedResults<Number> frac_to_bound_cache_; + //@} + + }; + + /** VectorSpace base class, corresponding to the Vector base class. + * For each Vector implementation, a corresponding VectorSpace has + * to be implemented. A VectorSpace is able to create new Vectors + * of a specific type. The VectorSpace should also store + * information that is common to all Vectors of that type. For + * example, the dimension of a Vector is stored in the VectorSpace + * base class. + */ + class VectorSpace : public ReferencedObject + { + public: + /** @name Constructors/Destructors */ + //@{ + /** Constructor, given the dimension of all vectors generated by + * this VectorSpace. + */ + VectorSpace(Index dim); + + /** Destructor */ + virtual ~VectorSpace() + {} + //@} + + /** Pure virtual method for creating a new Vector of the + * corresponding type. + */ + virtual Vector* MakeNew() const=0; + + /** Accessor function for the dimension of the vectors of this type.*/ + Index Dim() const + { + return dim_; + } + + private: + /**@name Default Compiler Generated Methods + * (Hidden to avoid implicit creation/calling). + * These methods are not implemented and + * we do not want the compiler to implement + * them for us, so we declare them private + * and do not define them. This ensures that + * they will not be implicitly created/called. */ + //@{ + /** default constructor */ + VectorSpace(); + + /** Copy constructor */ + VectorSpace(const VectorSpace&); + + /** Overloaded Equals Operator */ + VectorSpace& operator=(const VectorSpace&); + //@} + + /** Dimension of the vectors in this vector space. */ + const Index dim_; + }; + + /* inline methods */ + inline + Vector::~Vector() + {} + + inline + Vector::Vector(const VectorSpace* owner_space) + : + TaggedObject(), + owner_space_(owner_space), + dot_cache_(10), + nrm2_cache_tag_(0), + asum_cache_tag_(0), + amax_cache_tag_(0), + max_cache_tag_(0), + min_cache_tag_(0), + sum_cache_tag_(0), + sumlogs_cache_tag_(0), + cached_valid_(0) + { + DBG_ASSERT(IsValid(owner_space_)); + } + + inline + Vector* Vector::MakeNew() const + { + return owner_space_->MakeNew(); + } + + inline + Vector* Vector::MakeNewCopy() const + { + // ToDo: We can probably copy also the cached values for Norms etc here + Vector* copy = MakeNew(); + copy->Copy(*this); + return copy; + } + + inline + void Vector::Copy(const Vector& x) + { + CopyImpl(x); + ObjectChanged(); + // Also copy any cached scalar values from the original vector + // ToDo: Check if that is too much overhead + TaggedObject::Tag x_tag = x.GetTag(); + if (x_tag == x.nrm2_cache_tag_) { + nrm2_cache_tag_ = GetTag(); + cached_nrm2_ = x.cached_nrm2_; + } + if (x_tag == x.asum_cache_tag_) { + asum_cache_tag_ = GetTag(); + cached_asum_ = x.cached_asum_; + } + if (x_tag == x.amax_cache_tag_) { + amax_cache_tag_ = GetTag(); + cached_amax_ = x.cached_amax_; + } + if (x_tag == x.max_cache_tag_) { + max_cache_tag_ = GetTag(); + cached_max_ = x.cached_max_; + } + if (x_tag == x.min_cache_tag_) { + min_cache_tag_ = GetTag(); + cached_min_ = x.cached_min_; + } + if (x_tag == x.sum_cache_tag_) { + sum_cache_tag_ = GetTag(); + cached_sum_ = x.cached_sum_; + } + if (x_tag == x.sumlogs_cache_tag_) { + sumlogs_cache_tag_ = GetTag(); + cached_sumlogs_ = x.cached_sumlogs_; + } + } + + inline + void Vector::Axpy(Number alpha, const Vector &x) + { + AxpyImpl(alpha, x); + ObjectChanged(); + } + + inline + Number Vector::Dot(const Vector &x) const + { + // The current implementation of the caching doesn't allow to have + // a dependency of something with itself. Therefore, we use the + // Nrm2 method if the dot product is to be taken with the vector + // itself. Might be more efficient anyway. + if (this==&x) { + Number nrm2 = Nrm2(); + return nrm2*nrm2; + } + Number retValue; + if (!dot_cache_.GetCachedResult2Dep(retValue, this, &x)) { + retValue = DotImpl(x); + dot_cache_.AddCachedResult2Dep(retValue, this, &x); + } + return retValue; + } + + inline + Number Vector::Nrm2() const + { + if (nrm2_cache_tag_ != GetTag()) { + cached_nrm2_ = Nrm2Impl(); + nrm2_cache_tag_ = GetTag(); + } + return cached_nrm2_; + } + + inline + Number Vector::Asum() const + { + if (asum_cache_tag_ != GetTag()) { + cached_asum_ = AsumImpl(); + asum_cache_tag_ = GetTag(); + } + return cached_asum_; + } + + inline + Number Vector::Amax() const + { + if (amax_cache_tag_ != GetTag()) { + cached_amax_ = AmaxImpl(); + amax_cache_tag_ = GetTag(); + } + return cached_amax_; + } + + inline + Number Vector::Sum() const + { + if (sum_cache_tag_ != GetTag()) { + cached_sum_ = SumImpl(); + sum_cache_tag_ = GetTag(); + } + return cached_sum_; + } + + inline + Number Vector::SumLogs() const + { + if (sumlogs_cache_tag_ != GetTag()) { + cached_sumlogs_ = SumLogsImpl(); + sumlogs_cache_tag_ = GetTag(); + } + return cached_sumlogs_; + } + + inline + void Vector::ElementWiseSgn() + { + ElementWiseSgnImpl(); + ObjectChanged(); + } + + inline + void Vector::Set(Number alpha) + { + // Could initialize caches here + SetImpl(alpha); + ObjectChanged(); + } + + inline + void Vector::ElementWiseDivide(const Vector& x) + { + ElementWiseDivideImpl(x); + ObjectChanged(); + } + + inline + void Vector::ElementWiseMultiply(const Vector& x) + { + ElementWiseMultiplyImpl(x); + ObjectChanged(); + } + + inline + void Vector::ElementWiseReciprocal() + { + ElementWiseReciprocalImpl(); + ObjectChanged(); + } + + inline + void Vector::ElementWiseMax(const Vector& x) + { + // Could initialize some caches here + ElementWiseMaxImpl(x); + ObjectChanged(); + } + + inline + void Vector::ElementWiseMin(const Vector& x) + { + // Could initialize some caches here + ElementWiseMinImpl(x); + ObjectChanged(); + } + + inline + void Vector::ElementWiseAbs() + { + // Could initialize some caches here + ElementWiseAbsImpl(); + ObjectChanged(); + } + + inline + void Vector::ElementWiseSqrt() + { + ElementWiseSqrtImpl(); + ObjectChanged(); + } + + inline + void Vector::AddScalar(Number scalar) + { + // Could initialize some caches here + AddScalarImpl(scalar); + ObjectChanged(); + } + + inline + Number Vector::Max() const + { + if (max_cache_tag_ != GetTag()) { + cached_max_ = MaxImpl(); + max_cache_tag_ = GetTag(); + } + return cached_max_; + } + + inline + Number Vector::Min() const + { + if (min_cache_tag_ != GetTag()) { + cached_min_ = MinImpl(); + min_cache_tag_ = GetTag(); + } + return cached_min_; + } + + inline + void Vector::AddOneVector(Number a, const Vector& v1, Number c) + { + AddTwoVectors(a, v1, 0., v1, c); + } + + inline + void Vector::AddTwoVectors(Number a, const Vector& v1, + Number b, const Vector& v2, Number c) + { + AddTwoVectorsImpl(a, v1, b, v2, c); + ObjectChanged(); + } + + inline + Number Vector::FracToBound(const Vector& delta, Number tau) const + { + /* AW: I avoid the caching here, since it leads to overhead in the + quality function search. Caches for this are in + CalculatedQuantities. + Number retValue; + std::vector<const TaggedObject*> tdeps(1); + tdeps[0] = δ + std::vector<Number> sdeps(1); + sdeps[0] = tau; + if (!frac_to_bound_cache_.GetCachedResult(retValue, tdeps, sdeps)) { + retValue = FracToBoundImpl(delta, tau); + frac_to_bound_cache_.AddCachedResult(retValue, tdeps, sdeps); + } + return retValue; + */ + return FracToBoundImpl(delta, tau); + } + + inline + void Vector::AddVectorQuotient(Number a, const Vector& z, + const Vector& s, Number c) + { + AddVectorQuotientImpl(a, z, s, c); + ObjectChanged(); + } + + inline + bool Vector::HasValidNumbers() const + { + if (valid_cache_tag_ != GetTag()) { + cached_valid_ = HasValidNumbersImpl(); + valid_cache_tag_ = GetTag(); + } + return cached_valid_; + } + + inline + Index Vector::Dim() const + { + return owner_space_->Dim(); + } + + inline + SmartPtr<const VectorSpace> Vector::OwnerSpace() const + { + return owner_space_; + } + + inline + VectorSpace::VectorSpace(Index dim) + : + dim_(dim) + {} + +} // namespace Ipopt + +// Macro definitions for debugging vectors +#if COIN_IPOPT_VERBOSITY == 0 +# define DBG_PRINT_VECTOR(__verbose_level, __vec_name, __vec) +#else +# define DBG_PRINT_VECTOR(__verbose_level, __vec_name, __vec) \ + if (dbg_jrnl.Verbosity() >= (__verbose_level)) { \ + if (dbg_jrnl.Jnlst()!=NULL) { \ + (__vec).Print(dbg_jrnl.Jnlst(), \ + J_ERROR, J_DBG, \ + __vec_name, \ + dbg_jrnl.IndentationLevel()*2, \ + "# "); \ + } \ + } +#endif //if COIN_IPOPT_VERBOSITY == 0 + +#endif |