You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
121 lines
4.3 KiB
121 lines
4.3 KiB
1 year ago
|
// This file is part of libigl, a simple c++ geometry processing library.
|
||
|
|
//
|
||
|
|
// Copyright (C) 2013 Alec Jacobson <alecjacobson@gmail.com>
|
||
|
|
//
|
||
|
|
// This Source Code Form is subject to the terms of the Mozilla Public License
|
||
|
|
// v. 2.0. If a copy of the MPL was not distributed with this file, You can
|
||
|
|
// obtain one at http://mozilla.org/MPL/2.0/.
|
||
|
|
#ifndef IGL_ACTIVE_SET_H
|
||
|
|
#define IGL_ACTIVE_SET_H
|
||
|
|
|
||
|
|
#include "igl_inline.h"
|
||
|
|
#include "SolverStatus.h"
|
||
|
|
#include <Eigen/Core>
|
||
|
|
#include <Eigen/Sparse>
|
||
|
|
|
||
|
|
namespace igl
|
||
|
|
{
|
||
|
|
struct active_set_params;
|
||
|
|
///
|
||
|
|
/// Minimize convex quadratic energy subject to linear inequality constraints
|
||
|
|
///
|
||
|
|
/// min ½ Zᵀ A Z + Zᵀ B + constant
|
||
|
|
/// Z
|
||
|
|
/// subject to
|
||
|
|
/// Aeq Z = Beq
|
||
|
|
/// Aieq Z <= Bieq
|
||
|
|
/// lx <= Z <= ux
|
||
|
|
/// Z(known) = Y
|
||
|
|
///
|
||
|
|
/// that Z(known) = Y, optionally also subject to the constraints Aeq*Z = Beq,
|
||
|
|
/// and further optionally subject to the linear inequality constraints that
|
||
|
|
/// Aieq*Z <= Bieq and constant inequality constraints lx <= x <= ux
|
||
|
|
///
|
||
|
|
/// @param[in] A n by n matrix of quadratic coefficients
|
||
|
|
/// @param[in] B n by 1 column of linear coefficients
|
||
|
|
/// @param[in] known list of indices to known rows in Z
|
||
|
|
/// @param[in] Y list of fixed values corresponding to known rows in Z
|
||
|
|
/// @param[in] Aeq meq by n list of linear equality constraint coefficients
|
||
|
|
/// @param[in] Beq meq by 1 list of linear equality constraint constant values
|
||
|
|
/// @param[in] Aieq mieq by n list of linear inequality constraint coefficients
|
||
|
|
/// @param[in] Bieq mieq by 1 list of linear inequality constraint constant values
|
||
|
|
/// @param[in] lx n by 1 list of lower bounds [] implies -Inf
|
||
|
|
/// @param[in] ux n by 1 list of upper bounds [] implies Inf
|
||
|
|
/// @param[in] params struct of additional parameters (see below)
|
||
|
|
/// @param[in,out] Z if not empty, is taken to be an n by 1 list of initial guess values. Set to solution on output.
|
||
|
|
/// @return true on success, false on error
|
||
|
|
///
|
||
|
|
/// \note Benchmark: For a harmonic solve on a mesh with 325K facets, matlab 2.2
|
||
|
|
/// secs, igl/min_quad_with_fixed.h 7.1 secs
|
||
|
|
///
|
||
|
|
/// \pre rows of [Aeq;Aieq] **must** be linearly independent. Should be
|
||
|
|
/// using QR decomposition otherwise:
|
||
|
|
/// https://v8doc.sas.com/sashtml/ormp/chap5/sect32.htm
|
||
|
|
///
|
||
|
|
/// \warning This solver is fairly experimental. It works reasonably well for
|
||
|
|
/// bbw problems but doesn't generalize well to other problems. NASOQ and
|
||
|
|
/// OSQP are better general purpose solvers.
|
||
|
|
template <
|
||
|
|
typename AT,
|
||
|
|
typename DerivedB,
|
||
|
|
typename Derivedknown,
|
||
|
|
typename DerivedY,
|
||
|
|
typename AeqT,
|
||
|
|
typename DerivedBeq,
|
||
|
|
typename AieqT,
|
||
|
|
typename DerivedBieq,
|
||
|
|
typename Derivedlx,
|
||
|
|
typename Derivedux,
|
||
|
|
typename DerivedZ
|
||
|
|
>
|
||
|
|
IGL_INLINE igl::SolverStatus active_set(
|
||
|
|
const Eigen::SparseMatrix<AT>& A,
|
||
|
|
const Eigen::PlainObjectBase<DerivedB> & B,
|
||
|
|
const Eigen::PlainObjectBase<Derivedknown> & known,
|
||
|
|
const Eigen::PlainObjectBase<DerivedY> & Y,
|
||
|
|
const Eigen::SparseMatrix<AeqT>& Aeq,
|
||
|
|
const Eigen::PlainObjectBase<DerivedBeq> & Beq,
|
||
|
|
const Eigen::SparseMatrix<AieqT>& Aieq,
|
||
|
|
const Eigen::PlainObjectBase<DerivedBieq> & Bieq,
|
||
|
|
const Eigen::PlainObjectBase<Derivedlx> & lx,
|
||
|
|
const Eigen::PlainObjectBase<Derivedux> & ux,
|
||
|
|
const igl::active_set_params & params,
|
||
|
|
Eigen::PlainObjectBase<DerivedZ> & Z
|
||
|
|
);
|
||
|
|
};
|
||
|
|
|
||
|
|
#include "EPS.h"
|
||
|
|
/// Input parameters controling active_set
|
||
|
|
///
|
||
|
|
/// \fileinfo
|
||
|
|
struct igl::active_set_params
|
||
|
|
{
|
||
|
|
/// Auu_pd whether Auu is positive definite {false}
|
||
|
|
bool Auu_pd;
|
||
|
|
/// max_iter Maximum number of iterations (0 = Infinity, {100})
|
||
|
|
int max_iter;
|
||
|
|
/// inactive_threshold Threshold on Lagrange multiplier values to determine
|
||
|
|
/// whether to keep constraints active {EPS}
|
||
|
|
double inactive_threshold;
|
||
|
|
/// constraint_threshold Threshold on whether constraints are violated (0
|
||
|
|
/// is perfect) {EPS}
|
||
|
|
double constraint_threshold;
|
||
|
|
/// solution_diff_threshold Threshold on the squared norm of the difference
|
||
|
|
/// between two consecutive solutions {EPS}
|
||
|
|
double solution_diff_threshold;
|
||
|
|
/// @private
|
||
|
|
active_set_params():
|
||
|
|
Auu_pd(false),
|
||
|
|
max_iter(100),
|
||
|
|
inactive_threshold(igl::DOUBLE_EPS),
|
||
|
|
constraint_threshold(igl::DOUBLE_EPS),
|
||
|
|
solution_diff_threshold(igl::DOUBLE_EPS)
|
||
|
|
{};
|
||
|
|
};
|
||
|
|
|
||
|
|
#ifndef IGL_STATIC_LIBRARY
|
||
|
|
# include "active_set.cpp"
|
||
|
|
#endif
|
||
|
|
|
||
|
|
#endif
|