[PPL-devel] [GIT] ppl/ppl(master): Improved method Polyhedron::relation_with( const Congruence&) const.

[Enea Zaffanella]

Module: ppl/ppl
Branch: master
Commit: eca9369a8272645ba0db81f3244004691623bcba
URL:    http://www.cs.unipr.it/git/gitweb.cgi?p=ppl/ppl.git;a=commit;h=eca9369a8272645ba0db81f3244004691623bcba

Author: Enea Zaffanella <zaffanella at cs.unipr.it>
Date:   Fri May 15 15:40:35 2009 +0200

Improved method Polyhedron::relation_with(const Congruence&) const.
Avoided several temporary objects; added a FIXME asking for the rounding mode.

---

 src/Polyhedron_public.cc |   71 +++++++++++++++++++++------------------------
 1 files changed, 33 insertions(+), 38 deletions(-)

diff --git a/src/Polyhedron_public.cc b/src/Polyhedron_public.cc
index 1e8f2c5..a805773 100644
--- a/src/Polyhedron_public.cc
+++ b/src/Polyhedron_public.cc
@@ -303,69 +303,64 @@ PPL::Polyhedron::relation_with(const Congruence& cg) const {
       && Poly_Con_Relation::is_included()
       && Poly_Con_Relation::is_disjoint();
 
-  // Find the equality corresponding to the congruence (ignoring the modulus).
-  Linear_Expression expr;
-  for (dimension_type i = cg_space_dim; i-- > 0; ) {
-    const Variable v(i);
-    expr += cg.coefficient(v) * v;
-  }
-  expr += cg.inhomogeneous_term();
+  // Build the equality corresponding to the congruence (ignoring the modulus).
+  Linear_Expression expr = Linear_Expression(cg);
   Constraint c(expr == 0);
 
   // The polyhedron is non-empty so that there exists a point.
-  // For an arbitrary generator point find the scalar product with
+  // For an arbitrary generator point, compute the scalar product with
   // the equality.
   PPL_DIRTY_TEMP_COEFFICIENT(point_val);
 
-  for (Generator_System::const_iterator g = gen_sys.begin(),
-         gen_sys_end = gen_sys.end(); g != gen_sys_end; ++g) {
-    if (g->is_point()) {
-      Scalar_Products::assign(point_val, c, *g);
+  for (Generator_System::const_iterator gs_i = gen_sys.begin(),
+         gs_end = gen_sys.end(); gs_i != gs_end; ++gs_i) {
+    if (gs_i->is_point()) {
+      Scalar_Products::assign(point_val, c, *gs_i);
       break;
     }
   }
 
-  // Find the 2 hyperplanes that satisfies the congruence and are
-  // nearest to the point such that the point lies on or between these
-  // hyperplanes.
+  // Find 2 hyperplanes that satisfy the congruence and are near to
+  // the generating point (so that the point lies on or between these
+  // hyperplanes).
   // Then use the relations between the polyhedron and the corresponding
-  // half-spaces to determine its relation with the congruence.
+  // halfspaces to determine its relation with the congruence.
   const Coefficient& modulus = cg.modulus();
 
-  PPL_DIRTY_TEMP_COEFFICIENT(div);
-  div = modulus;
-
+  // FIXME: specify rounding mode.
   PPL_DIRTY_TEMP_COEFFICIENT(nearest);
-  nearest = (point_val / div) * div;
+  nearest = (point_val / modulus) * modulus;
 
   point_val -= nearest;
   expr -= nearest;
   if (point_val == 0)
      return relation_with(expr == 0);
 
-  Linear_Expression next_expr;
-  if (point_val > 0) {
-     next_expr = expr - modulus;
-  }
-  else {
-    expr = - (expr);
-    next_expr = expr - modulus;
-  }
+  // Build first halfspace.
+  const bool positive = (point_val > 0);
+  Constraint first_halfspace = positive ? (expr >= 0) : (expr <= 0);
 
-  Poly_Con_Relation relations = relation_with(expr >= 0);
-  assert(!relations.implies(Poly_Con_Relation::saturates())
-	 && !relations.implies(Poly_Con_Relation::is_disjoint()));
-  if (relations.implies(Poly_Con_Relation::strictly_intersects()))
+  Poly_Con_Relation first_rels = relation_with(first_halfspace);
+  assert(!first_rels.implies(Poly_Con_Relation::saturates())
+	 && !first_rels.implies(Poly_Con_Relation::is_disjoint()));
+  if (first_rels.implies(Poly_Con_Relation::strictly_intersects()))
     return Poly_Con_Relation::strictly_intersects();
 
-  assert(relations == Poly_Con_Relation::is_included());
-  Poly_Con_Relation next_relations = relation_with(next_expr <= 0);
-  assert(!next_relations.implies(Poly_Con_Relation::saturates())
-	 && !next_relations.implies(Poly_Con_Relation::is_disjoint()));
-  if (next_relations.implies(Poly_Con_Relation::strictly_intersects()))
+  // Build second halfspace.
+  if (positive)
+    expr -= modulus;
+  else
+    expr += modulus;
+  Constraint second_halfspace = positive ? (expr <= 0) : (expr >= 0);
+
+  assert(first_rels == Poly_Con_Relation::is_included());
+  Poly_Con_Relation second_rels = relation_with(second_halfspace);
+  assert(!second_rels.implies(Poly_Con_Relation::saturates())
+	 && !second_rels.implies(Poly_Con_Relation::is_disjoint()));
+  if (second_rels.implies(Poly_Con_Relation::strictly_intersects()))
     return Poly_Con_Relation::strictly_intersects();
 
-  assert(next_relations == Poly_Con_Relation::is_included());
+  assert(second_rels == Poly_Con_Relation::is_included());
   return Poly_Con_Relation::is_disjoint();
 }