quadratic-program.hh

Go to the documentation of this file.

// Copyright (c) 2018, Joseph Mirabel
// Authors: Joseph Mirabel (joseph.mirabel@laas.fr)
//
 
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// 1. Redistributions of source code must retain the above copyright
//    notice, this list of conditions and the following disclaimer.
//
// 2. Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
// HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
// SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
// LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
// THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
// OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
// DAMAGE.
 
#ifndef HPP_CORE_PATH_OPTIMIZATION_QUADRATIC_PROGRAM_HH
#define HPP_CORE_PATH_OPTIMIZATION_QUADRATIC_PROGRAM_HH
 
#include <hpp/core/fwd.hh>
#include <hpp/core/path-optimization/linear-constraint.hh>
 
namespace hpp {
namespace core {
namespace pathOptimization {
struct QuadraticProgram {
  typedef Eigen::JacobiSVD<matrix_t> Decomposition_t;
  typedef Eigen::LLT<matrix_t, Eigen::Lower> LLT_t;

  QuadraticProgram(size_type inputSize)
      : H(inputSize, inputSize),
        b(inputSize),
        dec(inputSize, inputSize, Eigen::ComputeThinU | Eigen::ComputeThinV),
        xStar(inputSize),
        accuracy_(1e-4) {
    H.setZero();
    b.setZero();
    bIsZero = true;
  }

  QuadraticProgram(const QuadraticProgram& QP, const LinearConstraint& lc)
      : H(lc.PK.cols(), lc.PK.cols()),
        b(lc.PK.cols()),
        bIsZero(false),
        dec(lc.PK.cols(), lc.PK.cols(),
            Eigen::ComputeThinU | Eigen::ComputeThinV),
        xStar(lc.PK.cols()),
        accuracy_(1e-4) {
    QP.reduced(lc, *this);
  }
 
  QuadraticProgram(const QuadraticProgram& QP)
      : H(QP.H),
        b(QP.b),
        bIsZero(QP.bIsZero),
        dec(QP.dec),
        xStar(QP.xStar),
        accuracy_(QP.accuracy_) {}
 
  ~QuadraticProgram();

  void accuracy(value_type acc) { accuracy_ = acc; }
  value_type accuracy() const { return accuracy_; }
  void addRows(const std::size_t& nbRows) {
    H.conservativeResize(H.rows() + nbRows, H.cols());
    b.conservativeResize(b.rows() + nbRows, b.cols());
 
    H.bottomRows(nbRows).setZero();
  }

 
  /*/ Compute the problem
   *  \f{eqnarray*}{
   *  \min & \frac{1}{2} * x^T H x + b^T x \\
   *      lc.J * x = lc.b
   *  \f}
   **/
  void reduced(const LinearConstraint& lc, QuadraticProgram& QPr) const {
    matrix_t H_PK(H * lc.PK);
    QPr.H.noalias() = lc.PK.transpose() * H_PK;
    QPr.b.noalias() = H_PK.transpose() * lc.xStar;
    if (!bIsZero) {
      QPr.b.noalias() += lc.PK.transpose() * b;
    }
    QPr.bIsZero = false;
  }
 
  void decompose();
 
  void solve() { xStar.noalias() = -dec.solve(b); }


 
  void computeLLT();

  double solve(const LinearConstraint& ce, const LinearConstraint& ci);


  matrix_t H;
  vector_t b;
  bool bIsZero;

  LLT_t llt;
  value_type trace;
  Eigen::VectorXi activeConstraint;
  int activeSetSize;

  Decomposition_t dec;
  vector_t xStar;
  value_type accuracy_;
};
}  // namespace pathOptimization
}  // namespace core
}  // namespace hpp
 
#endif  // HPP_CORE_PATH_OPTIMIZATION_QUADRATIC_PROGRAM_HH