graph_cl.cpp

Go to the documentation of this file.

#include <graph_cl.h>
#include <concurrency.h>

using namespace frechet;
using namespace reach;
using namespace app;

/*
 *  Factory Methods
 */
GraphPtr frechet::reach::newGraph(const GraphModel::ptr model)
{
    if (ConcurrencyContext::hasGpuSupport())
        return GraphPtr(new GraphCL(model));
    else
        return GraphPtr(new Graph(model));
}

GraphPtr frechet::reach::newGraph(const GraphModel::ptr model, IndexRange hmask)
{
    if (ConcurrencyContext::hasGpuSupport())
        return GraphPtr(new GraphCL(model,hmask));
    else
        return GraphPtr(new Graph(model, hmask));
}

GraphPtr frechet::reach::newGraph(const GraphModel::ptr model, Structure& str)
{
    if (ConcurrencyContext::hasGpuSupport())
        return GraphPtr(new GraphCL(model,str));
    else
        return GraphPtr(new Graph(model, str));
}


GraphCL::GraphCL(const GraphModel::ptr model)
    : Graph(model),
    clmtx{ {nullptr,nullptr},{nullptr,nullptr} },
    temps(), cond(), diagonalElementBuffer(nullptr)
{
    init_conditions(&cond);
}

GraphCL::GraphCL(const GraphModel::ptr model, IndexRange hmask)
    : Graph(model, hmask),
    clmtx{ { nullptr,nullptr },{ nullptr,nullptr } },
    temps(), cond(), diagonalElementBuffer(nullptr)
{
    init_conditions(&cond);
}

GraphCL::GraphCL(const GraphModel::ptr model, Structure& str)
    : Graph(model, str),
    clmtx{ {nullptr, nullptr}, { nullptr,nullptr } },
    temps(), cond(), diagonalElementBuffer(nullptr)
{
    init_conditions(&cond);
}

GraphCL::GraphCL(const GraphCL& that)
    : Graph(that)
{ 
    copy(that);
}

GraphCL::GraphCL(GraphCL&& that)
    : Graph(that)
{ 
    swap(that);
}

GraphCL::~GraphCL()
{
    for (Orientation o1 = HORIZONTAL; o1 <= VERTICAL; ++o1)
        for (Orientation o2 = HORIZONTAL; o2 <= VERTICAL; ++o2) {
            clm4rm_free(clmtx[o1][o2]);
            clmtx[o1][o2] = nullptr;
        }
    for (clmatrix_t* t : temps)
        clm4rm_free(t);
    temps.clear();

    if (diagonalElementBuffer)
        clReleaseMemObject(diagonalElementBuffer);
    release_conditions(&cond);
}

GraphCL& GraphCL::operator= (const GraphCL& that) {
    Graph::operator=(that);
    copy(that);
    return *this;
}

GraphCL& GraphCL::operator= (GraphCL&& that) {
    Graph::operator=(that);
    swap(that);
    return *this;
}

void GraphCL::copy(const GraphCL& that) 
{
    for (Orientation o1 = HORIZONTAL; o1 <= VERTICAL; ++o1)
        for (Orientation o2 = HORIZONTAL; o2 <= VERTICAL; ++o2)
            if (that.clmtx[o1][o2])
                clmtx[o1][o2] = clm4rm_copy(that.mtx[o1][o2], 32, true, ConcurrencyContext::clContext());
            else
                clmtx[o1][o2] = nullptr;
    //  don't copy diagonalElementBuffer
}

void GraphCL::swap(GraphCL& that)
{
    //  steal from 'that'
    for (Orientation o1 = HORIZONTAL; o1 <= VERTICAL; ++o1)
        for (Orientation o2 = HORIZONTAL; o2 <= VERTICAL; ++o2)
            std::swap(clmtx[o1][o2], that.clmtx[o1][o2]);
    std::swap(this->diagonalElementBuffer, that.diagonalElementBuffer);
}

void GraphCL::release(Orientation o1, Orientation o2)
{
    Graph::release(o1, o2);
    if (clmtx[o1][o2]) {
        clm4rm_free(clmtx[o1][o2]);
        clmtx[o1][o2] = nullptr;
    }
}

void GraphCL::finalize() {
    Graph::finalize();    //    releaseIfZero
    synchToGpu();
}
void GraphCL::resetConditions()
{
    release_conditions(&cond);
}

void GraphCL::synchToGpu() {
    //  copy to GPU memory
    for (Orientation o1 = HORIZONTAL; o1 <= VERTICAL; ++o1)
        for (Orientation o2 = HORIZONTAL; o2 <= VERTICAL; ++o2)
        {
            mzd_t* M = mtx[o1][o2];
            clmatrix_t* G = clmtx[o1][o2];

            if (!M && !G)
                continue;

            if (G) {
                clm4rm_free(G);
                G = clmtx[o1][o2] = nullptr;
            }

            if (M) {
                G = clmtx[o1][o2] = clm4rm_create(M->nrows, M->ncols, 32,
                                        true, ConcurrencyContext::clContext());
                clm4rm_write(G, M, ConcurrencyContext::clQueue(), &cond);
                //  Note: clm4rm_copy does the same but is always _blocking_
                join_conditions(&cond);
            }
        }
}

void GraphCL::synchFromGpu() {
    //  copy to GPU memory
    for (Orientation o1 = HORIZONTAL; o1 <= VERTICAL; ++o1)
        for (Orientation o2 = HORIZONTAL; o2 <= VERTICAL; ++o2)
        {
            if (!mtx[o1][o2] && !clmtx[o1][o2])
                continue;

            if (clmtx[o1][o2]) {
                mtx[o1][o2] = clm4rm_read(mtx[o1][o2], clmtx[o1][o2],
                                          ConcurrencyContext::clQueue(), &cond);
                join_conditions(&cond);
            }
        }
}

//  (*this)_VV &= that_VV
void GraphCL::combine(const Graph* Bg)
{
    const GraphCL* B = dynamic_cast<const GraphCL*>(Bg);
    Q_ASSERT(B);

    clmatrix_t* A_VV = this->clmtx[VERTICAL][VERTICAL];
    clmatrix_t* B_VV = B->clmtx[VERTICAL][VERTICAL];

    Q_ASSERT((this->mtx[VERTICAL][VERTICAL]==nullptr) || (A_VV!=nullptr));
    Q_ASSERT((B->mtx[VERTICAL][VERTICAL]==nullptr) || (B_VV!=nullptr));

    if (!A_VV) return;
    if (!B_VV) {
        release(VERTICAL, VERTICAL);
        return;
    }

    merge_conditions(&cond,&B->cond);
#if IMAGE2D
        //  can't do clm4rm_and in-place. Use temp matrix instead:
        clmatrix_t* T = tempMatrix(A_VV->nrows,A_VV->ncols);
        clm4rm_and(T, A_VV, B_VV, ConcurrencyContext::clQueue(),&cond);
        std::swap(A_VV->data,T->data);
        //  don't release T. Everything is (will be) asynchronous.
#else
        //  VV are upper triangular matrices
        clm4rm_and(A_VV, A_VV, B_VV, ConcurrencyContext::clQueue(),&cond);
        //  TODO clutri_and
#endif
    join_conditions(&cond);
}

void merge2(const Graph* A, const Graph* B, MatrixPool* pool);

void merge3(const Graph* A, const Graph* B, MatrixPool* pool);

//  @return (*this)_VV * that_VV
void GraphCL::merge2(const Graph* Ag, const Graph* Bg, MatrixPool* pool)
{
    const GraphCL* A = dynamic_cast<const GraphCL*>(Ag);
    const GraphCL* B = dynamic_cast<const GraphCL*>(Bg);
    Q_ASSERT(A && B);
    Q_ASSERT(is_adjacent_to(*A, *B));

    GraphCL* C = this;

    clmatrix_t* A_VV = A->clmtx[VERTICAL][VERTICAL];
    clmatrix_t* B_VV = B->clmtx[VERTICAL][VERTICAL];

    Q_ASSERT((A->mtx[VERTICAL][VERTICAL]==nullptr) || (A_VV!=nullptr));
    Q_ASSERT((B->mtx[VERTICAL][VERTICAL]==nullptr) || (B_VV!=nullptr));

    if (!A_VV || !B_VV)
        return;

    merge_conditions(&cond,&A->cond);
    merge_conditions(&cond,&B->cond);

    //  both matrices are upper triangular
    size2_t max_tile;
    ConcurrencyContext::maxMaxtrixTile(max_tile);
    clmatrix_t* C_VV = new_clmatrix(A_VV->nrows,B_VV->ncols,pool,&cond);
    join_conditions(&cond);

    clutri_mul(C_VV, A_VV, B_VV,  max_tile, ConcurrencyContext::clQueue(), &cond);
    join_conditions(&cond);

    this->clmtx[VERTICAL][VERTICAL] = C_VV;
    Q_ASSERT(C_VV->nrows == this->mask[VERTICAL].len());
    Q_ASSERT(C_VV->ncols == this->mask[VERTICAL].len());
}

//  @return (*this)_HV * that_VV * (*this)_VH
void GraphCL::merge3(const Graph* Ag, const Graph* Bg, MatrixPool* pool)
{
    const GraphCL* A = dynamic_cast<const GraphCL*>(Ag);
    const GraphCL* B = dynamic_cast<const GraphCL*>(Bg);
    Q_ASSERT(A && B);
    Q_ASSERT(is_adjacent_to(*A, *B));
    Q_ASSERT(is_adjacent_to(*B, *A));

    clmatrix_t* A_HV = A->clmtx[HORIZONTAL][VERTICAL];
    clmatrix_t* B_VV = B->clmtx[VERTICAL][VERTICAL];
    clmatrix_t* A_VH = A->clmtx[VERTICAL][HORIZONTAL];

    Q_ASSERT((A_HV!=nullptr)==(A->mtx[HORIZONTAL][VERTICAL]!=nullptr));
    Q_ASSERT((B->mtx[VERTICAL][VERTICAL]==nullptr) || (B_VV!=nullptr));
    Q_ASSERT((A_VH!=nullptr)==(A->mtx[VERTICAL][HORIZONTAL]!=nullptr));

    if (!A_HV || !B_VV || !A_VH)
        return;   //  shortcut for empty

    merge_conditions(&cond,&A->cond);
    merge_conditions(&cond,&B->cond);

    size2_t max_tile;
    ConcurrencyContext::maxMaxtrixTile(max_tile);
    clmatrix_t* temp = tempMatrix(A_HV->nrows,B_VV->ncols,pool);
    clmatrix_t* C_HH = new_clmatrix(temp->nrows,A_VH->ncols,pool,&cond);
    join_conditions(&cond);

    clcubic_mul(temp, A_HV, B_VV, max_tile, ConcurrencyContext::clQueue(), &cond);
    join_conditions(&cond);

    clcubic_mul(C_HH, temp, A_VH,  max_tile, ConcurrencyContext::clQueue(), &cond);
    join_conditions(&cond);

    this->clmtx[HORIZONTAL][HORIZONTAL] = C_HH;
    Q_ASSERT(C_HH->nrows == this->mask[HORIZONTAL].len());
    Q_ASSERT(C_HH->ncols == this->mask[HORIZONTAL].len());
}

void GraphCL::queryDiagonalElement() const
{
    diagonalElement=-1;
    clmatrix_t* M = clmtx[HORIZONTAL][HORIZONTAL];
    if (!M) return;

    Q_ASSERT(diagonalElementBuffer==nullptr);

    //  asynchronous call. returns buffer but no result
    diagonalElementBuffer = clm4rm_query_diagonal(M, ConcurrencyContext::clContext(), ConcurrencyContext::clQueue(), &cond);
    join_conditions(&cond);
}

int GraphCL::foundDiagonalElement() const
{
    if(diagonalElementBuffer) {
        //  copy result from GPU; blocking callAda!
        int i = clm4rm_query_result(diagonalElementBuffer, ConcurrencyContext::clQueue(), &cond);
        if (i >= 0)
            diagonalElement = i + mask[HORIZONTAL].lower;
        else
            diagonalElement = -1;

        diagonalElementBuffer=nullptr; // was released by clm4rm_query_result
    }
    return diagonalElement;
}


clmatrix_t* GraphCL::tempMatrix(int rows, int cols, MatrixPool* pool) const
{
    //  odd-sized temp matrices are not pooled
    clmatrix_t* t = new_clmatrix(rows, cols, pool, &cond);
    temps.push_back(t);
    return t;
}

void GraphCL::release(Orientation o1, Orientation o2, MatrixPool* pool) {
    Graph::release(o1,o2,pool);
    if (clmtx[o1][o2]) {
        reclaim(clmtx[o1][o2],pool);
        clmtx[o1][o2] = nullptr;
    }
    //  release temp matrices to pool
    if (pool) {
        for (clmatrix_t *T : temps)
            reclaim(T, pool);
        temps.clear();
    }
}