qp_solver.h

来自「CGAL is a collaborative effort of severa」· C头文件 代码 · 共 1,686 行 · 第 1/4 页

// Copyright (c) 1997-2001  ETH Zurich (Switzerland).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org); you may redistribute it under
// the terms of the Q Public License version 1.0.
// See the file LICENSE.QPL distributed with CGAL.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $Source: /CVSROOT/CGAL/Packages/_QP_solver/include/CGAL/_QP_solver/QP_solver.h,v $
// $Revision: 1.8 $ $Date: 2004/09/03 17:41:16 $
// $Name:  $
//
// Author(s)     : Sven Schoenherr <[e-mail]>
                                                                               
#ifndef CGAL_QP_SOLVER_H
#define CGAL_QP_SOLVER_H

// includes
#include <CGAL/Optimisation/basic.h>

#include <CGAL/_QP_solver/Basis_inverse.h>
//#include <CGAL/_QP_solver/Pricing_strategy_base.h>
//#include <CGAL/_QP_solver/Full_exact_pricing.h>

#include <CGAL/_QP_solver/Join_random_access_iterator.h>
#include <CGAL/_QP_solver/Sparse_vector_iterator.h>

#include <CGAL/_QP_solver/Compute_quotient.h>
#include <CGAL/_QP_solver/Access_by_index.h>

#include <CGAL/IO/Verbose_ostream.h>

#include <vector>
#include <functional>
#include <algorithm>
#include <iterator>
#include <numeric>


CGAL_BEGIN_NAMESPACE
                    

// Class declaration
// =================
template < class Rep_ >
class QP_solver;
                
template < class Rep_ >
class Pricing_strategy_base;
                
template < class Rep_ >
class Full_exact_pricing;
                

// Class interface
// ===============
template < class Rep_ >
class QP_solver {
  public:
    // self
    typedef  Rep_                       Rep;
    typedef  QP_solver<Rep>             Self;

    // types from the representation class
    typedef  typename Rep::NT           NT;
    typedef  typename Rep::ET           ET;

    typedef  typename Rep::A_iterator   A_iterator;
    typedef  typename Rep::B_iterator   B_iterator;
    typedef  typename Rep::C_iterator   C_iterator;
    typedef  typename Rep::D_iterator   D_iterator;

    typedef  typename Rep::Is_lp        Is_lp;

  private:
    // private types
    typedef  CGAL::Tag_true             Tag_true;
    typedef  CGAL::Tag_false            Tag_false;
    
    
    typedef  std::vector<int>           Indices;
    typedef  Indices::const_iterator    Index_iterator;
    
    typedef  std::vector<ET>            Values;
    typedef  typename Values::const_iterator
                                        Value_iterator;
    
    
    typedef  CGAL::Compute_quotient<ET> Compute_quotient;
    typedef  std::binder2nd< Compute_quotient >
                                        Make_quotient;
    
    
    typedef  CGAL::Access_by_index<Value_iterator,true,false>
                                        Access_value_checked;
    
    
    typedef  CGAL::Sparse_vector_iterator<NT>
                                        Artificial_column;
    
    typedef  std::vector<Artificial_column>
                                        A_artificial;
    
    
    typedef  std::vector<NT>            C_auxiliary;
    
    
    typedef  CGAL::Access_by_index< C_iterator,
                 false, false >         Access_c_B;
    typedef  CGAL::Join_random_access_iterator_1<
                 Index_iterator, Access_c_B >
                                        c_B_iterator;
    
    
    typedef  CGAL::Access_by_index<
                 typename std::iterator_traits<D_iterator>::value_type,
                 false, true >          Access_D_Bj;
    typedef  CGAL::Join_random_access_iterator_1<
                 Index_iterator, Access_D_Bj >
                                        D_Bj_iterator;
    
    

  public:
    
    // public types
    enum Status { UPDATE, INFEASIBLE, UNBOUNDED, OPTIMAL };
    
    
    typedef  Index_iterator             Basic_variable_index_iterator;
    typedef  Value_iterator             Basic_variable_numerator_iterator;
    typedef  CGAL::Join_random_access_iterator_1<
                 Basic_variable_numerator_iterator, Make_quotient >
                                        Basic_variable_value_iterator;
    
    
    typedef  CGAL::Join_random_access_iterator_1<
                 Index_iterator, Access_value_checked >
                                        Variable_numerator_iterator;
    typedef  CGAL::Join_random_access_iterator_1<
                 Variable_numerator_iterator, Make_quotient >
                                        Variable_value_iterator;
    
    
    typedef  Value_iterator             Lambda_numerator_iterator;
    typedef  CGAL::Join_random_access_iterator_1<
                 Lambda_numerator_iterator, Make_quotient >
                                        Lambda_value_iterator;
    
    
    typedef  CGAL::Pricing_strategy_base<Rep>
                                        Pricing_strategy;
    typedef  CGAL::Full_exact_pricing<Rep>
                                        Pricing_strategy_default;
    
    
    typedef  typename A_artificial::const_iterator
                                        A_artificial_iterator;
    
    
    typedef  typename C_auxiliary::const_iterator
                                        C_auxiliary_iterator;
    
    

    
    // creation
    QP_solver( int verbose = 0, std::ostream& stream = std::cout);
      /*
        : nt_0( 0), nt_1( 1), nt_minus_1( -nt_1),
          et_0( 0), et_1( 1), et_2( 2),
          vout1( verbose == 1, stream),
          vout2( verbose >= 2, stream),
          vout3( verbose == 3, stream),
          vout ( verbose >  0, stream),
          qp_n( 0), qp_m( 0), inv_M_B( vout3), m_phase( 0),
          d( inv_M_B.denominator()),
          
          strategyP( &strategy_default)
                                       
        {
            
            CGAL_optimisation_debug {
                vout2 << "======================================" << std::endl
                      << "The CGAL Solver for Quadratic Programs" << std::endl
                      << "======================================" << std::endl;
            }
             
        }
      */
    
    // set-up of QP
    void  set( int n, int m, int max_b,
               A_iterator A, B_iterator b, C_iterator c, D_iterator D);
    
    
    // initialization (of phase I)
    void  init( );
    
    
    // access to QP
    int  number_of_variables  ( ) const { return qp_n; }
    int  number_of_constraints( ) const { return qp_m; }
    
    A_iterator  a_begin( ) const { return qp_A;      }
    A_iterator  a_end  ( ) const { return qp_A+qp_n; }
    
    B_iterator  b_begin( ) const { return qp_b;      }
    B_iterator  b_end  ( ) const { return qp_b+qp_m; }
    
    C_iterator  c_begin( ) const { return qp_c;      }
    C_iterator  c_end  ( ) const { return qp_c+qp_n; }
    
    D_iterator  d_begin( ) const { return qp_D;      }
    D_iterator  d_end  ( ) const { return qp_D+qp_n; }
    
    
    // access to current status
    int     phase     ( ) const { return m_phase; }
    Status  status    ( ) const { return m_status; }
    int     iterations( ) const { return m_iterations; }
    bool    is_optimal( ) const { return status() == OPTIMAL; }
    
    
    // access to common denominator
    const ET&  variables_common_denominator( ) const { return d; }
    
    // access to basic variables
    int  number_of_basic_variables( ) const { return B.size(); }
    
    Basic_variable_index_iterator
    basic_variables_index_begin( ) const { return B.begin(); }
    
    Basic_variable_index_iterator
    basic_variables_index_end  ( ) const { return B.end(); }
    
    Basic_variable_numerator_iterator
    basic_variables_numerator_begin( ) const { return x_B.begin(); }
    
    Basic_variable_numerator_iterator
    basic_variables_numerator_end  ( ) const { return x_B.end(); }
    
    Basic_variable_value_iterator
    basic_variables_value_begin( ) const
        { return Basic_variable_value_iterator(
                     basic_variables_numerator_begin(),
                     Make_quotient( Compute_quotient(), d)); }
    
    Basic_variable_value_iterator
    basic_variables_value_end  ( ) const
        { return Basic_variable_value_iterator(
                     basic_variables_numerator_end  (),
                     Make_quotient( Compute_quotient(), d)); }
    
    
    bool  is_basic( int ii) const
        { CGAL_optimisation_precondition( ii >= 0);
          CGAL_optimisation_precondition( ii < ( phase() == 1
              ? number_of_variables()+number_of_constraints()
              : number_of_variables()));
          return ( in_B[ ii] >= 0); }
    
    
    // access to variables
    Variable_numerator_iterator
    variables_numerator_begin( ) const
        { return Variable_numerator_iterator(
                     in_B.begin(), Access_value_checked( x_B.begin(), et_0)); }
    
    Variable_numerator_iterator
    variables_numerator_end  ( ) const
        { return Variable_numerator_iterator(
                     in_B.end  (), Access_value_checked( x_B.begin(), et_0)); }
    
    Variable_value_iterator
    variables_value_begin( ) const
        { return Variable_value_iterator(
                     variables_numerator_begin(),
                     Make_quotient( Compute_quotient(), d)); }
    
    Variable_value_iterator
    variables_value_end  ( ) const
        { return Variable_value_iterator(
                     variables_numerator_end  (),
                     Make_quotient( Compute_quotient(), d)); }
    
    
    // access to current solution
    ET  solution_numerator  ( ) const;
    ET  solution_denominator( ) const { return d*d; }
    
    Quotient<ET>  solution( ) const
        { return Quotient<ET>( solution_numerator(), d*d); }
    
    
    // access to lambda
    Lambda_numerator_iterator
    lambda_numerator_begin( ) const { return lambda.begin(); }
    
    Lambda_numerator_iterator
    lambda_numerator_end  ( ) const { return lambda.end(); }
    
    Lambda_value_iterator
    lambda_value_begin( ) const
        { return Lambda_value_iterator(
                     lambda_numerator_begin(),
                     Make_quotient( Compute_quotient(), d)); }
    
    Lambda_value_iterator
    lambda_value_end  ( ) const
        { return Lambda_value_iterator(
                     lambda_numerator_end  (),
                     Make_quotient( Compute_quotient(), d)); }
    
    // access to auxiliary problem
    A_artificial_iterator  a_artificial_begin( ) const { return art_A.begin();}
    A_artificial_iterator  a_artificial_end  ( ) const { return art_A.end  ();}
    
    C_auxiliary_iterator   c_auxiliary_begin ( ) const { return aux_c.begin();}
    C_auxiliary_iterator   c_auxiliary_end   ( ) const { return aux_c.end  ();}
    
    
    // access to variables of dual LP
    ET  dual_variable( int i) const;
    
    
    // operations
    Status  pivot( )
        { CGAL_optimisation_precondition( phase() > 0);
          CGAL_optimisation_precondition( phase() < 3);
          pivot_step();
          return status(); }
    
    Status  solve( )
        { CGAL_optimisation_precondition( phase() > 0);
          while ( phase() < 3) { pivot_step(); }
          return status(); }
    
    
    // altering the pricing strategy
    const Pricing_strategy&  pricing_strategy() const;
    /*
      { return *strategyP; }
    */
    
    void  set_pricing_strategy( Pricing_strategy& pricing_strategy);
    /*
        {
            
            CGAL_optimisation_debug {
                vout2 << std::endl
                      << "-----------------------" << std::endl
                      << "Pricing Strategy Change" << std::endl
                      << "-----------------------" << std::endl;
            }
             
    
            strategyP = &pricing_strategy;
            strategyP->set( *this, vout2);
        }
    */
    
    // transition (to phase II)
    void  transition( );
    
    

  private:
    
    // some constants
    const NT                 nt_0, nt_1, nt_minus_1;
    const ET                 et_0, et_1, et_2;
                                              

    
    // verbose output streams
    CGAL::Verbose_ostream    vout1;     // used for some verbose output
    CGAL::Verbose_ostream    vout2;     // used for more verbose output
    CGAL::Verbose_ostream    vout3;     // used for very verbose output
    CGAL::Verbose_ostream    vout;      // used for any  verbose output
                                                                       
    
    
    // given QP
    int                      qp_n;      // number of variables
    int                      qp_m;      // number of constraints
    
    A_iterator               qp_A;      // constraint matrix
    B_iterator               qp_b;      // right-hand-side vector
    C_iterator               qp_c;      // objective vector
    D_iterator               qp_D;      // objective matrix
    
    // HACK
    unsigned int             max_basis; 
    
    // current status
    Indices                  B;         // basis
    
    Basis_inverse<ET,Is_lp>  inv_M_B;   // inverse of basis matrix
    
    Values                   lambda;
    Values                   x_B;       // solution restricted to basis
    
    int                      m_phase;   // phase of the Simplex method
    Status                   m_status;  // status of last pivot step
    int                      m_iterations;// number of pivot steps
    
    
    const ET&                d;         // reference to `inv_M_B.denominator()'
    
    
    Indices                  in_B;      // position in basis, -1 if non-basic
    
    
    // pricing strategy
    //Pricing_strategy_default strategy_default;
    Pricing_strategy*        strategyP;