IpVector.hpp

Go to the documentation of this file.

// Copyright (C) 2004, 2006 International Business Machines and others.
// All Rights Reserved.
// This code is published under the Common Public License.
//
// $Id: IpVector.hpp 759 2006-07-07 03:07:08Z andreasw $
//
// 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 <vector>

namespace Ipopt
{
  /* forward declarations */
  class VectorSpace;

  class Vector : public TaggedObject
  {
  public:
    Vector(const VectorSpace* owner_space);

    virtual ~Vector();

    Vector* MakeNew() const;

    Vector* MakeNewCopy() const;

    void Copy(const Vector& x);

    void Scal(Number alpha);

    void Axpy(Number alpha, const Vector &x);

    Number Dot(const Vector &x) const;

    Number Nrm2() const;

    Number Asum() const;

    Number Amax() const;

    void Set(Number alpha);

    void ElementWiseDivide(const Vector& x);

    void ElementWiseMultiply(const Vector& x);

    void ElementWiseMax(const Vector& x);

    void ElementWiseMin(const Vector& x);

    void ElementWiseReciprocal();

    void ElementWiseAbs();

    void ElementWiseSqrt();

    void ElementWiseSgn();

    void AddScalar(Number scalar);

    Number Max() const;

    Number Min() const;

    Number Sum() const;

    Number SumLogs() const;

    void AddOneVector(Number a, const Vector& v1, Number c);

    void AddTwoVectors(Number a, const Vector& v1,
                       Number b, const Vector& v2, Number c);
    Number FracToBound(const Vector& delta, Number tau) const;
    void AddVectorQuotient(Number a, const Vector& z, const Vector& s,
                           Number c);

    bool HasValidNumbers() const;

    Index Dim() const;

    SmartPtr<const VectorSpace> OwnerSpace() const;

    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:
    virtual void CopyImpl(const Vector& x)=0;

    virtual void ScalImpl(Number alpha)=0;

    virtual void AxpyImpl(Number alpha, const Vector &x)=0;

    virtual Number DotImpl(const Vector &x) const =0;

    virtual Number Nrm2Impl() const =0;

    virtual Number AsumImpl() const =0;

    virtual Number AmaxImpl() const =0;

    virtual void SetImpl(Number alpha)=0;

    virtual void ElementWiseDivideImpl(const Vector& x)=0;

    virtual void ElementWiseMultiplyImpl(const Vector& x)=0;

    virtual void ElementWiseMaxImpl(const Vector& x)=0;

    virtual void ElementWiseMinImpl(const Vector& x)=0;

    virtual void ElementWiseReciprocalImpl()=0;

    virtual void ElementWiseAbsImpl()=0;

    virtual void ElementWiseSqrtImpl()=0;

    virtual void ElementWiseSgnImpl()=0;

    virtual void AddScalarImpl(Number scalar)=0;

    virtual Number MaxImpl() const=0;

    virtual Number MinImpl() const=0;

    virtual Number SumImpl() const=0;

    virtual Number SumLogsImpl() const=0;

    virtual void AddTwoVectorsImpl(Number a, const Vector& v1,
                                   Number b, const Vector& v2, Number c);

    virtual Number FracToBoundImpl(const Vector& delta, Number tau) const;

    virtual void AddVectorQuotientImpl(Number a, const Vector& z,
                                       const Vector& s, Number c);

    virtual bool HasValidNumbersImpl() const;

    virtual void PrintImpl(const Journalist& jnlst,
                           EJournalLevel level,
                           EJournalCategory category,
                           const std::string& name,
                           Index indent,
                           const std::string& prefix) const =0;

  private:
    Vector();

    Vector(const Vector&);

    Vector& operator=(const Vector&);

    const SmartPtr<const VectorSpace> owner_space_;

    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_;

  };

  class VectorSpace : public ReferencedObject
  {
  public:
    VectorSpace(Index dim);

    virtual ~VectorSpace()
    {}

    virtual Vector* MakeNew() const=0;

    Index Dim() const
    {
      return dim_;
    }

  private:
    VectorSpace();

    VectorSpace(const VectorSpace&);

    VectorSpace& operator=(const VectorSpace&);

    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)
  {
    AddTwoVectorsImpl(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] = &delta;
    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
#ifndef IP_DEBUG
# 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 // IP_DEBUG

#endif