DGtal  1.1.0
BoundedRationalPolytope.h
1 
17 #pragma once
18 
31 #if defined(BoundedRationalPolytope_RECURSES)
32 #error Recursive header files inclusion detected in BoundedRationalPolytope.h
33 #else // defined(BoundedRationalPolytope_RECURSES)
34 
35 #define BoundedRationalPolytope_RECURSES
36 
37 #if !defined BoundedRationalPolytope_h
38 
39 #define BoundedRationalPolytope_h
40 
42 // Inclusions
43 #include <iostream>
44 #include <list>
45 #include <vector>
46 #include <string>
47 #include "DGtal/base/Common.h"
48 #include "DGtal/kernel/CSpace.h"
49 #include "DGtal/kernel/domains/HyperRectDomain.h"
50 #include "DGtal/arithmetic/IntegerComputer.h"
51 #include "DGtal/arithmetic/ClosedIntegerHalfPlane.h"
53 
54 namespace DGtal
55 {
56 
58  // template class BoundedRationalPolytope
73  template < typename TSpace >
75  {
77 
78  public:
80  typedef TSpace Space;
81  typedef typename Space::Integer Integer;
82  typedef typename Space::Point Point;
83  typedef typename Space::Vector Vector;
84  typedef std::vector<Vector> InequalityMatrix;
85  typedef std::vector<Integer> InequalityVector;
88 #ifdef WITH_BIGINTEGER
90 #else
91  typedef DGtal::int64_t BigInteger;
92 #endif
94 
99  struct UnitSegment {
101  UnitSegment( Dimension d ) : k( d ) {}
102  };
103 
109  struct UnitCell {
110  std::vector<Dimension> dims;
111  UnitCell( std::initializer_list<Dimension> l )
112  : dims( l.begin(), l.end() ) {}
113 
120  friend std::ostream&
121  operator<< ( std::ostream & out,
122  const UnitCell & object )
123  {
124  out << "{";
125  for ( Dimension i = 0; i < object.dims.size(); ++i ) out << object.dims[ i ];
126  out << "}";
127  return out;
128  }
129  };
130 
131 
134  struct Rational {
135  Integer p; // numerator
136  Integer q; // denominator
140  inline Rational( Integer a, Integer b ) : p( a ), q( b ) {}
141  };
142 
145 
150 
155 
160  BoundedRationalPolytope ( const Self & other ) = default;
161 
162 
173  BoundedRationalPolytope( std::initializer_list<Point> l );
174 
188  template <typename PointIterator>
189  BoundedRationalPolytope( Integer d, PointIterator itB, PointIterator itE );
190 
211  template <typename HalfSpaceIterator>
213  const Domain& domain,
214  HalfSpaceIterator itB, HalfSpaceIterator itE,
215  bool valid_edge_constraints = false );
216 
237  template <typename HalfSpaceIterator>
238  void init( Integer d, const Domain& domain,
239  HalfSpaceIterator itB, HalfSpaceIterator itE,
240  bool valid_edge_constraints = false );
241 
242 
260  template <typename PointIterator>
261  bool init( Integer d, PointIterator itB, PointIterator itE );
262 
268  Self & operator= ( const Self & other ) = default;
269 
271  void clear();
272 
274 
275  // ----------------------- Accessor services ------------------------------
276  public:
279 
281  const Domain& getDomain() const;
282 
285  const Domain& getLatticeDomain() const;
286 
289  const Domain& getRationalDomain() const;
290 
292  unsigned int nbHalfSpaces() const;
293 
298 
304  const Vector& getA( unsigned int i ) const;
305 
311  Integer getB( unsigned int i ) const;
312 
318  bool isLarge( unsigned int i ) const;
319 
321  const InequalityMatrix& getA() const;
322 
324  const InequalityVector& getB() const;
328  const std::vector<bool>& getI() const;
329 
334  bool canBeSummed() const;
335 
337 
338  // ----------------------- Check point services ------------------------------
339  public:
340 
343 
348  bool isInside( const Point& p ) const;
349 
356  bool isDomainPointInside( const Point& p ) const;
357 
362  bool isInterior( const Point& p ) const;
363 
368  bool isBoundary( const Point& p ) const;
369 
371 
372  // ----------------------- Modification services ------------------------------
373  public:
374 
377 
381 
394  unsigned int cut( Dimension k, bool pos, Integer b, bool large = true );
395 
413  unsigned int cut( const Vector& a, Integer b, bool large = true,
414  bool valid_edge_constraint = false );
415 
432  unsigned int cut( const HalfSpace & hs, bool large = true,
433  bool valid_edge_constraint = false );
434 
439  void swap( BoundedRationalPolytope & other );
440 
441 
448 
455 
463 
471 
473 
474  // ----------------------- Enumeration services ------------------------------
475  public:
476 
479 
487  Integer count() const;
488 
500 
512 
523  Integer countWithin( Point low, Point hi ) const;
524 
542 
551  void getPoints( std::vector<Point>& pts ) const;
552 
561  void getInteriorPoints( std::vector<Point>& pts ) const;
562 
571  void getBoundaryPoints( std::vector<Point>& pts ) const;
572 
583  template <typename PointSet>
584  void insertPoints( PointSet& pts_set ) const;
585 
587 
588 
589  // ----------------------- Interface --------------------------------------
590  public:
593 
598  void selfDisplay ( std::ostream & out ) const;
599 
606  bool isValid() const;
607 
611  std::string className() const;
612 
614 
615  // ------------------------- Protected Datas ------------------------------
616  protected:
617  // Denominator for constraints, i.e. \f$ q A x \le B \f$.
619  // The matrix A in the polytope representation \f$ q A x \le B \f$.
621  // The vector B in the polytope representation \f$ q A x \le B \f$.
623  // Tight bounded box
625  // Lattice bounded box (i.e. floor( rationalD/q ))
627  // Are inequalities large ?
628  std::vector<bool> I;
629  // Indicates if Minkowski sums with segments will be valid
631 
632  // ------------------------- Private Datas --------------------------------
633  private:
634 
635 
636  // ------------------------- Internals ------------------------------------
637  private:
645 
652 
659 
667 
674 
675  }; // end of class BoundedRationalPolytope
676 
677  namespace detail {
686  template <DGtal::Dimension N, typename TInteger>
688  typedef TInteger Integer;
690  typedef typename Space::Point Point;
691  typedef typename Space::Vector Vector;
694 
707  static void
708  addEdgeConstraint( Polytope& , unsigned int , unsigned int ,
709  const std::vector<Point>& )
710  {
711  trace.error() << "[BoundedRationalPolytopeHelper::addEdgeConstraint]"
712  << " this method is only implemented in 3D." << std::endl;
713  }
714 
717  static
718  Vector crossProduct( const Vector& , const Vector& )
719  {
720  trace.error() << "[BoundedRationalPolytopeHelper::crossProduct]"
721  << " this method is only implemented in 3D." << std::endl;
722  return Vector::zero;
723  }
724  };
725 
733  template <typename TInteger>
735  typedef TInteger Integer;
737  typedef typename Space::Point Point;
738  typedef typename Space::Vector Vector;
741 
751  static void
752  addEdgeConstraint( Polytope& P, unsigned int i, unsigned int j,
753  const std::vector<Point>& pts )
754  {
755  Vector ab = pts[ i ] - pts[ j ];
756  for ( int s = 0; s < 2; s++ )
757  for ( Dimension k = 0; k < dimension; ++k )
758  {
759  Vector n = ab.crossProduct( Point::base( k, (s == 0) ? 1 : -1 ) );
760  Integer b = n.dot( pts[ i ] );
761  int nb_in = 0;
762  for ( auto p : pts ) {
763  Integer v = n.dot( p );
764  if ( v < b ) nb_in++;
765  }
766  if ( nb_in == dimension - 1 ) {
767  P.cut( n, b, true, true );
768  }
769  }
770  }
775  static
776  Vector crossProduct( const Vector& v1, const Vector& v2 )
777  {
778  return v1.crossProduct( v2 );
779  }
780  };
781  }
782 
785 
792  template <typename TSpace>
793  std::ostream&
794  operator<< ( std::ostream & out,
795  const BoundedRationalPolytope<TSpace> & object );
796 
797 
803  template <typename TSpace>
807 
813  template <typename TSpace>
817 
818 
826  template <typename TSpace>
830 
838  template <typename TSpace>
842 
844 
845 } // namespace DGtal
846 
847 
849 // Includes inline functions.
850 #include "BoundedRationalPolytope.ih"
851 
852 // //
854 
855 #endif // !defined BoundedRationalPolytope_h
856 
857 #undef BoundedRationalPolytope_RECURSES
858 #endif // else defined(BoundedRationalPolytope_RECURSES)
DGtal::BigInteger
mpz_class BigInteger
Multi-precision integer with GMP implementation.
Definition: BasicTypes.h:79
DGtal::BoundedRationalPolytope::getA
const Vector & getA(unsigned int i) const
DGtal::BoundedRationalPolytope::isInterior
bool isInterior(const Point &p) const
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Space
SpaceND< 3, Integer > Space
Definition: BoundedRationalPolytope.h:736
DGtal::PointVector< dim, Integer >::zero
static Self zero
Static const for zero PointVector.
Definition: PointVector.h:1590
DGtal::detail::BoundedRationalPolytopeSpecializer
Aim: It is just a helper class for BoundedRationalPolytope to add dimension specific static methods.
Definition: BoundedRationalPolytope.h:687
DGtal::BoundedRationalPolytope::Rational::p
Integer p
Definition: BoundedRationalPolytope.h:135
DGtal::BoundedRationalPolytope::Rational::Rational
Rational(Integer a, Integer b)
Definition: BoundedRationalPolytope.h:140
DGtal::HyperRectDomain< Space >
DGtal::BoundedRationalPolytope::isLarge
bool isLarge(unsigned int i) const
DGtal::concepts::CSpace
Aim: Defines the concept describing a digital space, ie a cartesian product of integer lines.
Definition: CSpace.h:106
max
int max(int a, int b)
Definition: testArithmeticalDSS.cpp:1108
DGtal::BoundedRationalPolytope::operator=
Self & operator=(const Self &other)=default
DGtal::BoundedRationalPolytope::internalInitFromSegment3D
bool internalInitFromSegment3D(Point a, Point b)
DGtal::BoundedRationalPolytope::getInteriorPoints
void getInteriorPoints(std::vector< Point > &pts) const
DGtal::Trace::error
std::ostream & error()
DGtal::BoundedRationalPolytope::UnitCell::UnitCell
UnitCell(std::initializer_list< Dimension > l)
Definition: BoundedRationalPolytope.h:111
DGtal::BoundedRationalPolytope::init
bool init(Integer d, PointIterator itB, PointIterator itE)
DGtal::BoundedRationalPolytope::isDomainPointInside
bool isDomainPointInside(const Point &p) const
DGtal::BoundedRationalPolytope::Vector
Space::Vector Vector
Definition: BoundedRationalPolytope.h:83
DGtal::BoundedRationalPolytope::latticeD
Domain latticeD
Definition: BoundedRationalPolytope.h:626
DGtal::BoundedRationalPolytope::dimension
static const Dimension dimension
Definition: BoundedRationalPolytope.h:93
DGtal::BoundedRationalPolytope::BOOST_CONCEPT_ASSERT
BOOST_CONCEPT_ASSERT((concepts::CSpace< TSpace >))
DGtal::BoundedRationalPolytope::countUpTo
Integer countUpTo(Integer max) const
DGtal::BoundedRationalPolytope::I
std::vector< bool > I
Definition: BoundedRationalPolytope.h:628
DGtal::trace
Trace trace
Definition: Common.h:150
DGtal::BoundedRationalPolytope::Space
TSpace Space
Definition: BoundedRationalPolytope.h:80
DGtal::BoundedRationalPolytope::countInterior
Integer countInterior() const
DGtal::Dimension
DGtal::uint32_t Dimension
Definition: Common.h:133
DGtal::BoundedRationalPolytope::getA
const InequalityMatrix & getA() const
DGtal::BoundedRationalPolytope::isInside
bool isInside(const Point &p) const
DGtal::BoundedRationalPolytope::canBeSummed
bool canBeSummed() const
DGtal::detail::BoundedRationalPolytopeSpecializer::addEdgeConstraint
static void addEdgeConstraint(Polytope &, unsigned int, unsigned int, const std::vector< Point > &)
Definition: BoundedRationalPolytope.h:708
DGtal::BoundedRationalPolytope::Self
BoundedRationalPolytope< TSpace > Self
Definition: BoundedRationalPolytope.h:79
DGtal::detail::BoundedRationalPolytopeSpecializer::dimension
static const Dimension dimension
Definition: BoundedRationalPolytope.h:693
DGtal::BoundedRationalPolytope::BoundedRationalPolytope
BoundedRationalPolytope(const Self &other)=default
DGtal::BoundedRationalPolytope::nbHalfSpaces
unsigned int nbHalfSpaces() const
DGtal::detail::BoundedRationalPolytopeSpecializer::Integer
TInteger Integer
Definition: BoundedRationalPolytope.h:688
DGtal::BoundedRationalPolytope::BoundedRationalPolytope
BoundedRationalPolytope()
DGtal::BoundedRationalPolytope::Domain
HyperRectDomain< Space > Domain
Definition: BoundedRationalPolytope.h:86
DGtal::PointVector< dim, Integer >::base
static Self base(Dimension k, Component val=1)
DGtal::BoundedRationalPolytope::InequalityVector
std::vector< Integer > InequalityVector
Definition: BoundedRationalPolytope.h:85
DGtal::SpaceND
Definition: SpaceND.h:96
DGtal::PointVector::crossProduct
auto crossProduct(const PointVector< dim, OtherComponent, OtherStorage > &v) const -> decltype(DGtal::crossProduct(*this, v))
Cross product with a PointVector.
DGtal::BoundedRationalPolytope::UnitCell::operator<<
friend std::ostream & operator<<(std::ostream &out, const UnitCell &object)
Definition: BoundedRationalPolytope.h:121
DGtal::BoundedRationalPolytope::operator*=
Self & operator*=(Rational r)
DGtal::BoundedRationalPolytope::InequalityMatrix
std::vector< Vector > InequalityMatrix
Definition: BoundedRationalPolytope.h:84
DGtal::BoundedRationalPolytope::getLatticeDomain
const Domain & getLatticeDomain() const
DGtal::BoundedRationalPolytope::getB
Integer getB(unsigned int i) const
DGtal::BoundedRationalPolytope::selfDisplay
void selfDisplay(std::ostream &out) const
DGtal::detail::BoundedRationalPolytopeSpecializer::crossProduct
static Vector crossProduct(const Vector &, const Vector &)
Definition: BoundedRationalPolytope.h:718
DGtal::BoundedRationalPolytope::computeLatticeDomain
Domain computeLatticeDomain(const Domain &d)
DGtal::operator+
Circulator< TIterator > operator+(typename IteratorCirculatorTraits< TIterator >::Difference d, Circulator< TIterator > &object)
Definition: Circulator.h:453
DGtal::BoundedRationalPolytope::BoundedRationalPolytope
BoundedRationalPolytope(Integer d, const Domain &domain, HalfSpaceIterator itB, HalfSpaceIterator itE, bool valid_edge_constraints=false)
DGtal::detail::BoundedRationalPolytopeSpecializer::Space
SpaceND< N, Integer > Space
Definition: BoundedRationalPolytope.h:689
DGtal
DGtal is the top-level namespace which contains all DGtal functions and types.
Definition: ClosedIntegerHalfPlane.h:49
DGtal::BoundedRationalPolytope::swap
void swap(BoundedRationalPolytope &other)
DGtal::detail::BoundedRationalPolytopeSpecializer::Polytope
BoundedRationalPolytope< Space > Polytope
Definition: BoundedRationalPolytope.h:692
DGtal::BoundedRationalPolytope::Rational::q
Integer q
Definition: BoundedRationalPolytope.h:136
DGtal::BoundedRationalPolytope::getBoundaryPoints
void getBoundaryPoints(std::vector< Point > &pts) const
DGtal::BoundedRationalPolytope::cut
unsigned int cut(const Vector &a, Integer b, bool large=true, bool valid_edge_constraint=false)
DGtal::BoundedRationalPolytope::getDomain
const Domain & getDomain() const
DGtal::BoundedRationalPolytope::UnitSegment::k
Dimension k
Definition: BoundedRationalPolytope.h:100
DGtal::BoundedRationalPolytope::rationalD
Domain rationalD
Definition: BoundedRationalPolytope.h:624
DGtal::SpaceND::dimension
static const Dimension dimension
static constants to store the dimension.
Definition: SpaceND.h:132
DGtal::BoundedRationalPolytope::UnitCell::dims
std::vector< Dimension > dims
Definition: BoundedRationalPolytope.h:110
DGtal::BoundedRationalPolytope::Point
Space::Point Point
Definition: BoundedRationalPolytope.h:82
DGtal::BoundedRationalPolytope::UnitSegment::UnitSegment
UnitSegment(Dimension d)
Definition: BoundedRationalPolytope.h:101
DGtal::BoundedRationalPolytope::init
void init(Integer d, const Domain &domain, HalfSpaceIterator itB, HalfSpaceIterator itE, bool valid_edge_constraints=false)
DGtal::BoundedRationalPolytope::~BoundedRationalPolytope
~BoundedRationalPolytope()=default
DGtal::detail::BoundedRationalPolytopeSpecializer::Vector
Space::Vector Vector
Definition: BoundedRationalPolytope.h:691
DGtal::BoundedRationalPolytope::Integer
Space::Integer Integer
Definition: BoundedRationalPolytope.h:81
DGtal::BoundedRationalPolytope::B
InequalityVector B
Definition: BoundedRationalPolytope.h:622
DGtal::operator<<
std::ostream & operator<<(std::ostream &out, const ClosedIntegerHalfPlane< TSpace > &object)
DGtal::BoundedRationalPolytope::denominator
Integer denominator() const
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::crossProduct
static Vector crossProduct(const Vector &v1, const Vector &v2)
Definition: BoundedRationalPolytope.h:776
DGtal::BoundedRationalPolytope::count
Integer count() const
DGtal::BoundedRationalPolytope::Rational
Definition: BoundedRationalPolytope.h:134
DGtal::BoundedRationalPolytope::insertPoints
void insertPoints(PointSet &pts_set) const
DGtal::BoundedRationalPolytope::BoundedRationalPolytope
BoundedRationalPolytope(std::initializer_list< Point > l)
DGtal::ClosedIntegerHalfPlane
Aim: A half-space specified by a vector N and a constant c. The half-space is the set .
Definition: ClosedIntegerHalfPlane.h:64
DGtal::BoundedRationalPolytope
Aim: Represents an nD rational polytope, i.e. a convex polyhedron bounded by vertices with rational c...
Definition: BoundedRationalPolytope.h:75
DGtal::BoundedRationalPolytope::BigInteger
DGtal::BigInteger BigInteger
Definition: BoundedRationalPolytope.h:89
DGtal::BoundedRationalPolytope::isBoundary
bool isBoundary(const Point &p) const
DGtal::BoundedRationalPolytope::className
std::string className() const
DGtal::BoundedRationalPolytope::countBoundary
Integer countBoundary() const
DGtal::BoundedRationalPolytope::HalfSpace
ClosedIntegerHalfPlane< Space > HalfSpace
Definition: BoundedRationalPolytope.h:87
DGtal::BoundedRationalPolytope::operator*=
Self & operator*=(Integer t)
DGtal::detail::BoundedRationalPolytopeSpecializer::Point
Space::Point Point
Definition: BoundedRationalPolytope.h:690
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Integer
TInteger Integer
Definition: BoundedRationalPolytope.h:735
domain
Domain domain
Definition: testProjection.cpp:88
DGtal::BoundedRationalPolytope::operator+=
Self & operator+=(UnitCell c)
DGtal::PointVector< dim, Integer >
DGtal::BoundedRationalPolytope::internalInitFromTriangle3D
bool internalInitFromTriangle3D(Point a, Point b, Point c)
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Vector
Space::Vector Vector
Definition: BoundedRationalPolytope.h:738
DGtal::BoundedRationalPolytope::countWithin
Integer countWithin(Point low, Point hi) const
DGtal::BoundedRationalPolytope::getPoints
void getPoints(std::vector< Point > &pts) const
DGtal::BoundedRationalPolytope::operator+=
Self & operator+=(UnitSegment s)
DGtal::BoundedRationalPolytope::UnitCell
Definition: BoundedRationalPolytope.h:109
DGtal::operator*
KForm< Calculus, order, duality > operator*(const typename Calculus::Scalar &scalar, const KForm< Calculus, order, duality > &form)
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Point
Space::Point Point
Definition: BoundedRationalPolytope.h:737
DGtal::int64_t
boost::int64_t int64_t
signed 94-bit integer.
Definition: BasicTypes.h:74
DGtal::BoundedRationalPolytope::BoundedRationalPolytope
BoundedRationalPolytope(Integer d, PointIterator itB, PointIterator itE)
DGtal::BoundedRationalPolytope::isValid
bool isValid() const
DGtal::BoundedRationalPolytope::cut
unsigned int cut(const HalfSpace &hs, bool large=true, bool valid_edge_constraint=false)
DGtal::BoundedRationalPolytope::q
Integer q
Definition: BoundedRationalPolytope.h:618
DGtal::BoundedRationalPolytope::A
InequalityMatrix A
Definition: BoundedRationalPolytope.h:620
DGtal::BoundedRationalPolytope::getI
const std::vector< bool > & getI() const
DGtal::BoundedRationalPolytope::internalInitFromSegment2D
bool internalInitFromSegment2D(Point a, Point b)
DGtal::BoundedRationalPolytope::getRationalDomain
const Domain & getRationalDomain() const
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Polytope
BoundedRationalPolytope< Space > Polytope
Definition: BoundedRationalPolytope.h:739
DGtal::BoundedRationalPolytope::cut
unsigned int cut(Dimension k, bool pos, Integer b, bool large=true)
DGtal::BoundedRationalPolytope::computeRationalDomain
Domain computeRationalDomain(const Domain &d)
DGtal::BoundedRationalPolytope::myValidEdgeConstraints
bool myValidEdgeConstraints
Definition: BoundedRationalPolytope.h:630
DGtal::BoundedRationalPolytope::UnitSegment
Definition: BoundedRationalPolytope.h:99
DGtal::BoundedRationalPolytope::interiorPolytope
BoundedRationalPolytope interiorPolytope() const
DGtal::SpaceND::Integer
TInteger Integer
Arithmetic ring induced by (+,-,*) and Integer numbers.
Definition: SpaceND.h:102
DGtal::BoundedRationalPolytope::clear
void clear()
Clears the polytope.
DGtal::BoundedRationalPolytope::getB
const InequalityVector & getB() const
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::addEdgeConstraint
static void addEdgeConstraint(Polytope &P, unsigned int i, unsigned int j, const std::vector< Point > &pts)
Definition: BoundedRationalPolytope.h:752