1 /* BD_Shape class declaration.
2    Copyright (C) 2001-2010 Roberto Bagnara <bagnara@cs.unipr.it>
3    Copyright (C) 2010-2016 BUGSENG srl (http://bugseng.com)
4 
5 This file is part of the Parma Polyhedra Library (PPL).
6 
7 The PPL is free software; you can redistribute it and/or modify it
8 under the terms of the GNU General Public License as published by the
9 Free Software Foundation; either version 3 of the License, or (at your
10 option) any later version.
11 
12 The PPL is distributed in the hope that it will be useful, but WITHOUT
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
15 for more details.
16 
17 You should have received a copy of the GNU General Public License
18 along with this program; if not, write to the Free Software Foundation,
19 Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02111-1307, USA.
20 
21 For the most up-to-date information see the Parma Polyhedra Library
22 site: http://bugseng.com/products/ppl/ . */
23 
24 #ifndef PPL_BD_Shape_defs_hh
25 #define PPL_BD_Shape_defs_hh 1
26 
27 #include "BD_Shape_types.hh"
28 #include "globals_defs.hh"
29 #include "Constraint_types.hh"
30 #include "Generator_types.hh"
31 #include "Congruence_types.hh"
32 #include "Linear_Expression_types.hh"
33 #include "Constraint_System_types.hh"
34 #include "Generator_System_types.hh"
35 #include "Congruence_System_types.hh"
36 #include "Poly_Con_Relation_types.hh"
37 #include "Poly_Gen_Relation_types.hh"
38 #include "Polyhedron_types.hh"
39 #include "Box_types.hh"
40 #include "Grid_types.hh"
41 #include "Octagonal_Shape_types.hh"
42 #include "Variable_defs.hh"
43 #include "Variables_Set_types.hh"
44 #include "DB_Matrix_defs.hh"
45 #include "DB_Row_defs.hh"
46 #include "Checked_Number_defs.hh"
47 #include "WRD_coefficient_types_defs.hh"
48 #include "Bit_Matrix_defs.hh"
49 #include "Coefficient_defs.hh"
50 #include "Interval_types.hh"
51 #include "Linear_Form_types.hh"
52 #include <cstddef>
53 #include <iosfwd>
54 #include <vector>
55 
56 namespace Parma_Polyhedra_Library {
57 
58 namespace IO_Operators {
59 
60 //! Output operator.
61 /*! \relates Parma_Polyhedra_Library::BD_Shape
62   Writes a textual representation of \p bds on \p s:
63   <CODE>false</CODE> is written if \p bds is an empty polyhedron;
64   <CODE>true</CODE> is written if \p bds is the universe polyhedron;
65   a system of constraints defining \p bds is written otherwise,
66   all constraints separated by ", ".
67 */
68 template <typename T>
69 std::ostream&
70 operator<<(std::ostream& s, const BD_Shape<T>& bds);
71 
72 } // namespace IO_Operators
73 
74 //! Swaps \p x with \p y.
75 /*! \relates BD_Shape */
76 template <typename T>
77 void swap(BD_Shape<T>& x, BD_Shape<T>& y);
78 
79 //! Returns <CODE>true</CODE> if and only if \p x and \p y are the same BDS.
80 /*! \relates BD_Shape
81   Note that \p x and \p y may be dimension-incompatible shapes:
82   in this case, the value <CODE>false</CODE> is returned.
83 */
84 template <typename T>
85 bool operator==(const BD_Shape<T>& x, const BD_Shape<T>& y);
86 
87 //! Returns <CODE>true</CODE> if and only if \p x and \p y are not the same BDS.
88 /*! \relates BD_Shape
89   Note that \p x and \p y may be dimension-incompatible shapes:
90   in this case, the value <CODE>true</CODE> is returned.
91 */
92 template <typename T>
93 bool operator!=(const BD_Shape<T>& x, const BD_Shape<T>& y);
94 
95 //! Computes the rectilinear (or Manhattan) distance between \p x and \p y.
96 /*! \relates BD_Shape
97   If the rectilinear distance between \p x and \p y is defined,
98   stores an approximation of it into \p r and returns <CODE>true</CODE>;
99   returns <CODE>false</CODE> otherwise.
100 
101   The direction of the approximation is specified by \p dir.
102 
103   All computations are performed using variables of type
104   <CODE>Checked_Number\<To, Extended_Number_Policy\></CODE>.
105 */
106 template <typename To, typename T>
107 bool rectilinear_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
108                                  const BD_Shape<T>& x,
109                                  const BD_Shape<T>& y,
110                                  Rounding_Dir dir);
111 
112 //! Computes the rectilinear (or Manhattan) distance between \p x and \p y.
113 /*! \relates BD_Shape
114   If the rectilinear distance between \p x and \p y is defined,
115   stores an approximation of it into \p r and returns <CODE>true</CODE>;
116   returns <CODE>false</CODE> otherwise.
117 
118   The direction of the approximation is specified by \p dir.
119 
120   All computations are performed using variables of type
121   <CODE>Checked_Number\<Temp, Extended_Number_Policy\></CODE>.
122 */
123 template <typename Temp, typename To, typename T>
124 bool rectilinear_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
125                                  const BD_Shape<T>& x,
126                                  const BD_Shape<T>& y,
127                                  Rounding_Dir dir);
128 
129 //! Computes the rectilinear (or Manhattan) distance between \p x and \p y.
130 /*! \relates BD_Shape
131   If the rectilinear distance between \p x and \p y is defined,
132   stores an approximation of it into \p r and returns <CODE>true</CODE>;
133   returns <CODE>false</CODE> otherwise.
134 
135   The direction of the approximation is specified by \p dir.
136 
137   All computations are performed using the temporary variables
138   \p tmp0, \p tmp1 and \p tmp2.
139 */
140 template <typename Temp, typename To, typename T>
141 bool rectilinear_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
142                                  const BD_Shape<T>& x,
143                                  const BD_Shape<T>& y,
144                                  Rounding_Dir dir,
145                                  Temp& tmp0,
146                                  Temp& tmp1,
147                                  Temp& tmp2);
148 
149 //! Computes the euclidean distance between \p x and \p y.
150 /*! \relates BD_Shape
151   If the euclidean distance between \p x and \p y is defined,
152   stores an approximation of it into \p r and returns <CODE>true</CODE>;
153   returns <CODE>false</CODE> otherwise.
154 
155   The direction of the approximation is specified by \p dir.
156 
157   All computations are performed using variables of type
158   <CODE>Checked_Number\<To, Extended_Number_Policy\></CODE>.
159 */
160 template <typename To, typename T>
161 bool euclidean_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
162                                const BD_Shape<T>& x,
163                                const BD_Shape<T>& y,
164                                Rounding_Dir dir);
165 
166 //! Computes the euclidean distance between \p x and \p y.
167 /*! \relates BD_Shape
168   If the euclidean distance between \p x and \p y is defined,
169   stores an approximation of it into \p r and returns <CODE>true</CODE>;
170   returns <CODE>false</CODE> otherwise.
171 
172   The direction of the approximation is specified by \p dir.
173 
174   All computations are performed using variables of type
175   <CODE>Checked_Number\<Temp, Extended_Number_Policy\></CODE>.
176 */
177 template <typename Temp, typename To, typename T>
178 bool euclidean_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
179                                const BD_Shape<T>& x,
180                                const BD_Shape<T>& y,
181                                Rounding_Dir dir);
182 
183 //! Computes the euclidean distance between \p x and \p y.
184 /*! \relates BD_Shape
185   If the euclidean distance between \p x and \p y is defined,
186   stores an approximation of it into \p r and returns <CODE>true</CODE>;
187   returns <CODE>false</CODE> otherwise.
188 
189   The direction of the approximation is specified by \p dir.
190 
191   All computations are performed using the temporary variables
192   \p tmp0, \p tmp1 and \p tmp2.
193 */
194 template <typename Temp, typename To, typename T>
195 bool euclidean_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
196                                const BD_Shape<T>& x,
197                                const BD_Shape<T>& y,
198                                Rounding_Dir dir,
199                                Temp& tmp0,
200                                Temp& tmp1,
201                                Temp& tmp2);
202 
203 //! Computes the \f$L_\infty\f$ distance between \p x and \p y.
204 /*! \relates BD_Shape
205   If the \f$L_\infty\f$ distance between \p x and \p y is defined,
206   stores an approximation of it into \p r and returns <CODE>true</CODE>;
207   returns <CODE>false</CODE> otherwise.
208 
209   The direction of the approximation is specified by \p dir.
210 
211   All computations are performed using variables of type
212   <CODE>Checked_Number\<To, Extended_Number_Policy\></CODE>.
213 */
214 template <typename To, typename T>
215 bool l_infinity_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
216                                 const BD_Shape<T>& x,
217                                 const BD_Shape<T>& y,
218                                 Rounding_Dir dir);
219 
220 //! Computes the \f$L_\infty\f$ distance between \p x and \p y.
221 /*! \relates BD_Shape
222   If the \f$L_\infty\f$ distance between \p x and \p y is defined,
223   stores an approximation of it into \p r and returns <CODE>true</CODE>;
224   returns <CODE>false</CODE> otherwise.
225 
226   The direction of the approximation is specified by \p dir.
227 
228   All computations are performed using variables of type
229   <CODE>Checked_Number\<Temp, Extended_Number_Policy\></CODE>.
230 */
231 template <typename Temp, typename To, typename T>
232 bool l_infinity_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
233                                 const BD_Shape<T>& x,
234                                 const BD_Shape<T>& y,
235                                 Rounding_Dir dir);
236 
237 //! Computes the \f$L_\infty\f$ distance between \p x and \p y.
238 /*! \relates BD_Shape
239   If the \f$L_\infty\f$ distance between \p x and \p y is defined,
240   stores an approximation of it into \p r and returns <CODE>true</CODE>;
241   returns <CODE>false</CODE> otherwise.
242 
243   The direction of the approximation is specified by \p dir.
244 
245   All computations are performed using the temporary variables
246   \p tmp0, \p tmp1 and \p tmp2.
247 */
248 template <typename Temp, typename To, typename T>
249 bool l_infinity_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
250                                 const BD_Shape<T>& x,
251                                 const BD_Shape<T>& y,
252                                 Rounding_Dir dir,
253                                 Temp& tmp0,
254                                 Temp& tmp1,
255                                 Temp& tmp2);
256 
257 // This class contains some helper functions that need to be friends of
258 // Linear_Expression.
259 class BD_Shape_Helpers {
260 public:
261   #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS
262   //! Decodes the constraint \p c as a bounded difference.
263   /*! \relates BD_Shape
264     \return
265     <CODE>true</CODE> if the constraint \p c is a
266     \ref Bounded_Difference_Shapes "bounded difference";
267     <CODE>false</CODE> otherwise.
268 
269     \param c
270     The constraint to be decoded.
271 
272     \param c_num_vars
273     If <CODE>true</CODE> is returned, then it will be set to the number
274     of variables having a non-zero coefficient. The only legal values
275     will therefore be 0, 1 and 2.
276 
277     \param c_first_var
278     If <CODE>true</CODE> is returned and if \p c_num_vars is not set to 0,
279     then it will be set to the index of the first variable having
280     a non-zero coefficient in \p c.
281 
282     \param c_second_var
283     If <CODE>true</CODE> is returned and if \p c_num_vars is set to 2,
284     then it will be set to the index of the second variable having
285     a non-zero coefficient in \p c. If \p c_num_vars is set to 1, this must be
286     0.
287 
288     \param c_coeff
289     If <CODE>true</CODE> is returned and if \p c_num_vars is not set to 0,
290     then it will be set to the value of the first non-zero coefficient
291     in \p c.
292   */
293   #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS)
294   static bool extract_bounded_difference(const Constraint& c,
295                                          dimension_type& c_num_vars,
296                                          dimension_type& c_first_var,
297                                          dimension_type& c_second_var,
298                                          Coefficient& c_coeff);
299 };
300 
301 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS
302 //! Extracts leader indices from the predecessor relation.
303 /*! \relates BD_Shape */
304 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS)
305 void compute_leader_indices(const std::vector<dimension_type>& predecessor,
306                             std::vector<dimension_type>& indices);
307 
308 } // namespace Parma_Polyhedra_Library
309 
310 //! A bounded difference shape.
311 /*! \ingroup PPL_CXX_interface
312   The class template BD_Shape<T> allows for the efficient representation
313   of a restricted kind of <EM>topologically closed</EM> convex polyhedra
314   called <EM>bounded difference shapes</EM> (BDSs, for short).
315   The name comes from the fact that the closed affine half-spaces that
316   characterize the polyhedron can be expressed by constraints of the form
317   \f$\pm x_i \leq k\f$ or \f$x_i - x_j \leq k\f$, where the inhomogeneous
318   term \f$k\f$ is a rational number.
319 
320   Based on the class template type parameter \p T, a family of extended
321   numbers is built and used to approximate the inhomogeneous term of
322   bounded differences. These extended numbers provide a representation
323   for the value \f$+\infty\f$, as well as <EM>rounding-aware</EM>
324   implementations for several arithmetic functions.
325   The value of the type parameter \p T may be one of the following:
326     - a bounded precision integer type (e.g., \c int32_t or \c int64_t);
327     - a bounded precision floating point type (e.g., \c float or \c double);
328     - an unbounded integer or rational type, as provided by GMP
329       (i.e., \c mpz_class or \c mpq_class).
330 
331   The user interface for BDSs is meant to be as similar as possible to
332   the one developed for the polyhedron class C_Polyhedron.
333 
334   The domain of BD shapes <EM>optimally supports</EM>:
335     - tautological and inconsistent constraints and congruences;
336     - bounded difference constraints;
337     - non-proper congruences (i.e., equalities) that are expressible
338       as bounded-difference constraints.
339 
340   Depending on the method, using a constraint or congruence that is not
341   optimally supported by the domain will either raise an exception or
342   result in a (possibly non-optimal) upward approximation.
343 
344   A constraint is a bounded difference if it has the form
345     \f[
346       a_i x_i - a_j x_j \relsym b
347     \f]
348   where \f$\mathord{\relsym} \in \{ \leq, =, \geq \}\f$ and
349   \f$a_i\f$, \f$a_j\f$, \f$b\f$ are integer coefficients such that
350   \f$a_i = 0\f$, or \f$a_j = 0\f$, or \f$a_i = a_j\f$.
351   The user is warned that the above bounded difference Constraint object
352   will be mapped into a \e correct and \e optimal approximation that,
353   depending on the expressive power of the chosen template argument \p T,
354   may loose some precision. Also note that strict constraints are not
355   bounded differences.
356 
357   For instance, a Constraint object encoding \f$3x - 3y \leq 1\f$ will be
358   approximated by:
359     - \f$x - y \leq 1\f$,
360       if \p T is a (bounded or unbounded) integer type;
361     - \f$x - y \leq \frac{1}{3}\f$,
362       if \p T is the unbounded rational type \c mpq_class;
363     - \f$x - y \leq k\f$, where \f$k > \frac{1}{3}\f$,
364       if \p T is a floating point type (having no exact representation
365       for \f$\frac{1}{3}\f$).
366 
367   On the other hand, depending from the context, a Constraint object
368   encoding \f$3x - y \leq 1\f$ will be either upward approximated
369   (e.g., by safely ignoring it) or it will cause an exception.
370 
371   In the following examples it is assumed that the type argument \p T
372   is one of the possible instances listed above and that variables
373   <CODE>x</CODE>, <CODE>y</CODE> and <CODE>z</CODE> are defined
374   (where they are used) as follows:
375   \code
376     Variable x(0);
377     Variable y(1);
378     Variable z(2);
379   \endcode
380 
381   \par Example 1
382   The following code builds a BDS corresponding to a cube in \f$\Rset^3\f$,
383   given as a system of constraints:
384   \code
385     Constraint_System cs;
386     cs.insert(x >= 0);
387     cs.insert(x <= 1);
388     cs.insert(y >= 0);
389     cs.insert(y <= 1);
390     cs.insert(z >= 0);
391     cs.insert(z <= 1);
392     BD_Shape<T> bd(cs);
393   \endcode
394   Since only those constraints having the syntactic form of a
395   <EM>bounded difference</EM> are optimally supported, the following code
396   will throw an exception (caused by constraints 7, 8 and 9):
397   \code
398     Constraint_System cs;
399     cs.insert(x >= 0);
400     cs.insert(x <= 1);
401     cs.insert(y >= 0);
402     cs.insert(y <= 1);
403     cs.insert(z >= 0);
404     cs.insert(z <= 1);
405     cs.insert(x + y <= 0);      // 7
406     cs.insert(x - z + x >= 0);  // 8
407     cs.insert(3*z - y <= 1);    // 9
408     BD_Shape<T> bd(cs);
409   \endcode
410 */
411 template <typename T>
412 class Parma_Polyhedra_Library::BD_Shape {
413 private:
414   /*! \brief
415     The (extended) numeric type of the inhomogeneous term of
416     the inequalities defining a BDS.
417   */
418 #ifndef NDEBUG
419   typedef Checked_Number<T, Debug_WRD_Extended_Number_Policy> N;
420 #else
421   typedef Checked_Number<T, WRD_Extended_Number_Policy> N;
422 #endif
423 
424 public:
425   //! The numeric base type upon which bounded differences are built.
426   typedef T coefficient_type_base;
427 
428   /*! \brief
429     The (extended) numeric type of the inhomogeneous term of the
430     inequalities defining a BDS.
431   */
432   typedef N coefficient_type;
433 
434   //! Returns the maximum space dimension that a BDS can handle.
435   static dimension_type max_space_dimension();
436 
437   /*! \brief
438     Returns \c false indicating that this domain cannot recycle constraints.
439   */
440   static bool can_recycle_constraint_systems();
441 
442   /*! \brief
443     Returns \c false indicating that this domain cannot recycle congruences.
444   */
445   static bool can_recycle_congruence_systems();
446 
447   //! \name Constructors, Assignment, Swap and Destructor
448   //@{
449 
450   //! Builds a universe or empty BDS of the specified space dimension.
451   /*!
452     \param num_dimensions
453     The number of dimensions of the vector space enclosing the BDS;
454 
455     \param kind
456     Specifies whether the universe or the empty BDS has to be built.
457   */
458   explicit BD_Shape(dimension_type num_dimensions = 0,
459                     Degenerate_Element kind = UNIVERSE);
460 
461   //! Ordinary copy constructor.
462   /*!
463     The complexity argument is ignored.
464   */
465   BD_Shape(const BD_Shape& y,
466            Complexity_Class complexity = ANY_COMPLEXITY);
467 
468   //! Builds a conservative, upward approximation of \p y.
469   /*!
470     The complexity argument is ignored.
471   */
472   template <typename U>
473   explicit BD_Shape(const BD_Shape<U>& y,
474                     Complexity_Class complexity = ANY_COMPLEXITY);
475 
476   //! Builds a BDS from the system of constraints \p cs.
477   /*!
478     The BDS inherits the space dimension of \p cs.
479 
480     \param cs
481     A system of BD constraints.
482 
483     \exception std::invalid_argument
484     Thrown if \p cs contains a constraint which is not optimally supported
485     by the BD shape domain.
486   */
487   explicit BD_Shape(const Constraint_System& cs);
488 
489   //! Builds a BDS from a system of congruences.
490   /*!
491     The BDS inherits the space dimension of \p cgs
492 
493     \param cgs
494     A system of congruences.
495 
496     \exception std::invalid_argument
497     Thrown if \p cgs contains congruences which are not optimally
498     supported by the BD shape domain.
499   */
500   explicit BD_Shape(const Congruence_System& cgs);
501 
502   //! Builds a BDS from the system of generators \p gs.
503   /*!
504     Builds the smallest BDS containing the polyhedron defined by \p gs.
505     The BDS inherits the space dimension of \p gs.
506 
507     \exception std::invalid_argument
508     Thrown if the system of generators is not empty but has no points.
509   */
510   explicit BD_Shape(const Generator_System& gs);
511 
512   //! Builds a BDS from the polyhedron \p ph.
513   /*!
514     Builds a BDS containing \p ph using algorithms whose complexity
515     does not exceed the one specified by \p complexity.  If
516     \p complexity is \p ANY_COMPLEXITY, then the BDS built is the
517     smallest one containing \p ph.
518   */
519   explicit BD_Shape(const Polyhedron& ph,
520                     Complexity_Class complexity = ANY_COMPLEXITY);
521 
522   //! Builds a BDS out of a box.
523   /*!
524     The BDS inherits the space dimension of the box.
525     The built BDS is the most precise BDS that includes the box.
526 
527     \param box
528     The box representing the BDS to be built.
529 
530     \param complexity
531     This argument is ignored as the algorithm used has
532     polynomial complexity.
533 
534     \exception std::length_error
535     Thrown if the space dimension of \p box exceeds the maximum
536     allowed space dimension.
537   */
538   template <typename Interval>
539   explicit BD_Shape(const Box<Interval>& box,
540                     Complexity_Class complexity = ANY_COMPLEXITY);
541 
542   //! Builds a BDS out of a grid.
543   /*!
544     The BDS inherits the space dimension of the grid.
545     The built BDS is the most precise BDS that includes the grid.
546 
547     \param grid
548     The grid used to build the BDS.
549 
550     \param complexity
551     This argument is ignored as the algorithm used has
552     polynomial complexity.
553 
554     \exception std::length_error
555     Thrown if the space dimension of \p grid exceeds the maximum
556     allowed space dimension.
557   */
558   explicit BD_Shape(const Grid& grid,
559                     Complexity_Class complexity = ANY_COMPLEXITY);
560 
561   //! Builds a BDS from an octagonal shape.
562   /*!
563     The BDS inherits the space dimension of the octagonal shape.
564     The built BDS is the most precise BDS that includes the octagonal shape.
565 
566     \param os
567     The octagonal shape used to build the BDS.
568 
569     \param complexity
570     This argument is ignored as the algorithm used has
571     polynomial complexity.
572 
573     \exception std::length_error
574     Thrown if the space dimension of \p os exceeds the maximum
575     allowed space dimension.
576   */
577   template <typename U>
578   explicit BD_Shape(const Octagonal_Shape<U>& os,
579                     Complexity_Class complexity = ANY_COMPLEXITY);
580 
581   /*! \brief
582     The assignment operator
583     (\p *this and \p y can be dimension-incompatible).
584   */
585   BD_Shape& operator=(const BD_Shape& y);
586 
587   /*! \brief
588     Swaps \p *this with \p y
589     (\p *this and \p y can be dimension-incompatible).
590   */
591   void m_swap(BD_Shape& y);
592 
593   //! Destructor.
594   ~BD_Shape();
595 
596   //@} Constructors, Assignment, Swap and Destructor
597 
598   //! \name Member Functions that Do Not Modify the BD_Shape
599   //@{
600 
601   //! Returns the dimension of the vector space enclosing \p *this.
602   dimension_type space_dimension() const;
603 
604   /*! \brief
605     Returns \f$0\f$, if \p *this is empty; otherwise, returns the
606     \ref Affine_Independence_and_Affine_Dimension "affine dimension"
607     of \p *this.
608   */
609   dimension_type affine_dimension() const;
610 
611   //! Returns a system of constraints defining \p *this.
612   Constraint_System constraints() const;
613 
614   //! Returns a minimized system of constraints defining \p *this.
615   Constraint_System minimized_constraints() const;
616 
617   //! Returns a system of (equality) congruences satisfied by \p *this.
618   Congruence_System congruences() const;
619 
620   /*! \brief
621     Returns a minimal system of (equality) congruences
622     satisfied by \p *this with the same affine dimension as \p *this.
623   */
624   Congruence_System minimized_congruences() const;
625 
626   /*! \brief
627     Returns <CODE>true</CODE> if and only if \p expr is
628     bounded from above in \p *this.
629 
630     \exception std::invalid_argument
631     Thrown if \p expr and \p *this are dimension-incompatible.
632   */
633   bool bounds_from_above(const Linear_Expression& expr) const;
634 
635   /*! \brief
636     Returns <CODE>true</CODE> if and only if \p expr is
637     bounded from below in \p *this.
638 
639     \exception std::invalid_argument
640     Thrown if \p expr and \p *this are dimension-incompatible.
641   */
642   bool bounds_from_below(const Linear_Expression& expr) const;
643 
644   /*! \brief
645     Returns <CODE>true</CODE> if and only if \p *this is not empty
646     and \p expr is bounded from above in \p *this, in which case
647     the supremum value is computed.
648 
649     \param expr
650     The linear expression to be maximized subject to \p *this;
651 
652     \param sup_n
653     The numerator of the supremum value;
654 
655     \param sup_d
656     The denominator of the supremum value;
657 
658     \param maximum
659     <CODE>true</CODE> if and only if the supremum is also the maximum value.
660 
661     \exception std::invalid_argument
662     Thrown if \p expr and \p *this are dimension-incompatible.
663 
664     If \p *this is empty or \p expr is not bounded from above,
665     <CODE>false</CODE> is returned and \p sup_n, \p sup_d
666     and \p maximum are left untouched.
667   */
668   bool maximize(const Linear_Expression& expr,
669                 Coefficient& sup_n, Coefficient& sup_d, bool& maximum) const;
670 
671   /*! \brief
672     Returns <CODE>true</CODE> if and only if \p *this is not empty
673     and \p expr is bounded from above in \p *this, in which case
674     the supremum value and a point where \p expr reaches it are computed.
675 
676     \param expr
677     The linear expression to be maximized subject to \p *this;
678 
679     \param sup_n
680     The numerator of the supremum value;
681 
682     \param sup_d
683     The denominator of the supremum value;
684 
685     \param maximum
686     <CODE>true</CODE> if and only if the supremum is also the maximum value;
687 
688     \param g
689     When maximization succeeds, will be assigned the point or
690     closure point where \p expr reaches its supremum value.
691 
692     \exception std::invalid_argument
693     Thrown if \p expr and \p *this are dimension-incompatible.
694 
695     If \p *this is empty or \p expr is not bounded from above,
696     <CODE>false</CODE> is returned and \p sup_n, \p sup_d, \p maximum
697     and \p g are left untouched.
698   */
699   bool maximize(const Linear_Expression& expr,
700                 Coefficient& sup_n, Coefficient& sup_d, bool& maximum,
701                 Generator& g) const;
702 
703   /*! \brief
704     Returns <CODE>true</CODE> if and only if \p *this is not empty
705     and \p expr is bounded from below in \p *this, in which case
706     the infimum value is computed.
707 
708     \param expr
709     The linear expression to be minimized subject to \p *this;
710 
711     \param inf_n
712     The numerator of the infimum value;
713 
714     \param inf_d
715     The denominator of the infimum value;
716 
717     \param minimum
718     <CODE>true</CODE> if and only if the infimum is also the minimum value.
719 
720     \exception std::invalid_argument
721     Thrown if \p expr and \p *this are dimension-incompatible.
722 
723     If \p *this is empty or \p expr is not bounded from below,
724     <CODE>false</CODE> is returned and \p inf_n, \p inf_d
725     and \p minimum are left untouched.
726   */
727   bool minimize(const Linear_Expression& expr,
728                 Coefficient& inf_n, Coefficient& inf_d, bool& minimum) const;
729 
730   /*! \brief
731     Returns <CODE>true</CODE> if and only if \p *this is not empty
732     and \p expr is bounded from below in \p *this, in which case
733     the infimum value and a point where \p expr reaches it are computed.
734 
735     \param expr
736     The linear expression to be minimized subject to \p *this;
737 
738     \param inf_n
739     The numerator of the infimum value;
740 
741     \param inf_d
742     The denominator of the infimum value;
743 
744     \param minimum
745     <CODE>true</CODE> if and only if the infimum is also the minimum value;
746 
747     \param g
748     When minimization succeeds, will be assigned a point or
749     closure point where \p expr reaches its infimum value.
750 
751     \exception std::invalid_argument
752     Thrown if \p expr and \p *this are dimension-incompatible.
753 
754     If \p *this is empty or \p expr is not bounded from below,
755     <CODE>false</CODE> is returned and \p inf_n, \p inf_d, \p minimum
756     and \p g are left untouched.
757   */
758   bool minimize(const Linear_Expression& expr,
759                 Coefficient& inf_n, Coefficient& inf_d, bool& minimum,
760                 Generator& g) const;
761 
762   /*! \brief
763     Returns <CODE>true</CODE> if and only if there exist a
764     unique value \p val such that \p *this
765     saturates the equality <CODE>expr = val</CODE>.
766 
767     \param expr
768     The linear expression for which the frequency is needed;
769 
770     \param freq_n
771     If <CODE>true</CODE> is returned, the value is set to \f$0\f$;
772     Present for interface compatibility with class Grid, where
773     the \ref Grid_Frequency "frequency" can have a non-zero value;
774 
775     \param freq_d
776     If <CODE>true</CODE> is returned, the value is set to \f$1\f$;
777 
778     \param val_n
779     The numerator of \p val;
780 
781     \param val_d
782     The denominator of \p val;
783 
784     \exception std::invalid_argument
785     Thrown if \p expr and \p *this are dimension-incompatible.
786 
787     If <CODE>false</CODE> is returned, then \p freq_n, \p freq_d,
788     \p val_n and \p val_d are left untouched.
789   */
790   bool frequency(const Linear_Expression& expr,
791                  Coefficient& freq_n, Coefficient& freq_d,
792                  Coefficient& val_n, Coefficient& val_d) const;
793 
794   //! Returns <CODE>true</CODE> if and only if \p *this contains \p y.
795   /*!
796     \exception std::invalid_argument
797     Thrown if \p *this and \p y are dimension-incompatible.
798   */
799   bool contains(const BD_Shape& y) const;
800 
801   //! Returns <CODE>true</CODE> if and only if \p *this strictly contains \p y.
802   /*!
803     \exception std::invalid_argument
804     Thrown if \p *this and \p y are dimension-incompatible.
805   */
806   bool strictly_contains(const BD_Shape& y) const;
807 
808   //! Returns <CODE>true</CODE> if and only if \p *this and \p y are disjoint.
809   /*!
810     \exception std::invalid_argument
811     Thrown if \p x and \p y are topology-incompatible or
812     dimension-incompatible.
813   */
814   bool is_disjoint_from(const BD_Shape& y) const;
815 
816   //! Returns the relations holding between \p *this and the constraint \p c.
817   /*!
818     \exception std::invalid_argument
819     Thrown if \p *this and constraint \p c are dimension-incompatible.
820   */
821   Poly_Con_Relation relation_with(const Constraint& c) const;
822 
823   //! Returns the relations holding between \p *this and the congruence \p cg.
824   /*!
825     \exception std::invalid_argument
826     Thrown if \p *this and congruence \p cg are dimension-incompatible.
827   */
828   Poly_Con_Relation relation_with(const Congruence& cg) const;
829 
830   //! Returns the relations holding between \p *this and the generator \p g.
831   /*!
832     \exception std::invalid_argument
833     Thrown if \p *this and generator \p g are dimension-incompatible.
834   */
835   Poly_Gen_Relation relation_with(const Generator& g) const;
836 
837   //! Returns <CODE>true</CODE> if and only if \p *this is an empty BDS.
838   bool is_empty() const;
839 
840   //! Returns <CODE>true</CODE> if and only if \p *this is a universe BDS.
841   bool is_universe() const;
842 
843   //! Returns <CODE>true</CODE> if and only if \p *this is discrete.
844   bool is_discrete() const;
845 
846   /*! \brief
847     Returns <CODE>true</CODE> if and only if \p *this
848     is a topologically closed subset of the vector space.
849   */
850   bool is_topologically_closed() const;
851 
852   //! Returns <CODE>true</CODE> if and only if \p *this is a bounded BDS.
853   bool is_bounded() const;
854 
855   /*! \brief
856     Returns <CODE>true</CODE> if and only if \p *this
857     contains at least one integer point.
858   */
859   bool contains_integer_point() const;
860 
861   /*! \brief
862     Returns <CODE>true</CODE> if and only if \p var is constrained in
863     \p *this.
864 
865     \exception std::invalid_argument
866     Thrown if \p var is not a space dimension of \p *this.
867   */
868   bool constrains(Variable var) const;
869 
870   /*! \brief
871     Returns <CODE>true</CODE> if and only if \p *this satisfies
872     all its invariants.
873   */
874   bool OK() const;
875 
876   //@} Member Functions that Do Not Modify the BD_Shape
877 
878   //! \name Space-Dimension Preserving Member Functions that May Modify the BD_Shape
879   //@{
880 
881   /*! \brief
882     Adds a copy of constraint \p c to the system of bounded differences
883     defining \p *this.
884 
885     \param c
886     The constraint to be added.
887 
888     \exception std::invalid_argument
889     Thrown if \p *this and constraint \p c are dimension-incompatible,
890     or \p c is not optimally supported by the BD shape domain.
891   */
892   void add_constraint(const Constraint& c);
893 
894   /*! \brief
895     Adds a copy of congruence \p cg to the system of congruences of \p *this.
896 
897     \param cg
898     The congruence to be added.
899 
900     \exception std::invalid_argument
901     Thrown if \p *this and congruence \p cg are dimension-incompatible,
902     or \p cg is not optimally supported by the BD shape domain.
903   */
904   void add_congruence(const Congruence& cg);
905 
906   /*! \brief
907     Adds the constraints in \p cs to the system of bounded differences
908     defining \p *this.
909 
910     \param  cs
911     The constraints that will be added.
912 
913     \exception std::invalid_argument
914     Thrown if \p *this and \p cs are dimension-incompatible,
915     or \p cs contains a constraint which is not optimally supported
916     by the BD shape domain.
917   */
918   void add_constraints(const Constraint_System& cs);
919 
920   /*! \brief
921     Adds the constraints in \p cs to the system of constraints
922     of \p *this.
923 
924     \param cs
925     The constraint system to be added to \p *this.  The constraints in
926     \p cs may be recycled.
927 
928     \exception std::invalid_argument
929     Thrown if \p *this and \p cs are dimension-incompatible,
930     or \p cs contains a constraint which is not optimally supported
931     by the BD shape domain.
932 
933     \warning
934     The only assumption that can be made on \p cs upon successful or
935     exceptional return is that it can be safely destroyed.
936   */
937   void add_recycled_constraints(Constraint_System& cs);
938 
939   /*! \brief
940     Adds to \p *this constraints equivalent to the congruences in \p cgs.
941 
942     \param cgs
943     Contains the congruences that will be added to the system of
944     constraints of \p *this.
945 
946     \exception std::invalid_argument
947     Thrown if \p *this and \p cgs are dimension-incompatible,
948     or \p cgs contains a congruence which is not optimally supported
949     by the BD shape domain.
950   */
951   void add_congruences(const Congruence_System& cgs);
952 
953   /*! \brief
954     Adds to \p *this constraints equivalent to the congruences in \p cgs.
955 
956     \param cgs
957     Contains the congruences that will be added to the system of
958     constraints of \p *this. Its elements may be recycled.
959 
960     \exception std::invalid_argument
961     Thrown if \p *this and \p cgs are dimension-incompatible,
962     or \p cgs contains a congruence which is not optimally supported
963     by the BD shape domain.
964 
965     \warning
966     The only assumption that can be made on \p cgs upon successful or
967     exceptional return is that it can be safely destroyed.
968   */
969   void add_recycled_congruences(Congruence_System& cgs);
970 
971   /*! \brief
972     Uses a copy of constraint \p c to refine the system of bounded differences
973     defining \p *this.
974 
975     \param c
976     The constraint. If it is not a bounded difference, it will be ignored.
977 
978     \exception std::invalid_argument
979     Thrown if \p *this and constraint \p c are dimension-incompatible.
980   */
981   void refine_with_constraint(const Constraint& c);
982 
983   /*! \brief
984     Uses a copy of congruence \p cg to refine the system of
985     bounded differences  of \p *this.
986 
987     \param cg
988     The congruence. If it is not a bounded difference equality, it
989     will be ignored.
990 
991     \exception std::invalid_argument
992     Thrown if \p *this and congruence \p cg are dimension-incompatible.
993   */
994   void refine_with_congruence(const Congruence& cg);
995 
996   /*! \brief
997     Uses a copy of the constraints in \p cs to refine the system of
998     bounded differences defining \p *this.
999 
1000     \param  cs
1001     The constraint system to be used. Constraints that are not bounded
1002     differences are ignored.
1003 
1004     \exception std::invalid_argument
1005     Thrown if \p *this and \p cs are dimension-incompatible.
1006   */
1007   void refine_with_constraints(const Constraint_System& cs);
1008 
1009   /*! \brief
1010     Uses a copy of the congruences in \p cgs to refine the system of
1011     bounded differences defining \p *this.
1012 
1013     \param  cgs
1014     The congruence system to be used. Congruences that are not bounded
1015     difference equalities are ignored.
1016 
1017     \exception std::invalid_argument
1018     Thrown if \p *this and \p cgs are dimension-incompatible.
1019   */
1020   void refine_with_congruences(const Congruence_System& cgs);
1021 
1022   /*! \brief
1023     Refines the system of BD_Shape constraints defining \p *this using
1024     the constraint expressed by \p left \f$\leq\f$ \p right.
1025 
1026     \param left
1027     The linear form on intervals with floating point boundaries that
1028     is at the left of the comparison operator. All of its coefficients
1029     MUST be bounded.
1030 
1031     \param right
1032     The linear form on intervals with floating point boundaries that
1033     is at the right of the comparison operator. All of its coefficients
1034     MUST be bounded.
1035 
1036     \exception std::invalid_argument
1037     Thrown if \p left (or \p right) is dimension-incompatible with \p *this.
1038 
1039     This function is used in abstract interpretation to model a filter
1040     that is generated by a comparison of two expressions that are correctly
1041     approximated by \p left and \p right respectively.
1042   */
1043   template <typename Interval_Info>
1044   void refine_with_linear_form_inequality(
1045                    const Linear_Form<Interval<T, Interval_Info> >& left,
1046                    const Linear_Form<Interval<T, Interval_Info> >& right);
1047 
1048   /*! \brief
1049     Refines the system of BD_Shape constraints defining \p *this using
1050     the constraint expressed by \p left \f$\relsym\f$ \p right, where
1051     \f$\relsym\f$ is the relation symbol specified by \p relsym.
1052 
1053     \param left
1054     The linear form on intervals with floating point boundaries that
1055     is at the left of the comparison operator. All of its coefficients
1056     MUST be bounded.
1057 
1058     \param right
1059     The linear form on intervals with floating point boundaries that
1060     is at the right of the comparison operator. All of its coefficients
1061     MUST be bounded.
1062 
1063     \param relsym
1064     The relation symbol.
1065 
1066     \exception std::invalid_argument
1067     Thrown if \p left (or \p right) is dimension-incompatible with \p *this.
1068 
1069     \exception std::runtime_error
1070     Thrown if \p relsym is not a valid relation symbol.
1071 
1072     This function is used in abstract interpretation to model a filter
1073     that is generated by a comparison of two expressions that are correctly
1074     approximated by \p left and \p right respectively.
1075   */
1076   template <typename Interval_Info>
1077   void generalized_refine_with_linear_form_inequality(
1078                    const Linear_Form<Interval<T, Interval_Info> >& left,
1079                    const Linear_Form<Interval<T, Interval_Info> >& right,
1080                    Relation_Symbol relsym);
1081 
1082   //! Applies to \p dest the interval constraints embedded in \p *this.
1083   /*!
1084     \param dest
1085     The object to which the constraints will be added.
1086 
1087     \exception std::invalid_argument
1088     Thrown if \p *this is dimension-incompatible with \p dest.
1089 
1090     The template type parameter U must provide the following methods.
1091     \code
1092       dimension_type space_dimension() const
1093     \endcode
1094     returns the space dimension of the object.
1095     \code
1096       void set_empty()
1097     \endcode
1098     sets the object to an empty object.
1099     \code
1100       bool restrict_lower(dimension_type dim, const T& lb)
1101     \endcode
1102     restricts the object by applying the lower bound \p lb to the space
1103     dimension \p dim and returns <CODE>false</CODE> if and only if the
1104     object becomes empty.
1105     \code
1106       bool restrict_upper(dimension_type dim, const T& ub)
1107     \endcode
1108     restricts the object by applying the upper bound \p ub to the space
1109     dimension \p dim and returns <CODE>false</CODE> if and only if the
1110     object becomes empty.
1111   */
1112   template <typename U>
1113   void export_interval_constraints(U& dest) const;
1114 
1115   /*! \brief
1116     Computes the \ref Cylindrification "cylindrification" of \p *this with
1117     respect to space dimension \p var, assigning the result to \p *this.
1118 
1119     \param var
1120     The space dimension that will be unconstrained.
1121 
1122     \exception std::invalid_argument
1123     Thrown if \p var is not a space dimension of \p *this.
1124   */
1125   void unconstrain(Variable var);
1126 
1127   /*! \brief
1128     Computes the \ref Cylindrification "cylindrification" of \p *this with
1129     respect to the set of space dimensions \p vars,
1130     assigning the result to \p *this.
1131 
1132     \param vars
1133     The set of space dimension that will be unconstrained.
1134 
1135     \exception std::invalid_argument
1136     Thrown if \p *this is dimension-incompatible with one of the
1137     Variable objects contained in \p vars.
1138   */
1139   void unconstrain(const Variables_Set& vars);
1140 
1141   //! Assigns to \p *this the intersection of \p *this and \p y.
1142   /*!
1143     \exception std::invalid_argument
1144     Thrown if \p *this and \p y are dimension-incompatible.
1145   */
1146   void intersection_assign(const BD_Shape& y);
1147 
1148   /*! \brief
1149     Assigns to \p *this the smallest BDS containing the union
1150     of \p *this and \p y.
1151 
1152     \exception std::invalid_argument
1153     Thrown if \p *this and \p y are dimension-incompatible.
1154   */
1155   void upper_bound_assign(const BD_Shape& y);
1156 
1157   /*! \brief
1158     If the upper bound of \p *this and \p y is exact, it is assigned
1159     to \p *this and <CODE>true</CODE> is returned,
1160     otherwise <CODE>false</CODE> is returned.
1161 
1162     \exception std::invalid_argument
1163     Thrown if \p *this and \p y are dimension-incompatible.
1164   */
1165   bool upper_bound_assign_if_exact(const BD_Shape& y);
1166 
1167   /*! \brief
1168     If the \e integer upper bound of \p *this and \p y is exact,
1169     it is assigned to \p *this and <CODE>true</CODE> is returned;
1170     otherwise <CODE>false</CODE> is returned.
1171 
1172     \exception std::invalid_argument
1173     Thrown if \p *this and \p y are dimension-incompatible.
1174 
1175     \note
1176     The integer upper bound of two rational BDS is the smallest rational
1177     BDS containing all the integral points of the two arguments.
1178     This method requires that the coefficient type parameter \c T is
1179     an integral type.
1180   */
1181   bool integer_upper_bound_assign_if_exact(const BD_Shape& y);
1182 
1183   /*! \brief
1184     Assigns to \p *this the smallest BD shape containing
1185     the set difference of \p *this and \p y.
1186 
1187     \exception std::invalid_argument
1188     Thrown if \p *this and \p y are dimension-incompatible.
1189   */
1190   void difference_assign(const BD_Shape& y);
1191 
1192   /*! \brief
1193     Assigns to \p *this a \ref Meet_Preserving_Simplification
1194     "meet-preserving simplification" of \p *this with respect to \p y.
1195     If \c false is returned, then the intersection is empty.
1196 
1197     \exception std::invalid_argument
1198     Thrown if \p *this and \p y are topology-incompatible or
1199     dimension-incompatible.
1200   */
1201   bool simplify_using_context_assign(const BD_Shape& y);
1202 
1203   /*! \brief
1204     Assigns to \p *this the
1205     \ref Single_Update_Affine_Functions "affine image"
1206     of \p *this under the function mapping variable \p var into the
1207     affine expression specified by \p expr and \p denominator.
1208 
1209     \param var
1210     The variable to which the affine expression is assigned.
1211 
1212     \param expr
1213     The numerator of the affine expression.
1214 
1215     \param denominator
1216     The denominator of the affine expression.
1217 
1218     \exception std::invalid_argument
1219     Thrown if \p denominator is zero or if \p expr and \p *this
1220     are dimension-incompatible or if \p var is not a dimension of \p *this.
1221   */
1222   void affine_image(Variable var,
1223                     const Linear_Expression& expr,
1224                     Coefficient_traits::const_reference denominator
1225                     = Coefficient_one());
1226 
1227   // FIXME: To be completed.
1228   /*! \brief
1229     Assigns to \p *this the \ref affine_form_relation "affine form image"
1230     of \p *this under the function mapping variable \p var into the
1231     affine expression(s) specified by \p lf.
1232 
1233     \param var
1234     The variable to which the affine expression is assigned.
1235 
1236     \param lf
1237     The linear form on intervals with floating point coefficients that
1238     defines the affine expression. ALL of its coefficients MUST be bounded.
1239 
1240     \exception std::invalid_argument
1241     Thrown if \p lf and \p *this are dimension-incompatible or if \p var
1242     is not a dimension of \p *this.
1243   */
1244   template <typename Interval_Info>
1245   void affine_form_image(Variable var,
1246                         const Linear_Form<Interval<T, Interval_Info> >& lf);
1247 
1248   /*! \brief
1249     Assigns to \p *this the
1250     \ref Single_Update_Affine_Functions "affine preimage"
1251     of \p *this under the function mapping variable \p var into the
1252     affine expression specified by \p expr and \p denominator.
1253 
1254     \param var
1255     The variable to which the affine expression is substituted.
1256 
1257     \param expr
1258     The numerator of the affine expression.
1259 
1260     \param denominator
1261     The denominator of the affine expression.
1262 
1263     \exception std::invalid_argument
1264     Thrown if \p denominator is zero or if \p expr and \p *this
1265     are dimension-incompatible or if \p var is not a dimension of \p *this.
1266   */
1267   void affine_preimage(Variable var,
1268                        const Linear_Expression& expr,
1269                        Coefficient_traits::const_reference denominator
1270                        = Coefficient_one());
1271 
1272   /*! \brief
1273     Assigns to \p *this the image of \p *this with respect to the
1274     \ref Generalized_Affine_Relations "affine relation"
1275     \f$\mathrm{var}' \relsym \frac{\mathrm{expr}}{\mathrm{denominator}}\f$,
1276     where \f$\mathord{\relsym}\f$ is the relation symbol encoded
1277     by \p relsym.
1278 
1279     \param var
1280     The left hand side variable of the generalized affine transfer function.
1281 
1282     \param relsym
1283     The relation symbol.
1284 
1285     \param expr
1286     The numerator of the right hand side affine expression.
1287 
1288     \param denominator
1289     The denominator of the right hand side affine expression.
1290 
1291     \exception std::invalid_argument
1292     Thrown if \p denominator is zero or if \p expr and \p *this
1293     are dimension-incompatible or if \p var is not a dimension
1294     of \p *this or if \p relsym is a strict relation symbol.
1295   */
1296   void generalized_affine_image(Variable var,
1297                                 Relation_Symbol relsym,
1298                                 const Linear_Expression& expr,
1299                                 Coefficient_traits::const_reference denominator
1300                                 = Coefficient_one());
1301 
1302   /*! \brief
1303     Assigns to \p *this the image of \p *this with respect to the
1304     \ref Generalized_Affine_Relations "affine relation"
1305     \f$\mathrm{lhs}' \relsym \mathrm{rhs}\f$, where
1306     \f$\mathord{\relsym}\f$ is the relation symbol encoded by \p relsym.
1307 
1308     \param lhs
1309     The left hand side affine expression.
1310 
1311     \param relsym
1312     The relation symbol.
1313 
1314     \param rhs
1315     The right hand side affine expression.
1316 
1317     \exception std::invalid_argument
1318     Thrown if \p *this is dimension-incompatible with \p lhs or \p rhs
1319     or if \p relsym is a strict relation symbol.
1320   */
1321   void generalized_affine_image(const Linear_Expression& lhs,
1322                                 Relation_Symbol relsym,
1323                                 const Linear_Expression& rhs);
1324 
1325   /*! \brief
1326     Assigns to \p *this the preimage of \p *this with respect to the
1327     \ref Generalized_Affine_Relations "affine relation"
1328     \f$\mathrm{var}' \relsym \frac{\mathrm{expr}}{\mathrm{denominator}}\f$,
1329     where \f$\mathord{\relsym}\f$ is the relation symbol encoded
1330     by \p relsym.
1331 
1332     \param var
1333     The left hand side variable of the generalized affine transfer function.
1334 
1335     \param relsym
1336     The relation symbol.
1337 
1338     \param expr
1339     The numerator of the right hand side affine expression.
1340 
1341     \param denominator
1342     The denominator of the right hand side affine expression.
1343 
1344     \exception std::invalid_argument
1345     Thrown if \p denominator is zero or if \p expr and \p *this
1346     are dimension-incompatible or if \p var is not a dimension
1347     of \p *this or if \p relsym is a strict relation symbol.
1348   */
1349   void generalized_affine_preimage(Variable var,
1350                                    Relation_Symbol relsym,
1351                                    const Linear_Expression& expr,
1352                                    Coefficient_traits::const_reference
1353                                    denominator = Coefficient_one());
1354 
1355   /*! \brief
1356     Assigns to \p *this the preimage of \p *this with respect to the
1357     \ref Generalized_Affine_Relations "affine relation"
1358     \f$\mathrm{lhs}' \relsym \mathrm{rhs}\f$, where
1359     \f$\mathord{\relsym}\f$ is the relation symbol encoded by \p relsym.
1360 
1361     \param lhs
1362     The left hand side affine expression.
1363 
1364     \param relsym
1365     The relation symbol.
1366 
1367     \param rhs
1368     The right hand side affine expression.
1369 
1370     \exception std::invalid_argument
1371     Thrown if \p *this is dimension-incompatible with \p lhs or \p rhs
1372     or if \p relsym is a strict relation symbol.
1373   */
1374   void generalized_affine_preimage(const Linear_Expression& lhs,
1375                                    Relation_Symbol relsym,
1376                                    const Linear_Expression& rhs);
1377 
1378   /*!
1379     \brief
1380     Assigns to \p *this the image of \p *this with respect to the
1381     \ref Single_Update_Bounded_Affine_Relations "bounded affine relation"
1382     \f$\frac{\mathrm{lb\_expr}}{\mathrm{denominator}}
1383          \leq \mathrm{var}'
1384            \leq \frac{\mathrm{ub\_expr}}{\mathrm{denominator}}\f$.
1385 
1386     \param var
1387     The variable updated by the affine relation;
1388 
1389     \param lb_expr
1390     The numerator of the lower bounding affine expression;
1391 
1392     \param ub_expr
1393     The numerator of the upper bounding affine expression;
1394 
1395     \param denominator
1396     The (common) denominator for the lower and upper bounding
1397     affine expressions (optional argument with default value 1).
1398 
1399     \exception std::invalid_argument
1400     Thrown if \p denominator is zero or if \p lb_expr (resp., \p ub_expr)
1401     and \p *this are dimension-incompatible or if \p var is not a space
1402     dimension of \p *this.
1403   */
1404   void bounded_affine_image(Variable var,
1405                             const Linear_Expression& lb_expr,
1406                             const Linear_Expression& ub_expr,
1407                             Coefficient_traits::const_reference denominator
1408                             = Coefficient_one());
1409 
1410   /*!
1411     \brief
1412     Assigns to \p *this the preimage of \p *this with respect to the
1413     \ref Single_Update_Bounded_Affine_Relations "bounded affine relation"
1414     \f$\frac{\mathrm{lb\_expr}}{\mathrm{denominator}}
1415          \leq \mathrm{var}'
1416            \leq \frac{\mathrm{ub\_expr}}{\mathrm{denominator}}\f$.
1417 
1418     \param var
1419     The variable updated by the affine relation;
1420 
1421     \param lb_expr
1422     The numerator of the lower bounding affine expression;
1423 
1424     \param ub_expr
1425     The numerator of the upper bounding affine expression;
1426 
1427     \param denominator
1428     The (common) denominator for the lower and upper bounding
1429     affine expressions (optional argument with default value 1).
1430 
1431     \exception std::invalid_argument
1432     Thrown if \p denominator is zero or if \p lb_expr (resp., \p ub_expr)
1433     and \p *this are dimension-incompatible or if \p var is not a space
1434     dimension of \p *this.
1435   */
1436   void bounded_affine_preimage(Variable var,
1437                                const Linear_Expression& lb_expr,
1438                                const Linear_Expression& ub_expr,
1439                                Coefficient_traits::const_reference denominator
1440                                = Coefficient_one());
1441   /*! \brief
1442     Assigns to \p *this the result of computing the
1443     \ref Time_Elapse_Operator "time-elapse" between \p *this and \p y.
1444 
1445     \exception std::invalid_argument
1446     Thrown if \p *this and \p y are dimension-incompatible.
1447   */
1448   void time_elapse_assign(const BD_Shape& y);
1449 
1450   /*! \brief
1451     \ref Wrapping_Operator "Wraps" the specified dimensions of the
1452     vector space.
1453 
1454     \param vars
1455     The set of Variable objects corresponding to the space dimensions
1456     to be wrapped.
1457 
1458     \param w
1459     The width of the bounded integer type corresponding to
1460     all the dimensions to be wrapped.
1461 
1462     \param r
1463     The representation of the bounded integer type corresponding to
1464     all the dimensions to be wrapped.
1465 
1466     \param o
1467     The overflow behavior of the bounded integer type corresponding to
1468     all the dimensions to be wrapped.
1469 
1470     \param cs_p
1471     Possibly null pointer to a constraint system whose variables
1472     are contained in \p vars.  If <CODE>*cs_p</CODE> depends on
1473     variables not in \p vars, the behavior is undefined.
1474     When non-null, the pointed-to constraint system is assumed to
1475     represent the conditional or looping construct guard with respect
1476     to which wrapping is performed.  Since wrapping requires the
1477     computation of upper bounds and due to non-distributivity of
1478     constraint refinement over upper bounds, passing a constraint
1479     system in this way can be more precise than refining the result of
1480     the wrapping operation with the constraints in <CODE>*cs_p</CODE>.
1481 
1482     \param complexity_threshold
1483     A precision parameter of the \ref Wrapping_Operator "wrapping operator":
1484     higher values result in possibly improved precision.
1485 
1486     \param wrap_individually
1487     <CODE>true</CODE> if the dimensions should be wrapped individually
1488     (something that results in much greater efficiency to the detriment of
1489     precision).
1490 
1491     \exception std::invalid_argument
1492     Thrown if <CODE>*cs_p</CODE> is dimension-incompatible with
1493     \p vars, or if \p *this is dimension-incompatible \p vars or with
1494     <CODE>*cs_p</CODE>.
1495   */
1496   void wrap_assign(const Variables_Set& vars,
1497                    Bounded_Integer_Type_Width w,
1498                    Bounded_Integer_Type_Representation r,
1499                    Bounded_Integer_Type_Overflow o,
1500                    const Constraint_System* cs_p = 0,
1501                    unsigned complexity_threshold = 16,
1502                    bool wrap_individually = true);
1503 
1504   /*! \brief
1505     Possibly tightens \p *this by dropping some points with non-integer
1506     coordinates.
1507 
1508     \param complexity
1509     The maximal complexity of any algorithms used.
1510 
1511     \note
1512     Currently there is no optimality guarantee, not even if
1513     \p complexity is <CODE>ANY_COMPLEXITY</CODE>.
1514   */
1515   void drop_some_non_integer_points(Complexity_Class complexity
1516                                     = ANY_COMPLEXITY);
1517 
1518   /*! \brief
1519     Possibly tightens \p *this by dropping some points with non-integer
1520     coordinates for the space dimensions corresponding to \p vars.
1521 
1522     \param vars
1523     Points with non-integer coordinates for these variables/space-dimensions
1524     can be discarded.
1525 
1526     \param complexity
1527     The maximal complexity of any algorithms used.
1528 
1529     \note
1530     Currently there is no optimality guarantee, not even if
1531     \p complexity is <CODE>ANY_COMPLEXITY</CODE>.
1532   */
1533   void drop_some_non_integer_points(const Variables_Set& vars,
1534                                     Complexity_Class complexity
1535                                     = ANY_COMPLEXITY);
1536 
1537   //! Assigns to \p *this its topological closure.
1538   void topological_closure_assign();
1539 
1540   /*! \brief
1541     Assigns to \p *this the result of computing the
1542     \ref CC76_extrapolation "CC76-extrapolation" between \p *this and \p y.
1543 
1544     \param y
1545     A BDS that <EM>must</EM> be contained in \p *this.
1546 
1547     \param tp
1548     An optional pointer to an unsigned variable storing the number of
1549     available tokens (to be used when applying the
1550     \ref Widening_with_Tokens "widening with tokens" delay technique).
1551 
1552     \exception std::invalid_argument
1553     Thrown if \p *this and \p y are dimension-incompatible.
1554   */
1555   void CC76_extrapolation_assign(const BD_Shape& y, unsigned* tp = 0);
1556 
1557   /*! \brief
1558     Assigns to \p *this the result of computing the
1559     \ref CC76_extrapolation "CC76-extrapolation" between \p *this and \p y.
1560 
1561     \param y
1562     A BDS that <EM>must</EM> be contained in \p *this.
1563 
1564     \param first
1565     An iterator referencing the first stop-point.
1566 
1567     \param last
1568     An iterator referencing one past the last stop-point.
1569 
1570     \param tp
1571     An optional pointer to an unsigned variable storing the number of
1572     available tokens (to be used when applying the
1573     \ref Widening_with_Tokens "widening with tokens" delay technique).
1574 
1575     \exception std::invalid_argument
1576     Thrown if \p *this and \p y are dimension-incompatible.
1577   */
1578   template <typename Iterator>
1579   void CC76_extrapolation_assign(const BD_Shape& y,
1580                                  Iterator first, Iterator last,
1581                                  unsigned* tp = 0);
1582 
1583   /*! \brief
1584     Assigns to \p *this the result of computing the
1585     \ref BHMZ05_widening "BHMZ05-widening" of \p *this and \p y.
1586 
1587     \param y
1588     A BDS that <EM>must</EM> be contained in \p *this.
1589 
1590     \param tp
1591     An optional pointer to an unsigned variable storing the number of
1592     available tokens (to be used when applying the
1593     \ref Widening_with_Tokens "widening with tokens" delay technique).
1594 
1595     \exception std::invalid_argument
1596     Thrown if \p *this and \p y are dimension-incompatible.
1597   */
1598   void BHMZ05_widening_assign(const BD_Shape& y, unsigned* tp = 0);
1599 
1600   /*! \brief
1601     Improves the result of the \ref BHMZ05_widening "BHMZ05-widening"
1602     computation by also enforcing those constraints in \p cs that are
1603     satisfied by all the points of \p *this.
1604 
1605     \param y
1606     A BDS that <EM>must</EM> be contained in \p *this.
1607 
1608     \param cs
1609     The system of constraints used to improve the widened BDS.
1610 
1611     \param tp
1612     An optional pointer to an unsigned variable storing the number of
1613     available tokens (to be used when applying the
1614     \ref Widening_with_Tokens "widening with tokens" delay technique).
1615 
1616     \exception std::invalid_argument
1617     Thrown if \p *this, \p y and \p cs are dimension-incompatible or
1618     if \p cs contains a strict inequality.
1619   */
1620   void limited_BHMZ05_extrapolation_assign(const BD_Shape& y,
1621                                            const Constraint_System& cs,
1622                                            unsigned* tp = 0);
1623 
1624   /*! \brief
1625     Assigns to \p *this the result of restoring in \p y the constraints
1626     of \p *this that were lost by
1627     \ref CC76_extrapolation "CC76-extrapolation" applications.
1628 
1629     \param y
1630     A BDS that <EM>must</EM> contain \p *this.
1631 
1632     \exception std::invalid_argument
1633     Thrown if \p *this and \p y are dimension-incompatible.
1634 
1635     \note
1636     As was the case for widening operators, the argument \p y is meant to
1637     denote the value computed in the previous iteration step, whereas
1638     \p *this denotes the value computed in the current iteration step
1639     (in the <EM>decreasing</EM> iteration sequence). Hence, the call
1640     <CODE>x.CC76_narrowing_assign(y)</CODE> will assign to \p x
1641     the result of the computation \f$\mathtt{y} \Delta \mathtt{x}\f$.
1642   */
1643   void CC76_narrowing_assign(const BD_Shape& y);
1644 
1645   /*! \brief
1646     Improves the result of the \ref CC76_extrapolation "CC76-extrapolation"
1647     computation by also enforcing those constraints in \p cs that are
1648     satisfied by all the points of \p *this.
1649 
1650     \param y
1651     A BDS that <EM>must</EM> be contained in \p *this.
1652 
1653     \param cs
1654     The system of constraints used to improve the widened BDS.
1655 
1656     \param tp
1657     An optional pointer to an unsigned variable storing the number of
1658     available tokens (to be used when applying the
1659     \ref Widening_with_Tokens "widening with tokens" delay technique).
1660 
1661     \exception std::invalid_argument
1662     Thrown if \p *this, \p y and \p cs are dimension-incompatible or
1663     if \p cs contains a strict inequality.
1664   */
1665   void limited_CC76_extrapolation_assign(const BD_Shape& y,
1666                                          const Constraint_System& cs,
1667                                          unsigned* tp = 0);
1668 
1669   /*! \brief
1670     Assigns to \p *this the result of computing the
1671     \ref H79_widening "H79-widening" between \p *this and \p y.
1672 
1673     \param y
1674     A BDS that <EM>must</EM> be contained in \p *this.
1675 
1676     \param tp
1677     An optional pointer to an unsigned variable storing the number of
1678     available tokens (to be used when applying the
1679     \ref Widening_with_Tokens "widening with tokens" delay technique).
1680 
1681     \exception std::invalid_argument
1682     Thrown if \p *this and \p y are dimension-incompatible.
1683   */
1684   void H79_widening_assign(const BD_Shape& y, unsigned* tp = 0);
1685 
1686   //! Same as H79_widening_assign(y, tp).
1687   void widening_assign(const BD_Shape& y, unsigned* tp = 0);
1688 
1689   /*! \brief
1690     Improves the result of the \ref H79_widening "H79-widening"
1691     computation by also enforcing those constraints in \p cs that are
1692     satisfied by all the points of \p *this.
1693 
1694     \param y
1695     A BDS that <EM>must</EM> be contained in \p *this.
1696 
1697     \param cs
1698     The system of constraints used to improve the widened BDS.
1699 
1700     \param tp
1701     An optional pointer to an unsigned variable storing the number of
1702     available tokens (to be used when applying the
1703     \ref Widening_with_Tokens "widening with tokens" delay technique).
1704 
1705     \exception std::invalid_argument
1706     Thrown if \p *this, \p y and \p cs are dimension-incompatible.
1707   */
1708   void limited_H79_extrapolation_assign(const BD_Shape& y,
1709                                         const Constraint_System& cs,
1710                                         unsigned* tp = 0);
1711 
1712   //@} Space-Dimension Preserving Member Functions that May Modify [...]
1713 
1714   //! \name Member Functions that May Modify the Dimension of the Vector Space
1715   //@{
1716 
1717   //! Adds \p m new dimensions and embeds the old BDS into the new space.
1718   /*!
1719     \param m
1720     The number of dimensions to add.
1721 
1722     The new dimensions will be those having the highest indexes in the new
1723     BDS, which is defined by a system of bounded differences in which the
1724     variables running through the new dimensions are unconstrained.
1725     For instance, when starting from the BDS \f$\cB \sseq \Rset^2\f$
1726     and adding a third dimension, the result will be the BDS
1727     \f[
1728       \bigl\{\,
1729         (x, y, z)^\transpose \in \Rset^3
1730       \bigm|
1731         (x, y)^\transpose \in \cB
1732       \,\bigr\}.
1733     \f]
1734   */
1735   void add_space_dimensions_and_embed(dimension_type m);
1736 
1737   /*! \brief
1738     Adds \p m new dimensions to the BDS and does not embed it in
1739     the new vector space.
1740 
1741     \param m
1742     The number of dimensions to add.
1743 
1744     The new dimensions will be those having the highest indexes in the
1745     new BDS, which is defined by a system of bounded differences in
1746     which the variables running through the new dimensions are all
1747     constrained to be equal to 0.
1748     For instance, when starting from the BDS \f$\cB \sseq \Rset^2\f$
1749     and adding a third dimension, the result will be the BDS
1750     \f[
1751       \bigl\{\,
1752         (x, y, 0)^\transpose \in \Rset^3
1753       \bigm|
1754         (x, y)^\transpose \in \cB
1755       \,\bigr\}.
1756     \f]
1757   */
1758   void add_space_dimensions_and_project(dimension_type m);
1759 
1760   /*! \brief
1761     Assigns to \p *this the \ref Concatenating_Polyhedra "concatenation"
1762     of \p *this and \p y, taken in this order.
1763 
1764     \exception std::length_error
1765     Thrown if the concatenation would cause the vector space
1766     to exceed dimension <CODE>max_space_dimension()</CODE>.
1767   */
1768   void concatenate_assign(const BD_Shape& y);
1769 
1770   //! Removes all the specified dimensions.
1771   /*!
1772     \param vars
1773     The set of Variable objects corresponding to the dimensions to be removed.
1774 
1775     \exception std::invalid_argument
1776     Thrown if \p *this is dimension-incompatible with one of the Variable
1777     objects contained in \p vars.
1778   */
1779   void remove_space_dimensions(const Variables_Set& vars);
1780 
1781   /*! \brief
1782     Removes the higher dimensions so that the resulting space
1783     will have dimension \p new_dimension.
1784 
1785     \exception std::invalid_argument
1786     Thrown if \p new_dimension is greater than the space dimension
1787     of \p *this.
1788   */
1789   void remove_higher_space_dimensions(dimension_type new_dimension);
1790 
1791   /*! \brief
1792     Remaps the dimensions of the vector space according to
1793     a \ref Mapping_the_Dimensions_of_the_Vector_Space "partial function".
1794 
1795     \param pfunc
1796     The partial function specifying the destiny of each dimension.
1797 
1798     The template type parameter Partial_Function must provide
1799     the following methods.
1800     \code
1801       bool has_empty_codomain() const
1802     \endcode
1803     returns <CODE>true</CODE> if and only if the represented partial
1804     function has an empty co-domain (i.e., it is always undefined).
1805     The <CODE>has_empty_codomain()</CODE> method will always be called
1806     before the methods below.  However, if
1807     <CODE>has_empty_codomain()</CODE> returns <CODE>true</CODE>, none
1808     of the functions below will be called.
1809     \code
1810       dimension_type max_in_codomain() const
1811     \endcode
1812     returns the maximum value that belongs to the co-domain
1813     of the partial function.
1814     \code
1815       bool maps(dimension_type i, dimension_type& j) const
1816     \endcode
1817     Let \f$f\f$ be the represented function and \f$k\f$ be the value
1818     of \p i.  If \f$f\f$ is defined in \f$k\f$, then \f$f(k)\f$ is
1819     assigned to \p j and <CODE>true</CODE> is returned.
1820     If \f$f\f$ is undefined in \f$k\f$, then <CODE>false</CODE> is
1821     returned.
1822 
1823     The result is undefined if \p pfunc does not encode a partial
1824     function with the properties described in the
1825     \ref Mapping_the_Dimensions_of_the_Vector_Space
1826     "specification of the mapping operator".
1827   */
1828   template <typename Partial_Function>
1829   void map_space_dimensions(const Partial_Function& pfunc);
1830 
1831   //! Creates \p m copies of the space dimension corresponding to \p var.
1832   /*!
1833     \param var
1834     The variable corresponding to the space dimension to be replicated;
1835 
1836     \param m
1837     The number of replicas to be created.
1838 
1839     \exception std::invalid_argument
1840     Thrown if \p var does not correspond to a dimension of the vector space.
1841 
1842     \exception std::length_error
1843     Thrown if adding \p m new space dimensions would cause the
1844     vector space to exceed dimension <CODE>max_space_dimension()</CODE>.
1845 
1846     If \p *this has space dimension \f$n\f$, with \f$n > 0\f$,
1847     and <CODE>var</CODE> has space dimension \f$k \leq n\f$,
1848     then the \f$k\f$-th space dimension is
1849     \ref expand_space_dimension "expanded" to \p m new space dimensions
1850     \f$n\f$, \f$n+1\f$, \f$\dots\f$, \f$n+m-1\f$.
1851   */
1852   void expand_space_dimension(Variable var, dimension_type m);
1853 
1854   //! Folds the space dimensions in \p vars into \p dest.
1855   /*!
1856     \param vars
1857     The set of Variable objects corresponding to the space dimensions
1858     to be folded;
1859 
1860     \param dest
1861     The variable corresponding to the space dimension that is the
1862     destination of the folding operation.
1863 
1864     \exception std::invalid_argument
1865     Thrown if \p *this is dimension-incompatible with \p dest or with
1866     one of the Variable objects contained in \p vars.
1867     Also thrown if \p dest is contained in \p vars.
1868 
1869     If \p *this has space dimension \f$n\f$, with \f$n > 0\f$,
1870     <CODE>dest</CODE> has space dimension \f$k \leq n\f$,
1871     \p vars is a set of variables whose maximum space dimension
1872     is also less than or equal to \f$n\f$, and \p dest is not a member
1873     of \p vars, then the space dimensions corresponding to
1874     variables in \p vars are \ref fold_space_dimensions "folded"
1875     into the \f$k\f$-th space dimension.
1876   */
1877   void fold_space_dimensions(const Variables_Set& vars, Variable dest);
1878 
1879   //! Refines \p store with the constraints defining \p *this.
1880   /*!
1881     \param store
1882     The interval floating point abstract store to refine.
1883   */
1884   template <typename Interval_Info>
1885   void refine_fp_interval_abstract_store(Box<Interval<T, Interval_Info> >&
1886                                          store) const;
1887 
1888 
1889   //@} // Member Functions that May Modify the Dimension of the Vector Space
1890 
1891   PPL_OUTPUT_DECLARATIONS
1892 
1893   /*! \brief
1894     Loads from \p s an ASCII representation (as produced by
1895     ascii_dump(std::ostream&) const) and sets \p *this accordingly.
1896     Returns <CODE>true</CODE> if successful, <CODE>false</CODE> otherwise.
1897   */
1898   bool ascii_load(std::istream& s);
1899 
1900   //! Returns the total size in bytes of the memory occupied by \p *this.
1901   memory_size_type total_memory_in_bytes() const;
1902 
1903   //! Returns the size in bytes of the memory managed by \p *this.
1904   memory_size_type external_memory_in_bytes() const;
1905 
1906   /*! \brief
1907     Returns a 32-bit hash code for \p *this.
1908 
1909     If \p x and \p y are such that <CODE>x == y</CODE>,
1910     then <CODE>x.hash_code() == y.hash_code()</CODE>.
1911   */
1912   int32_t hash_code() const;
1913 
1914   friend bool operator==<T>(const BD_Shape<T>& x, const BD_Shape<T>& y);
1915 
1916   template <typename Temp, typename To, typename U>
1917   friend bool Parma_Polyhedra_Library
1918   ::rectilinear_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
1919                                 const BD_Shape<U>& x, const BD_Shape<U>& y,
1920                                 const Rounding_Dir dir,
1921                                 Temp& tmp0, Temp& tmp1, Temp& tmp2);
1922   template <typename Temp, typename To, typename U>
1923   friend bool Parma_Polyhedra_Library
1924   ::euclidean_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
1925                               const BD_Shape<U>& x, const BD_Shape<U>& y,
1926                               const Rounding_Dir dir,
1927                               Temp& tmp0, Temp& tmp1, Temp& tmp2);
1928   template <typename Temp, typename To, typename U>
1929   friend bool Parma_Polyhedra_Library
1930   ::l_infinity_distance_assign(Checked_Number<To, Extended_Number_Policy>& r,
1931                                const BD_Shape<U>& x, const BD_Shape<U>& y,
1932                                const Rounding_Dir dir,
1933                                Temp& tmp0, Temp& tmp1, Temp& tmp2);
1934 
1935 private:
1936   template <typename U> friend class Parma_Polyhedra_Library::BD_Shape;
1937   template <typename Interval> friend class Parma_Polyhedra_Library::Box;
1938 
1939   //! The matrix representing the system of bounded differences.
1940   DB_Matrix<N> dbm;
1941 
1942 #define PPL_IN_BD_Shape_CLASS
1943 #include "BDS_Status_idefs.hh"
1944 #undef PPL_IN_BD_Shape_CLASS
1945 
1946   //! The status flags to keep track of the internal state.
1947   Status status;
1948 
1949   //! A matrix indicating which constraints are redundant.
1950   Bit_Matrix redundancy_dbm;
1951 
1952   //! Returns <CODE>true</CODE> if the BDS is the zero-dimensional universe.
1953   bool marked_zero_dim_univ() const;
1954 
1955   /*! \brief
1956     Returns <CODE>true</CODE> if the BDS is known to be empty.
1957 
1958     The return value <CODE>false</CODE> does not necessarily
1959     implies that \p *this is non-empty.
1960   */
1961   bool marked_empty() const;
1962 
1963   /*! \brief
1964     Returns <CODE>true</CODE> if the system of bounded differences
1965     is known to be shortest-path closed.
1966 
1967     The return value <CODE>false</CODE> does not necessarily
1968     implies that <CODE>this->dbm</CODE> is not shortest-path closed.
1969   */
1970   bool marked_shortest_path_closed() const;
1971 
1972   /*! \brief
1973     Returns <CODE>true</CODE> if the system of bounded differences
1974     is known to be shortest-path reduced.
1975 
1976     The return value <CODE>false</CODE> does not necessarily
1977     implies that <CODE>this->dbm</CODE> is not shortest-path reduced.
1978   */
1979   bool marked_shortest_path_reduced() const;
1980 
1981   //! Turns \p *this into an empty BDS.
1982   void set_empty();
1983 
1984   //! Turns \p *this into an zero-dimensional universe BDS.
1985   void set_zero_dim_univ();
1986 
1987   //! Marks \p *this as shortest-path closed.
1988   void set_shortest_path_closed();
1989 
1990   //! Marks \p *this as shortest-path closed.
1991   void set_shortest_path_reduced();
1992 
1993   //! Marks \p *this as possibly not shortest-path closed.
1994   void reset_shortest_path_closed();
1995 
1996   //! Marks \p *this as possibly not shortest-path reduced.
1997   void reset_shortest_path_reduced();
1998 
1999   //! Assigns to <CODE>this->dbm</CODE> its shortest-path closure.
2000   void shortest_path_closure_assign() const;
2001 
2002   /*! \brief
2003     Assigns to <CODE>this->dbm</CODE> its shortest-path closure and
2004     records into <CODE>this->redundancy_dbm</CODE> which of the entries
2005     in <CODE>this->dbm</CODE> are redundant.
2006   */
2007   void shortest_path_reduction_assign() const;
2008 
2009   /*! \brief
2010     Returns <CODE>true</CODE> if and only if <CODE>this->dbm</CODE>
2011     is shortest-path closed and <CODE>this->redundancy_dbm</CODE>
2012     correctly flags the redundant entries in <CODE>this->dbm</CODE>.
2013   */
2014   bool is_shortest_path_reduced() const;
2015 
2016   /*! \brief
2017     Incrementally computes shortest-path closure, assuming that only
2018     constraints affecting variable \p var need to be considered.
2019 
2020     \note
2021     It is assumed that \c *this, which was shortest-path closed,
2022     has only been modified by adding constraints affecting variable
2023     \p var. If this assumption is not satisfied, i.e., if a non-redundant
2024     constraint not affecting variable \p var has been added, the behavior
2025     is undefined.
2026   */
2027   void incremental_shortest_path_closure_assign(Variable var) const;
2028 
2029   //! Checks if and how \p expr is bounded in \p *this.
2030   /*!
2031     Returns <CODE>true</CODE> if and only if \p from_above is
2032     <CODE>true</CODE> and \p expr is bounded from above in \p *this,
2033     or \p from_above is <CODE>false</CODE> and \p expr is bounded
2034     from below in \p *this.
2035 
2036     \param expr
2037     The linear expression to test;
2038 
2039     \param from_above
2040     <CODE>true</CODE> if and only if the boundedness of interest is
2041     "from above".
2042 
2043     \exception std::invalid_argument
2044     Thrown if \p expr and \p *this are dimension-incompatible.
2045   */
2046   bool bounds(const Linear_Expression& expr, bool from_above) const;
2047 
2048   //! Maximizes or minimizes \p expr subject to \p *this.
2049   /*!
2050     \param expr
2051     The linear expression to be maximized or minimized subject to \p
2052     *this;
2053 
2054     \param maximize
2055     <CODE>true</CODE> if maximization is what is wanted;
2056 
2057     \param ext_n
2058     The numerator of the extremum value;
2059 
2060     \param ext_d
2061     The denominator of the extremum value;
2062 
2063     \param included
2064     <CODE>true</CODE> if and only if the extremum of \p expr can
2065     actually be reached in \p * this;
2066 
2067     \param g
2068     When maximization or minimization succeeds, will be assigned
2069     a point or closure point where \p expr reaches the
2070     corresponding extremum value.
2071 
2072     \exception std::invalid_argument
2073     Thrown if \p expr and \p *this are dimension-incompatible.
2074 
2075     If \p *this is empty or \p expr is not bounded in the appropriate
2076     direction, <CODE>false</CODE> is returned and \p ext_n, \p ext_d,
2077     \p included and \p g are left untouched.
2078   */
2079   bool max_min(const Linear_Expression& expr,
2080                bool maximize,
2081                Coefficient& ext_n, Coefficient& ext_d, bool& included,
2082                Generator& g) const;
2083 
2084   //! Maximizes or minimizes \p expr subject to \p *this.
2085   /*!
2086     \param expr
2087     The linear expression to be maximized or minimized subject to \p
2088     *this;
2089 
2090     \param maximize
2091     <CODE>true</CODE> if maximization is what is wanted;
2092 
2093     \param ext_n
2094     The numerator of the extremum value;
2095 
2096     \param ext_d
2097     The denominator of the extremum value;
2098 
2099     \param included
2100     <CODE>true</CODE> if and only if the extremum of \p expr can
2101     actually be reached in \p * this;
2102 
2103     \exception std::invalid_argument
2104     Thrown if \p expr and \p *this are dimension-incompatible.
2105 
2106     If \p *this is empty or \p expr is not bounded in the appropriate
2107     direction, <CODE>false</CODE> is returned and \p ext_n, \p ext_d,
2108     \p included and \p point are left untouched.
2109   */
2110   bool max_min(const Linear_Expression& expr,
2111                bool maximize,
2112                Coefficient& ext_n, Coefficient& ext_d, bool& included) const;
2113 
2114   /*! \brief
2115     If the upper bound of \p *this and \p y is exact it is assigned
2116     to \p *this and \c true is returned, otherwise \c false is returned.
2117 
2118     Current implementation is based on a variant of Algorithm 4.1 in
2119       A. Bemporad, K. Fukuda, and F. D. Torrisi
2120       <em>Convexity Recognition of the Union of Polyhedra</em>
2121       Technical Report AUT00-13, ETH Zurich, 2000
2122     tailored to the special case of BD shapes.
2123 
2124     \note
2125     It is assumed that \p *this and \p y are dimension-compatible;
2126     if the assumption does not hold, the behavior is undefined.
2127   */
2128   bool BFT00_upper_bound_assign_if_exact(const BD_Shape& y);
2129 
2130   /*! \brief
2131     If the upper bound of \p *this and \p y is exact it is assigned
2132     to \p *this and \c true is returned, otherwise \c false is returned.
2133 
2134     Implementation for the rational (resp., integer) case is based on
2135     Theorem 5.2 (resp. Theorem 5.3) of \ref BHZ09b "[BHZ09b]".
2136     The Boolean template parameter \c integer_upper_bound allows for
2137     choosing between the rational and integer upper bound algorithms.
2138 
2139     \note
2140     It is assumed that \p *this and \p y are dimension-compatible;
2141     if the assumption does not hold, the behavior is undefined.
2142 
2143     \note
2144     The integer case is only enabled if T is an integer data type.
2145   */
2146   template <bool integer_upper_bound>
2147   bool BHZ09_upper_bound_assign_if_exact(const BD_Shape& y);
2148 
2149   /*! \brief
2150     Uses the constraint \p c to refine \p *this.
2151 
2152     \param c
2153     The constraint to be added. Non BD constraints are ignored.
2154 
2155     \warning
2156     If \p c and \p *this are dimension-incompatible,
2157     the behavior is undefined.
2158   */
2159   void refine_no_check(const Constraint& c);
2160 
2161   /*! \brief
2162     Uses the congruence \p cg to refine \p *this.
2163 
2164     \param cg
2165     The congruence to be added.
2166     Nontrivial proper congruences are ignored.
2167     Non BD equalities are ignored.
2168 
2169     \warning
2170     If \p cg and \p *this are dimension-incompatible,
2171     the behavior is undefined.
2172   */
2173   void refine_no_check(const Congruence& cg);
2174 
2175   //! Adds the constraint <CODE>dbm[i][j] \<= k</CODE>.
2176   void add_dbm_constraint(dimension_type i, dimension_type j, const N& k);
2177 
2178   //! Adds the constraint <CODE>dbm[i][j] \<= numer/denom</CODE>.
2179   void add_dbm_constraint(dimension_type i, dimension_type j,
2180                           Coefficient_traits::const_reference numer,
2181                           Coefficient_traits::const_reference denom);
2182 
2183   /*! \brief
2184     Adds to the BDS the constraint
2185     \f$\mathrm{var} \relsym \frac{\mathrm{expr}}{\mathrm{denominator}}\f$.
2186 
2187     Note that the coefficient of \p var in \p expr is null.
2188   */
2189   void refine(Variable var, Relation_Symbol relsym,
2190               const Linear_Expression& expr,
2191               Coefficient_traits::const_reference denominator
2192               = Coefficient_one());
2193 
2194   //! Removes all the constraints on row/column \p v.
2195   void forget_all_dbm_constraints(dimension_type v);
2196   //! Removes all binary constraints on row/column \p v.
2197   void forget_binary_dbm_constraints(dimension_type v);
2198 
2199   //! An helper function for the computation of affine relations.
2200   /*!
2201     For each dbm index \p u (less than or equal to \p last_v and different
2202     from \p v), deduce constraints of the form <CODE>v - u \<= c</CODE>,
2203     starting from \p ub_v which is an upper bound for \p v.
2204 
2205     The shortest-path closure is able to deduce the constraint
2206     <CODE>v - u \<= ub_v - lb_u</CODE>. We can be more precise if variable
2207     \p u played an active role in the computation of the upper bound for
2208     \p v, i.e., if the corresponding coefficient
2209     <CODE>q == sc_expr[u]/sc_denom</CODE> is greater than zero. In particular:
2210       - if <CODE>q \>= 1</CODE>, then <CODE>v - u \<= ub_v - ub_u</CODE>;
2211       - if <CODE>0 \< q \< 1</CODE>, then
2212         <CODE>v - u \<= ub_v - (q*ub_u + (1-q)*lb_u)</CODE>.
2213   */
2214   void deduce_v_minus_u_bounds(dimension_type v,
2215                                dimension_type last_v,
2216                                const Linear_Expression& sc_expr,
2217                                Coefficient_traits::const_reference sc_denom,
2218                                const N& ub_v);
2219 
2220   /*! \brief
2221     Auxiliary function for \ref affine_form_relation "affine form image" that
2222     handle the general case: \f$l = c\f$
2223   */
2224   template <typename Interval_Info>
2225   void inhomogeneous_affine_form_image(const dimension_type& var_id,
2226                                        const Interval<T, Interval_Info>& b);
2227 
2228   /*! \brief
2229     Auxiliary function for \ref affine_form_relation "affine form
2230     image" that handle the general case: \f$l = ax + c\f$
2231   */
2232   template <typename Interval_Info>
2233   void
2234   one_variable_affine_form_image(const dimension_type& var_id,
2235                                  const Interval<T, Interval_Info>& b,
2236                                  const Interval<T, Interval_Info>& w_coeff,
2237                                  const dimension_type& w_id,
2238                                  const dimension_type& space_dim);
2239 
2240   /*! \brief
2241     Auxiliary function for \ref affine_form_relation "affine form image" that
2242     handle the general case: \f$l = ax + by + c\f$
2243   */
2244   template <typename Interval_Info>
2245   void
2246   two_variables_affine_form_image(const dimension_type& var_id,
2247                                   const Linear_Form<Interval<T,Interval_Info> >& lf,
2248                                   const dimension_type& space_dim);
2249 
2250   /*! \brief
2251     Auxiliary function for refine with linear form that handle
2252     the general case: \f$l = ax + c\f$
2253   */
2254   template <typename Interval_Info>
2255   void
2256   left_inhomogeneous_refine(const dimension_type& right_t,
2257                             const dimension_type& right_w_id,
2258                             const Linear_Form<Interval<T, Interval_Info> >& left,
2259                             const Linear_Form<Interval<T, Interval_Info> >& right);
2260 
2261   /*! \brief
2262     Auxiliary function for refine with linear form that handle
2263     the general case: \f$ax + b = cy + d\f$
2264   */
2265   template <typename Interval_Info>
2266   void
2267   left_one_var_refine(const dimension_type& left_w_id,
2268                       const dimension_type& right_t,
2269                       const dimension_type& right_w_id,
2270                       const Linear_Form<Interval<T, Interval_Info> >& left,
2271                       const Linear_Form<Interval<T, Interval_Info> >& right);
2272 
2273   /*! \brief
2274     Auxiliary function for refine with linear form that handle
2275     the general case.
2276   */
2277   template <typename Interval_Info>
2278   void general_refine(const dimension_type& left_w_id,
2279                       const dimension_type& right_w_id,
2280                       const Linear_Form<Interval<T, Interval_Info> >& left,
2281                       const Linear_Form<Interval<T, Interval_Info> >& right);
2282 
2283   template <typename Interval_Info>
2284   void linear_form_upper_bound(const Linear_Form<Interval<T, Interval_Info> >&
2285                                lf,
2286                                N& result) const;
2287 
2288   //! An helper function for the computation of affine relations.
2289   /*!
2290     For each dbm index \p u (less than or equal to \p last_v and different
2291     from \p v), deduce constraints of the form <CODE>u - v \<= c</CODE>,
2292     starting from \p minus_lb_v which is a lower bound for \p v.
2293 
2294     The shortest-path closure is able to deduce the constraint
2295     <CODE>u - v \<= ub_u - lb_v</CODE>. We can be more precise if variable
2296     \p u played an active role in the computation of the lower bound for
2297     \p v, i.e., if the corresponding coefficient
2298     <CODE>q == sc_expr[u]/sc_denom</CODE> is greater than zero.
2299     In particular:
2300       - if <CODE>q \>= 1</CODE>, then <CODE>u - v \<= lb_u - lb_v</CODE>;
2301       - if <CODE>0 \< q \< 1</CODE>, then
2302         <CODE>u - v \<= (q*lb_u + (1-q)*ub_u) - lb_v</CODE>.
2303   */
2304   void deduce_u_minus_v_bounds(dimension_type v,
2305                                dimension_type last_v,
2306                                const Linear_Expression& sc_expr,
2307                                Coefficient_traits::const_reference sc_denom,
2308                                const N& minus_lb_v);
2309 
2310   /*! \brief
2311     Adds to \p limiting_shape the bounded differences in \p cs
2312     that are satisfied by \p *this.
2313   */
2314   void get_limiting_shape(const Constraint_System& cs,
2315                           BD_Shape& limiting_shape) const;
2316 
2317   //! Compute the (zero-equivalence classes) predecessor relation.
2318   /*!
2319     It is assumed that the BDS is not empty and shortest-path closed.
2320   */
2321   void compute_predecessors(std::vector<dimension_type>& predecessor) const;
2322 
2323   //! Compute the leaders of zero-equivalence classes.
2324   /*!
2325     It is assumed that the BDS is not empty and shortest-path closed.
2326   */
2327   void compute_leaders(std::vector<dimension_type>& leaders) const;
2328 
2329   void drop_some_non_integer_points_helper(N& elem);
2330 
2331   friend std::ostream&
2332   Parma_Polyhedra_Library::IO_Operators
2333   ::operator<<<>(std::ostream& s, const BD_Shape<T>& c);
2334 
2335   //! \name Exception Throwers
2336   //@{
2337   void throw_dimension_incompatible(const char* method,
2338                                     const BD_Shape& y) const;
2339 
2340   void throw_dimension_incompatible(const char* method,
2341                                     dimension_type required_dim) const;
2342 
2343   void throw_dimension_incompatible(const char* method,
2344                                     const Constraint& c) const;
2345 
2346   void throw_dimension_incompatible(const char* method,
2347                                     const Congruence& cg) const;
2348 
2349   void throw_dimension_incompatible(const char* method,
2350                                     const Generator& g) const;
2351 
2352   void throw_dimension_incompatible(const char* method,
2353                                     const char* le_name,
2354                                     const Linear_Expression& le) const;
2355 
2356   template<typename Interval_Info>
2357   void
2358   throw_dimension_incompatible(const char* method,
2359                                const char* lf_name,
2360                                const Linear_Form<Interval<T, Interval_Info> >&
2361                                lf) const;
2362 
2363   static void throw_expression_too_complex(const char* method,
2364                                            const Linear_Expression& le);
2365 
2366   static void throw_invalid_argument(const char* method, const char* reason);
2367   //@} // Exception Throwers
2368 };
2369 
2370 #include "BDS_Status_inlines.hh"
2371 #include "BD_Shape_inlines.hh"
2372 #include "BD_Shape_templates.hh"
2373 
2374 #endif // !defined(PPL_BD_Shape_defs_hh)
2375