5 #if defined(__cplusplus) && ! defined(_MSC_VER) 8 # include <m4ri/m4ri.h> 9 # include <m4ri/mzd_bool.h> 10 # include <m4ri/russian_bool.h> 11 #if defined(__cplusplus) && !defined(_MSC_VER) 16 #include <boost/cstdint.hpp> 17 #include <boost/multiprecision/integer.hpp> 20 using namespace reach;
25 origin{
Origin::RG,-1,0,-1,-1,
nullptr,
nullptr,
nullptr},
29 mask[o1] =
IndexRange(o1, 0, model->count2(o1));
32 mtx[o1][o2] =
nullptr;
86 Q_ASSERT(
mtx[o1][o2]==
nullptr);
87 mtx[o1][o2] = mzd_init(
mask[o1].len(),
mask[o2].len());
89 Q_ASSERT(
mtx[o1][o2]->nrows==
mask[o1].len());
90 Q_ASSERT(
mtx[o1][o2]->ncols==
mask[o2].len());
134 origin.
j = A->origin.j + (B->origin.j - B->origin.i);
138 A->origin.succesors++;
139 B->origin.succesors++;
145 Q_ASSERT(A->origin.blevel == 0);
154 B->origin.succesors++;
172 Q_ASSERT(
origin.
A.get() !=
this);
173 Q_ASSERT(
origin.
B.get() !=
this);
182 Q_ASSERT(
origin.
A.get() !=
this);
183 Q_ASSERT(
origin.
B.get() !=
this);
221 mzd_free(
mtx[o1][o2]);
222 mtx[o1][o2] =
nullptr;
249 mzd_t* m1 =
mtx[o1][o2];
250 mzd_t* m2 = that.
mtx[o1][o2];
251 if (!mzd_equal(m1,m2))
265 mzd_t* m2 = that.
mtx[o1][o2];
266 mzd_t*& m1 =
mtx[o1][o2];
269 && (m1->nrows == m2->nrows)
270 && (m1->ncols == m2->ncols))
279 m1 = mzd_copy(
nullptr,m2);
305 for( ; i1; i1=i1->
next(), i2=i2->
next())
317 Q_ASSERT(lh->
ori==r2.
ori);
329 C->merge2(
this,B,pool);
340 if (! A_VV || ! B_VV)
343 mzd_t* C_VV =
new_mzd(A_VV->nrows,B_VV->ncols,pool);
347 int threads = ConcurrencyContext::countThreads();
348 mzd_optimal_parameters(A_VV, B_VV,
349 ConcurrencyContext::cacheSize(2),
350 ConcurrencyContext::cacheSize(3),
353 mzd_bool_mul_uptri(C_VV, A_VV, B_VV, k,blocksize,threads);
355 if (! mzd_is_zero(C_VV)) {
357 Q_ASSERT(C_VV->nrows == this->mask[
VERTICAL].len());
358 Q_ASSERT(C_VV->ncols == this->mask[
VERTICAL].len());
369 C->merge3(
this,B,pool);
383 if (! A_HV || ! B_VV || ! A_VH)
387 int threads = ConcurrencyContext::countThreads();
388 mzd_optimal_parameters(A_HV, B_VV,
389 ConcurrencyContext::cacheSize(2),
390 ConcurrencyContext::cacheSize(3),
393 mzd_t* temp =
new_mzd(A_HV->nrows,B_VV->ncols,
nullptr);
394 mzd_bool_mul_m4rm(temp, A_HV, B_VV, k,blocksize,threads);
396 mzd_optimal_parameters(temp, A_VH,
397 ConcurrencyContext::cacheSize(2),
398 ConcurrencyContext::cacheSize(3),
401 mzd_t* C_HH =
new_mzd(A_HV->nrows,A_VH->ncols,pool);
403 Q_ASSERT(C_HH->nrows==A_HV->nrows);
404 Q_ASSERT(C_HH->ncols==A_VH->ncols);
405 mzd_bool_mul_m4rm(C_HH, temp, A_VH, k,blocksize,threads);
408 if (! mzd_is_zero(C_HH)) {
430 int threads = ConcurrencyContext::countThreads();
431 mzd_bool_and_uptri(A_VV, A_VV, B_VV, threads);
433 if (mzd_is_zero(A_VV))
443 for (
int k = 0; k <
std::min(M->ncols,M->nrows); ++k)
444 if (mzd_read_bit(M, k, k)) {
464 mzd_t* M =
mtx[ori1][ori2];
465 if (!M)
return false;
467 if (!
mask[ori1].contains(i1) || !
mask[ori2].contains(i2))
472 return mzd_read_bit(M,row,col);
478 mzd_t* M =
mtx[qfrom.ori][qto.ori];
479 if (!M)
return false;
481 qfrom &=
mask[qfrom.ori];
482 qto &=
mask[qto.ori];
487 Q_ASSERT(qfrom.lower >= 0 && qfrom.upper <=
mask[qfrom.ori].
len());
488 Q_ASSERT(qto.lower >= 0 && qto.upper <=
mask[qto.ori].
len());
490 for(
int row=qfrom.lower; row < qfrom.upper; ++row)
491 if (
has_bits(M,row,qto.lower,qto.upper))
498 return !
mtx[o1][o2];
503 return mtx[o1][o2] && mzd_is_zero(
mtx[o1][o2]);
533 mzd_t* M =
mtx[ri.ori][rj.
ori];
540 for(
int i=ri.lower; i < ri.upper; ++i)
550 for(
int i=ri.lower; i < ri.upper; ++i)
558 Q_ASSERT(i >= 0 && i <= M->ncols);
559 Q_ASSERT(j >= i && j <= M->ncols);
560 Q_ASSERT(r >= 0 && r < M->nrows);
562 wi_t w1 = i / m4ri_radix;
563 wi_t w2 = j / m4ri_radix;
565 word mask1 = ((word)-1) << (i % m4ri_radix);
566 word mask2 = (((word)1) << (j % m4ri_radix)) - 1;
567 word* row = M->rows[r];
568 Q_ASSERT(mzd_row(M,r)==row);
569 Q_ASSERT(w1>=0 && w1<M->width);
570 Q_ASSERT(w2>=0 && w2<M->width);
573 row[w1] |= mask1 & mask2;
580 for(wi_t w=w1+1; w < w2; ++w)
587 Q_ASSERT(i >= 0 && i <= M->ncols);
588 Q_ASSERT(j >= i && j <= M->ncols);
589 Q_ASSERT(r >= 0 && r < M->nrows);
591 wi_t w1 = i / m4ri_radix;
592 wi_t w2 = j / m4ri_radix;
594 word mask1 = ((word)-1) << (i % m4ri_radix);
595 word mask2 = (((word)1) << (j % m4ri_radix)) - 1;
596 word* row = M->rows[r];
597 Q_ASSERT(mzd_row(M,r)==row);
600 return (row[w1] & mask1 & mask2);
602 if (row[w1] & mask1)
return true;
603 if (row[w2] & mask2)
return true;
605 for(wi_t w=w1+1; w < w2; ++w)
606 if (row[w])
return true;
614 mzd_t* m =
mtx[ori1][ori2];
616 return m->nrows * m->rowstride *
sizeof(word);
621 mzd_t* m =
mtx[ori1][ori2];
623 return mzd_density(m,0);
642 mzd_t* m =
mtx[o1][o2];
643 result +=
density(o1,o2) * m->nrows*m->ncols;
646 return result / (nrows*nrows);
659 mzd_set_ui(
mtx[o1][o2],0);
669 if (m->nrows != m->ncols)
672 for(
int i=1; i < m->nrows; ++i)
673 for (
int j=0; j < i; ++j)
674 if (mzd_read_bit(m,i,j))
681 mzd_t* M =
mtx[o1][o2];
682 if (!M)
return INT_MAX;
684 if ((i >=
mask[o1].upper) || (j >=
mask[o2].upper))
697 mzd_t* M =
mtx[o1][o2];
700 if (!
mask[o1].contains(i) || !
mask[o2].contains(j))
716 Q_ASSERT(j2 <=
mask[o2].upper);
721 #define INVERT(w) ((w) ^ ((uint64_t)-1)) 730 return boost::multiprecision::lsb(x);
735 Q_ASSERT(r>=0 && r<M->nrows);
738 word* row = M->rows[r];
739 wi_t w = i / m4ri_radix;
743 word v = row[w] >> (i % m4ri_radix);
747 while(++w < M->width)
749 return w*m4ri_radix +
lsb_pos(row[w]);
756 Q_ASSERT(r>=0 && r<M->nrows);
759 word* row = M->rows[r];
760 wi_t w = i / m4ri_radix;
764 word v = (
INVERT(row[w]) >> (i % m4ri_radix));
768 while(++w < M->width)
769 if (row[w]!=(uint64_t)-1)
782 const mzd_t* M =
mtx[o1][o2];
788 for(
int row=0; row < M->nrows; ++row) {
789 for (
int col = 0; col < M->ncols; col += 64) {
790 word w = mzd_read_bits(M, row, col, 64);
791 out << std::hex << std::setw(16) << w <<
" ";
801 mtx[o1][o2] =
nullptr;
int find_next0(Orientation o1, int i, Orientation o2, int j) const
find the next missing edge in a sub-graph
bool contains(int i) const
containment test
void swap(gpuword **A, gpuword **B)
int count2(Orientation ori) const
void print(std::ostream &out, Orientation o1, Orientation o2) const
print a sub-adjacancy matrix
Graph::ptr P
valid placement applied to the graph
mzd_t * allocate(Orientation o1, Orientation o2)
allocate an adjacancy sub-matrix (one of four parts)
int i
if this is a reachability graph: column range in free-space. this = RG(i,j), or CRG(i,...
virtual void combine(const Graph *P)
apply the COMBINE operation, filtering edges with valid placements. Effectively performs a Boolean AN...
global definitions for all algorithms.
mzd_t * mtx[2][2]
adjacency matrix (M4RI structure) split into four parts to allow for memory savings.
const Orientation ori
horizontal or vertical
non-accesible. no interval is reachable from this one
int map_lower(Pointer p) const
map the lower bound of a reachability structure interval
boost::shared_ptr< GraphModel > ptr
smart pointer to a GraphModel object
int diagonalElement
result of call to searchDiagonalElement
The StructureIterator class.
void releaseIfZero()
release memory that is not neeeded (empty sub-graphs)
bool zero(Orientation o1, Orientation o2) const
empty test for one of the four sub-graphs
void setOriginMerge3(Graph::ptr A, Graph::ptr B)
this graph is the transitive closure of two other graphs. It is the result of a final MERGE operation...
boundary interval in the reachability structure. Represents an interval on the boundary of the FreeSp...
mzd_t * new_mzd(int rows, int cols, MatrixPool *pool)
allocate a new mzd_t structure (a matrix for the M4RI algorithms)
bool empty(Orientation o1, Orientation o2) const
empty test for one of the four sub-graphs
bool contains_edge(Orientation ori_from, int i_from, Orientation ori_to, int i_to) const
check if an edge is present
bool is_upper_triangular(Orientation ori1, Orientation ori2) const
int upper
upper index (exclusive)
int lower
lower index (inclusive)
GraphPtr newGraph(const GraphModel::ptr model)
Type type
reachability label
virtual void release()
release memory for all parts of the adjacancy matrix
a very simple Rectangle structure, with integer boundaries.
Graph::ptr merge2(const Graph *B, MatrixPool *pool=nullptr) const
apply the MERGE operation, computing the transitive closure of two graphs. Effectively performs a mat...
Graph::ptr merge3(const Graph *B, MatrixPool *pool=nullptr) const
apply the final MERGE operation, computing the transitive closure of two graphs. Effectively performs...
void allocateAll()
allocate data for adjacenty matrixes (all four parts)
void setOriginRG(int i)
this graph represents a reachability graph, it was created from a reachability structure
void setOriginPlacement(int di, int dj)
this graph represents a set of valid placements
a range of node indices in a Reachability Graph
void setOriginCombine(Graph::ptr P)
a sub-set of valid placements is applied to the graph
Orientation
Segment Orientation.
IndexRange map(Orientation ori, const Interval &, bool upperBoundInclusive) const
map a free-space interval to a range of graph nodes
boost::shared_ptr< Graph > ptr
this is a rechability graph, constructed from a reachability structure
memory pool for matrix objects (M4RI matrices mzd_t* and OpenCL matrices clm4rm_t*)
double min(double a, double b)
minimum function with checks for NAN
void reclaim(mzd_t *m, MatrixPool *pool)
reclaim an object (i.e. put it into the recycling list)
Graph(const GraphModel &model, graph_t *ag)
Represents a Reachability Graph. Vertices correspond to intervals in the reachability structure,...
void add_range(IndexRange r_from, IndexRange r_to)
add a set of edges to the graph, covering a range of nodes
this graph is the result of a MERGE operation (transitive closure) on two graphs
IndexRange mergedHMask(const IndexRange &m1, const IndexRange &m2) const
merged two node ranges representing the horizontal masks of reachability graph (see frechet::reach::G...
Orientation ori
horizontal or vertical part of the reachability structure?
IndexRange find_next_range(Orientation o1, int i, Orientation o2, int j) const
find a range of consecutive edges in a sub-graph
void setOriginMerge2(Graph::ptr A, Graph::ptr B)
this graph is the transitive closure of two other graphs. It is the result of a MERGE operation.
void add_edge(int i, int j)
this graph represents a set of valid placements
struct frechet::reach::Graph::Origin origin
int find_next1(Orientation o1, int i, Orientation o2, int j) const
find the next edge in a sub-graph
static bool has_bits(mzd_t *M, int row, int i, int j)
query bits in a Boolean matrix
Graph::ptr A
if this graph was constructed as transitive closure: original graphs this = MERGE(A,...
int blevel
bottom-level (for topologial sorting; leaves = blevel 0
virtual void queryDiagonalElement() const
find an edge on the diagonal of the adjacancy matrix. Does not return a result. To query the result o...
static void set_bits(mzd_t *M, int row, int i, int j)
set some bits in a Boolean matrix
const GraphModel::ptr graphModel() const
The Reachability Structure; maintains a list of intervals on the border of Free Space,...
const IndexRange & hmask() const
horizontal range covered by this graph
int succesors
number of successors
void printOrigin(std::ostream &out) const
debug output
virtual void finalize()
release memory that is not neeeded (empty sub-graphs)
Rect rect(Orientation o1, Orientation o2) const
range covered by this graph
static bool is_adjacent_to(const Graph &A, const Graph &B)
Graph & operator=(const Graph &)
this graph is the result of a MERGE operation (transitive closure) on three graphs
double max(double a, double b)
maximum function with checks for NAN
enum frechet::reach::Graph::Origin::@7 operation
virtual int foundDiagonalElement() const
bool operator==(const Graph &that) const
comparison operator