sequential_dgetrf.c

Go to the documentation of this file.

/**
 *
 * @file sequential_dgetrf.c
 *
 * @copyright 2012-2023 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
 *                      Univ. Bordeaux. All rights reserved.
 *
 * @version 6.3.2
 * @author Pascal Henon
 * @author Xavier Lacoste
 * @author Pierre Ramet
 * @author Mathieu Faverge
 * @author Esragul Korkmaz
 * @author Gregoire Pichon
 * @author Tony Delarue
 * @date 2023-10-23
 *
 * @generated from /builds/solverstack/pastix/sopalin/sequential_zgetrf.c, normal z -> d, Wed Dec 13 12:09:47 2023
 *
 **/
#include "common.h"
#include "isched.h"
#include "blend/solver.h"
#include "sopalin/sopalin_data.h"
#include "sopalin/coeftab_d.h"
#include "pastix_dcores.h"
 
#if defined(PASTIX_WITH_PARSEC)
#include "parsec/pastix_dparsec.h"
#endif
 
#if defined(PASTIX_WITH_STARPU)
#include "starpu/pastix_dstarpu.h"
#endif
 
/**
 *******************************************************************************
 *
 * @brief TODO
 *
 *******************************************************************************
 *
 * @param[in] pastix_data
 *          TODO
 *
 * @param[in] sopalin_data
 *          TODO
 *
 *******************************************************************************/
void
sequential_dgetrf( pastix_data_t  *pastix_data,
                   sopalin_data_t *sopalin_data )
{
    SolverMatrix       *datacode = pastix_data->solvmatr;
    SolverCblk         *cblk;
    double *work;
    pastix_int_t  i, lwork;
    (void)sopalin_data;
 
    lwork = datacode->gemmmax;
    if ( (datacode->lowrank.compress_when != PastixCompressNever) &&
         (datacode->lowrank.ilu_lvl < INT_MAX) )
    {
        lwork = pastix_imax( lwork, 2 * datacode->blokmax );
    }
    MALLOC_INTERN( work, lwork, double );
 
    cblk = datacode->cblktab;
    for (i=0; i<datacode->cblknbr; i++, cblk++){
        if ( cblk->cblktype & CBLK_IN_SCHUR ) {
            break;
        }
 
        /* Wait for incoming dependencies */
        if ( cpucblk_dincoming_deps( 0, PastixLUCoef,
                                     datacode, cblk ) )
        {
            continue;
        }
 
        /* Compute */
        cpucblk_dgetrfsp1d( datacode, cblk,
                            work, lwork );
    }
 
    memFree_null( work );
}
 
/**
 *******************************************************************************
 *
 * @brief TODO
 *
 *******************************************************************************
 *
 * @param[in] ctx
 *          TODO
 *
 * @param[in] args
 *          TODO
 *
 *******************************************************************************/
void
thread_dgetrf_static( isched_thread_t *ctx,
                      void            *args )
{
    sopalin_data_t     *sopalin_data = (sopalin_data_t*)args;
    SolverMatrix       *datacode = sopalin_data->solvmtx;
    SolverCblk         *cblk;
    Task               *t;
    double *work;
    pastix_int_t i, ii, lwork;
    pastix_int_t tasknbr, *tasktab;
    int rank = ctx->rank;
 
    lwork = datacode->gemmmax;
    if ( (datacode->lowrank.compress_when != PastixCompressNever) &&
         (datacode->lowrank.ilu_lvl < INT_MAX) )
    {
        lwork = pastix_imax( lwork, 2 * datacode->blokmax );
    }
    MALLOC_INTERN( work, lwork, double );
 
    tasknbr = datacode->ttsknbr[rank];
    tasktab = datacode->ttsktab[rank];
 
    for (ii=0; ii<tasknbr; ii++) {
        i = tasktab[ii];
        t = datacode->tasktab + i;
        cblk = datacode->cblktab + t->cblknum;
 
        if ( cblk->cblktype & CBLK_IN_SCHUR ) {
            continue;
        }
 
        /* Wait for incoming dependencies */
        if ( cpucblk_dincoming_deps( rank, PastixLUCoef,
                                     datacode, cblk ) )
        {
            continue;
        }
 
        /* Compute */
        cpucblk_dgetrfsp1d( datacode, cblk,
                            work, lwork );
    }
 
    memFree_null( work );
}
 
/**
 *******************************************************************************
 *
 * @brief TODO
 *
 *******************************************************************************
 *
 * @param[in] pastix_data
 *          TODO
 *
 * @param[in] sopalin_data
 *          TODO
 *
 *******************************************************************************/
void
static_dgetrf( pastix_data_t  *pastix_data,
               sopalin_data_t *sopalin_data )
{
    isched_parallel_call( pastix_data->isched, thread_dgetrf_static, sopalin_data );
}
 
/**
 * @brief TODO
 */
struct args_dgetrf_t
{
    sopalin_data_t   *sopalin_data;
    volatile int32_t  taskcnt;
};
 
/**
 *******************************************************************************
 *
 * @brief TODO
 *
 *******************************************************************************
 *
 * @param[in] ctx
 *          TODO
 *
 * @param[in] args
 *          TODO
 *
 *******************************************************************************/
void
thread_dgetrf_dynamic( isched_thread_t *ctx,
                       void            *args )
{
    struct args_dgetrf_t *arg = (struct args_dgetrf_t*)args;
    sopalin_data_t       *sopalin_data = arg->sopalin_data;
    SolverMatrix         *datacode = sopalin_data->solvmtx;
    SolverCblk           *cblk;
    SolverBlok           *blok;
    Task                 *t;
    pastix_queue_t       *computeQueue;
    double   *work;
    pastix_int_t          i, ii, lwork;
    pastix_int_t          tasknbr, *tasktab, cblknum, bloknum;
    int32_t               local_taskcnt = 0;
    int                   rank = ctx->rank;
 
    lwork = datacode->gemmmax;
    if ( (datacode->lowrank.compress_when != PastixCompressNever) &&
         (datacode->lowrank.ilu_lvl < INT_MAX) )
    {
        lwork = pastix_imax( lwork, 2 * datacode->blokmax );
    }
    MALLOC_INTERN( work, lwork, double );
    MALLOC_INTERN( datacode->computeQueue[rank], 1, pastix_queue_t );
 
    tasknbr      = datacode->ttsknbr[rank];
    tasktab      = datacode->ttsktab[rank];
    computeQueue = datacode->computeQueue[rank];
    pqueueInit( computeQueue, tasknbr );
 
    /* Initialize the local task queue with available cblks */
    for (ii=0; ii<tasknbr; ii++) {
        i = tasktab[ii];
        t = datacode->tasktab + i;
 
        if ( !(t->ctrbcnt) ) {
            cblk = datacode->cblktab + t->cblknum;
            pqueuePush1( computeQueue, t->cblknum, cblk->priority );
        }
    }
 
    /* Make sure that all computeQueues are allocated */
    isched_barrier_wait( &(ctx->global_ctx->barrier) );
 
    while( arg->taskcnt > 0 )
    {
        cblknum = pqueuePop( computeQueue );
 
#if defined(PASTIX_WITH_MPI)
        /* Nothing to do, let's make progress on communications */
        if( cblknum == -1 ) {
            cpucblk_dmpi_progress( PastixLUCoef, datacode, rank );
            cblknum = pqueuePop( computeQueue );
        }
#endif
 
        /* No more local job, let's steal our neighbors */
        if( cblknum == -1 ) {
            if ( local_taskcnt ) {
                pastix_atomic_sub_32b( &(arg->taskcnt), local_taskcnt );
                local_taskcnt = 0;
            }
            cblknum = stealQueue( datacode, rank,
                                  ctx->global_ctx->world_size );
        }
 
        /* Still no job, let's loop again */
        if ( cblknum == -1 ) {
            continue;
        }
 
        if ( cblknum >= 0 ) {
            cblk = datacode->cblktab + cblknum;
            if ( cblk->cblktype & CBLK_IN_SCHUR ) {
                continue;
            }
            cblk->threadid = rank;
 
            /* Compute */
            if ( cblk->cblktype & CBLK_TASKS_2D ) {
                cpucblk_dgetrfsp1dplus( datacode, cblk );
            }
            else {
                cpucblk_dgetrfsp1d( datacode, cblk,
                                    work, lwork );
            }
        }
        else {
            bloknum = - cblknum - 1;
            blok    = datacode->bloktab + bloknum;
            cpucblk_dgetrfsp1dplus_update( datacode, blok, work, lwork );
        }
        local_taskcnt++;
    }
    memFree_null( work );
 
    /* Make sure that everyone is done before freeing */
    isched_barrier_wait( &(ctx->global_ctx->barrier) );
    pqueueExit( computeQueue );
    memFree_null( computeQueue );
}
 
/**
 *******************************************************************************
 *
 * @brief TODO
 *
 *******************************************************************************
 *
 * @param[in] pastix_data
 *          TODO
 *
 * @param[in] sopalin_data
 *          TODO
 *
 *******************************************************************************/
void
dynamic_dgetrf( pastix_data_t  *pastix_data,
                sopalin_data_t *sopalin_data )
{
    SolverMatrix         *datacode    = sopalin_data->solvmtx;
    int32_t               taskcnt     = datacode->tasknbr_1dp;
    struct args_dgetrf_t  args_dgetrf = { sopalin_data, taskcnt };
 
    /* Allocate the computeQueue */
    MALLOC_INTERN( datacode->computeQueue,
                   pastix_data->isched->world_size, pastix_queue_t * );
 
    isched_parallel_call( pastix_data->isched, thread_dgetrf_dynamic, &args_dgetrf );
 
    memFree_null( datacode->computeQueue );
}
 
#ifndef DOXYGEN_SHOULD_SKIP_THIS
static void (*dgetrf_table[5])(pastix_data_t *, sopalin_data_t *) = {
    sequential_dgetrf,
    static_dgetrf,
#if defined(PASTIX_WITH_PARSEC)
    parsec_dgetrf,
#else
    NULL,
#endif
#if defined(PASTIX_WITH_STARPU)
    starpu_dgetrf,
#else
    NULL,
#endif
    dynamic_dgetrf
};
#endif /* DOXYGEN_SHOULD_SKIP_THIS */
 
/**
 *******************************************************************************
 *
 * @brief TODO
 *
 *******************************************************************************
 *
 * @param[in] pastix_data
 *          TODO
 *
 * @param[in] sopalin_data
 *          TODO
 *
 *******************************************************************************/
void
sopalin_dgetrf( pastix_data_t  *pastix_data,
                sopalin_data_t *sopalin_data )
{
    int sched = pastix_data->iparm[IPARM_SCHEDULER];
    void (*dgetrf)(pastix_data_t *, sopalin_data_t *) = dgetrf_table[ sched ];
 
    if (dgetrf == NULL) {
        sched  = PastixSchedDynamic;
        dgetrf = dynamic_dgetrf;
    }
 
    if ( (sched == PastixSchedSequential) ||
         (sched == PastixSchedStatic)     ||
         (sched == PastixSchedDynamic) )
    {
        solverRequestInit( PastixFacto, sopalin_data->solvmtx );
        solverRecvInit( PastixLUCoef, sopalin_data->solvmtx, PastixDouble );
    }
 
    dgetrf( pastix_data, sopalin_data );
 
    if ( (sched == PastixSchedSequential) ||
         (sched == PastixSchedStatic)     ||
         (sched == PastixSchedDynamic) )
    {
        cpucblk_drequest_cleanup( PastixLUCoef, sched, sopalin_data->solvmtx );
        solverRequestExit( sopalin_data->solvmtx );
        solverRecvExit( sopalin_data->solvmtx );
    }
 
#if defined(PASTIX_DEBUG_FACTO)
    coeftab_ddump( pastix_data, sopalin_data->solvmtx, "getrf" );
#endif
}