important variables doc

.gitattributes

@@ -1642,12 +1642,13 @@ QP_solver/examples/QP_solver/first_lp.cpp -text
 QP_solver/examples/QP_solver/first_nonnegative_lp.cpp -text
 QP_solver/examples/QP_solver/first_nonnegative_qp.cpp -text
 QP_solver/examples/QP_solver/first_qp.cpp -text
+QP_solver/examples/QP_solver/important_variables.cpp -text
 QP_solver/examples/QP_solver/integer_qp_solver.cin -text
 QP_solver/examples/QP_solver/integer_qp_solver.data -text
-QP_solver/examples/QP_solver/is_in_convex_hull.h -text
-QP_solver/examples/QP_solver/is_in_convex_hull2.h -text
 QP_solver/examples/QP_solver/rational_qp_solver.cin -text
 QP_solver/examples/QP_solver/rational_qp_solver.data -text
+QP_solver/examples/QP_solver/solve_convex_hull_lp.h -text
+QP_solver/examples/QP_solver/solve_convex_hull_lp2.h -text
 QP_solver/include/CGAL/QP_functions.h -text
 QP_solver/include/CGAL/QP_models.h -text
 QP_solver/include/CGAL/QP_solution.h -text

QP_solver/doc_tex/QP_solver/main.tex

@@ -330,14 +330,19 @@ equivalently, if there are $\mu_1,\ldots,\mu_n$
 \]
 
 The linear program now tests for the existence of nonnegative $\mu_j$
-that satisfy the latter equation. Below is the code; it defines a function
+that satisfy the latter equation. Below is the code; it defines a
-that performs the test, given $p$ and $p_1,\ldots,p_n$ (through an iterator
+function that solves the linear program, given $p$ and
-range), and then calls it with $p_1=(0,0), p_2=(10,0), p_3=(0,10)$ fixed 
+$p_1,\ldots,p_n$ (through an iterator range).
-(they define a triangle) and all integral points $p$ in $[0,10]^2$. We know
-that $p$ is in the convex hull of $\{p_1,p_2,p_3\}$ if and only if its two
-coordinates sum up to $10$ at most.
 
-\ccIncludeExampleCode{QP_solver/is_in_convex_hull.cpp}
+\ccIncludeExampleCode{QP_solver/solve_convex_hull_lp.h}
+
+To see this in action, let us call it with $p_1=(0,0), p_2=(10,0),
+p_3=(0,10)$ fixed (they define a triangle) and all integral points $p$
+in $[0,10]^2$. We know that $p$ is in the convex hull of
+$\{p_1,p_2,p_3\}$ if and only if its two coordinates sum up to $10$ at
+most.
+
+\ccIncludeExampleCode{QP_solver/convex_hull_containment.cpp}
 
 \subsection{Using Makers}
 You already noticed in the previous example that the actual 
@@ -354,9 +359,25 @@ You can avoid the explicit construction of the type
 if you only need an expression of it, e.g.\ to pass it directly 
 as an argument to the solving function. Here is how this works.
 
-\ccIncludeExampleCode{QP_solver/is_in_convex_hull2.cpp}
+\ccIncludeExampleCode{QP_solver/solve_convex_hull_lp2.h}
 
 \section{The Important Variables}
+If we have an optimal solution $x^*$ of a linear or quadratic program,
+the ``important'' variables are typically the ones that are not on
+their bounds. In case of a nonnegative program, these are the nonzero
+variables.  Going back to the example of the previous Section
+\ref{sec:QP-working_from_iterators}, we can easily interpret their
+importance: the nonzero variables correspond to points $p_j$ that
+actually contribute to the convex combination that yields $p$. 
+
+The following example shows how we can access the important variables.
+We generate a set of points on a line and then find the ones that 
+contribute to the convex combinations of intermediate points.
+
+\ccIncludeExampleCode{QP_solver/important_variables.cpp}
+
+Not surprisingly, only two points contribute to any convex
+combination.
 
 \section{Working from MPS Files}
 

QP_solver/examples/QP_solver/convex_hull_containment.cpp

@@ -3,11 +3,20 @@
 #include <vector>
 #include <CGAL/Homogeneous_d.h>
 #include <CGAL/MP_Float.h>
-#include "is_in_convex_hull.h"
+#include "solve_convex_hull_lp.h"
 
 typedef CGAL::Homogeneous_d<double> Kernel_d;
 typedef Kernel_d::Point_d Point_d;
 
+bool is_in_convex_hull (const Point_d& p, 
+			std::vector<Point_d>::const_iterator begin,
+			std::vector<Point_d>::const_iterator end)
+{
+  CGAL::Quadratic_program_solution<CGAL::MP_Float> s =
+    solve_convex_hull_lp (p, begin, end, CGAL::MP_Float());
+  return s.status() != CGAL::QP_INFEASIBLE;
+}
+
 int main()
 {
   std::vector<Point_d> points;
@@ -19,7 +28,7 @@ int main()
     for (int j=0; j<=10; ++j) {
       // (i,j) is in the simplex iff i+j <= 10 
       bool contained = is_in_convex_hull 
-	(Point_d (i, j), points.begin(), points.end(), CGAL::MP_Float());
+	(Point_d (i, j), points.begin(), points.end());
       assert (contained == (i+j<=10));
     }
 	

QP_solver/examples/QP_solver/convex_hull_containment2.cpp

@@ -3,11 +3,20 @@
 #include <vector>
 #include <CGAL/Homogeneous_d.h>
 #include <CGAL/MP_Float.h>
-#include "is_in_convex_hull2.h"
+#include "solve_convex_hull_lp2.h"
 
 typedef CGAL::Homogeneous_d<double> Kernel_d;
 typedef Kernel_d::Point_d Point_d;
 
+bool is_in_convex_hull (const Point_d& p, 
+			std::vector<Point_d>::const_iterator begin,
+			std::vector<Point_d>::const_iterator end)
+{
+  CGAL::Quadratic_program_solution<CGAL::MP_Float> s =
+    solve_convex_hull_lp (p, begin, end, CGAL::MP_Float());
+  return s.status() != CGAL::QP_INFEASIBLE;
+}
+
 int main()
 {
   std::vector<Point_d> points;
@@ -19,7 +28,7 @@ int main()
     for (int j=0; j<=10; ++j) {
       // (i,j) is in the simplex iff i+j <= 10 
       bool contained = is_in_convex_hull 
-	(Point_d (i, j), points.begin(), points.end(), CGAL::MP_Float());
+	(Point_d (i, j), points.begin(), points.end());
       assert (contained == (i+j<=10));
     }
 	

QP_solver/examples/QP_solver/important_variables.cpp

@@ -0,0 +1,37 @@
+// Example: find the points that contribute to a convex combination
+#include <cassert>
+#include <vector>
+#include <CGAL/Cartesian_d.h>
+#include <CGAL/MP_Float.h>
+#include "solve_convex_hull_lp2.h"
+
+typedef CGAL::Cartesian_d<double> Kernel_d;
+typedef Kernel_d::Point_d Point_d;
+typedef CGAL::Quadratic_program_solution<CGAL::MP_Float> Solution;
+
+int main()
+{
+  std::vector<Point_d> points;
+  // convex hull: line spanned by {(0,0), (1,0),..., (9,0)} 
+  for (int j=0; j<10; ++j) 
+    points.push_back (Point_d (j, 0));
+   
+  for (double f=0.5; f<10; ++f) {
+    Point_d p (f, 0.0);
+    Solution s = solve_convex_hull_lp 
+      (p, points.begin(), points.end(), CGAL::MP_Float());
+    std::cout << p; 
+    if (s.status() == CGAL::QP_INFEASIBLE) 
+      std::cout << " is not in the convex hull." << std::endl;
+    else {
+      std::cout << " is a convex combination of the points ";
+      Solution::Index_iterator it = s.basic_variable_indices_begin();
+      Solution::Index_iterator end = s.basic_variable_indices_end();
+      for (; it != end; ++it) 
+	std::cout << points[*it] << " ";
+      std::cout << std::endl;
+    }
+  }
+	
+  return 0;
+}

QP_solver/examples/QP_solver/solve_convex_hull_lp.h

@@ -36,8 +36,6 @@ solve_convex_hull_lp (const Point_d& p,
   typedef
     CGAL::Nonnegative_linear_program_from_iterators<A_it, B_it, R_it, C_it>
     Program;
-  // the solution type
-  typedef CGAL::Quadratic_program_solution<ET> Solution;
 
   // ok, we are prepared now: construct program and solve it
   Program lp (end-begin,        // number of variables 
@@ -46,7 +44,3 @@ solve_convex_hull_lp (const Point_d& p,
 	      R_it (CGAL::EQUAL), C_it (0));
   return CGAL::solve_nonnegative_linear_program (lp, ET(0));
 }
-  
-  // p is in the convex hull if and only if the program is feasible
-  return (s.status() != CGAL::QP_INFEASIBLE);
-}

QP_solver/examples/QP_solver/solve_convex_hull_lp2.h

@@ -15,18 +15,16 @@ struct Homogeneous_begin  {
 // function to test whether point is in the convex hull of other points;
 // the type ET is an exact type used for the computations
 template <class Point_d, class RandomAccessIterator, class ET>
-bool is_in_convex_hull
+CGAL::Quadratic_program_solution<ET> 
-(const Point_d& p, 
+solve_convex_hull_lp (const Point_d& p, 
- RandomAccessIterator begin, RandomAccessIterator end, const ET& dummy) 
+		      RandomAccessIterator begin, 
+		      RandomAccessIterator end, const ET& dummy) 
 {
   // the right-hand side type
   typedef typename Point_d::Homogeneous_const_iterator B_it;
 
-  // the solution type
-  typedef CGAL::Quadratic_program_solution<ET> Solution;
-
   // construct program and solve it
-  Solution s = CGAL::solve_nonnegative_linear_program 
+  return CGAL::solve_nonnegative_linear_program 
     (CGAL::make_nonnegative_linear_program_from_iterators
      (end-begin,                                                         // n
       p.dimension()+1,                                                   // m
@@ -36,7 +34,5 @@ bool is_in_convex_hull
       CGAL::Const_oneset_iterator<CGAL::Comparison_result>(CGAL::EQUAL), // ~
       CGAL::Const_oneset_iterator
       <typename std::iterator_traits<B_it>::value_type> (0)), ET(0));    // c
-
-  // p is in the convex hull if and only if the program is feasible
-  return (s.status() != CGAL::QP_INFEASIBLE);
 }
+

QP_solver/include/CGAL/QP_solver.h

@@ -1277,9 +1277,9 @@ private:
   bool  check_basis_inverse( Tag_false is_linear);
 
   // diagnostic output
-  void  print_program ( );
+  void  print_program ( ) const;
-  void  print_basis   ( );
+  void  print_basis   ( ) const;
-  void  print_solution( );
+  void  print_solution( ) const;
   void  print_ratio_1_original(int k, const ET& x_k, const ET& q_k);
   void  print_ratio_1_slack(int k, const ET& x_k, const ET& q_k);
 

QP_solver/include/CGAL/QP_solver/QP_solver_impl.h

@@ -2972,7 +2972,7 @@ set_verbosity( int verbose, std::ostream& stream)
 
 template < typename Q, typename ET, typename Tags >
 void  QP_solver<Q, ET, Tags>::
-print_program( )
+print_program( ) const
 {
     int  row, i;
 
@@ -3050,7 +3050,7 @@ print_program( )
 
 template < typename Q, typename ET, typename Tags >
 void  QP_solver<Q, ET, Tags>::
-print_basis( )
+print_basis( ) const
 {
     char label;
     vout1 << "  basis: ";
@@ -3069,7 +3069,8 @@ print_basis( )
 	if ( has_ineq) {
 	    vout2.out() << std::endl
 			<< "basic constraints:  ";
-	    for (Index_iterator i_it = C.begin(); i_it != C.end(); ++i_it) {
+	    for (Index_const_iterator i_it = 
+		   C.begin(); i_it != C.end(); ++i_it) {
 	        label = (qp_r[*i_it] == CGAL::EQUAL) ? 'e' : 'i';
 		    vout2.out() << *i_it << ":" << label << " ";
 	    }
@@ -3103,7 +3104,7 @@ print_basis( )
 
 template < typename Q, typename ET, typename Tags >
 void  QP_solver<Q, ET, Tags>::
-print_solution( )
+print_solution( ) const
 {
     if ( vout3.verbose()) {
 	   vout3.out() << std::endl