Commit 86ee80ff authored by Mathieu Faverge's avatar Mathieu Faverge
Browse files

Update the 3 other amigos to be able to run throught parsec or sequential

parent 2535591d
......@@ -25,28 +25,45 @@
#include "d_bcsc.h"
#include "s_bcsc.h"
#if defined(PASTIX_DEBUG_SOLVE)
static inline void dump_rhs( char *name, int n, double *b )
{
int i;
fprintf(stderr,"%s :", name );
for (i=0; i<n; i++) {
if (i%10 == 0)
fprintf(stderr, "\n");
fprintf(stderr,"%e ", b[i]);
}
fprintf(stderr,"\n");
}
#else
#define dump_rhs(...) do {} while(0)
#endif
void
pastix_trsm( int flttype, int side, int uplo, int trans, int diag,
pastix_trsm( pastix_data_t *pastix_data,
int flttype, int side, int uplo, int trans, int diag,
sopalin_data_t *sopalin_data,
int nrhs, void *b, int ldb )
{
switch (flttype) {
case PastixComplex64:
sequential_ztrsm( side, uplo, trans, diag,
sequential_ztrsm( pastix_data, side, uplo, trans, diag,
sopalin_data, nrhs, (pastix_complex64_t *)b, ldb );
break;
case PastixComplex32:
sequential_ctrsm( side, uplo, trans, diag,
sequential_ctrsm( pastix_data, side, uplo, trans, diag,
sopalin_data, nrhs, (pastix_complex32_t *)b, ldb );
break;
case PastixDouble:
trans = (trans == PastixConjTrans) ? PastixTrans : trans;
sequential_dtrsm( side, uplo, trans, diag,
sequential_dtrsm( pastix_data, side, uplo, trans, diag,
sopalin_data, nrhs, (double *)b, ldb );
break;
case PastixFloat:
trans = (trans == PastixConjTrans) ? PastixTrans : trans;
sequential_strsm( side, uplo, trans, diag,
sequential_strsm( pastix_data, side, uplo, trans, diag,
sopalin_data, nrhs, (float *)b, ldb );
break;
default:
......@@ -55,22 +72,22 @@ pastix_trsm( int flttype, int side, int uplo, int trans, int diag,
}
void
pastix_diag( int flttype,
pastix_diag( pastix_data_t *pastix_data, int flttype,
sopalin_data_t *sopalin_data,
int nrhs, void *b, int ldb )
{
switch (flttype) {
case PastixComplex64:
sequential_zdiag( sopalin_data, nrhs, (pastix_complex64_t *)b, ldb );
sequential_zdiag( pastix_data, sopalin_data, nrhs, (pastix_complex64_t *)b, ldb );
break;
case PastixComplex32:
sequential_cdiag( sopalin_data, nrhs, (pastix_complex32_t *)b, ldb );
sequential_cdiag( pastix_data, sopalin_data, nrhs, (pastix_complex32_t *)b, ldb );
break;
case PastixDouble:
sequential_ddiag( sopalin_data, nrhs, (double *)b, ldb );
sequential_ddiag( pastix_data, sopalin_data, nrhs, (double *)b, ldb );
break;
case PastixFloat:
sequential_sdiag( sopalin_data, nrhs, (float *)b, ldb );
sequential_sdiag( pastix_data, sopalin_data, nrhs, (float *)b, ldb );
break;
default:
fprintf(stderr, "Unknown floating point arithmetic\n" );
......@@ -176,26 +193,33 @@ pastix_task_solve( pastix_data_t *pastix_data,
clockStart(timer);
switch ( pastix_data->iparm[IPARM_FACTORIZATION] ){
case PastixFactLLT:
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixNonUnit, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixConjTrans, PastixNonUnit, &sopalin_data, nrhs, b, ldb );
dump_rhs( "AfterPerm", csc->gN, b );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixNonUnit, &sopalin_data, nrhs, b, ldb );
dump_rhs( "AfterDown", csc->gN, b );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixConjTrans, PastixNonUnit, &sopalin_data, nrhs, b, ldb );
dump_rhs( "AfterUp", csc->gN, b );
break;
case PastixFactLDLT:
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
pastix_diag( pastix_data->bcsc->flttype, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
dump_rhs( "AfterPerm", csc->gN, b );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
dump_rhs( "AfterDown", csc->gN, b );
pastix_diag( pastix_data, pastix_data->bcsc->flttype, &sopalin_data, nrhs, b, ldb );
dump_rhs( "AfterDiag", csc->gN, b );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
dump_rhs( "AfterUp", csc->gN, b );
break;
case PastixFactLDLH:
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
pastix_diag( pastix_data->bcsc->flttype, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixConjTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
pastix_diag( pastix_data, pastix_data->bcsc->flttype, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixConjTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
break;
case PastixFactLU:
default:
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data->bcsc->flttype, PastixLeft, PastixUpper, PastixNoTrans, PastixNonUnit, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixLower, PastixNoTrans, PastixUnit, &sopalin_data, nrhs, b, ldb );
pastix_trsm( pastix_data, pastix_data->bcsc->flttype, PastixLeft, PastixUpper, PastixNoTrans, PastixNonUnit, &sopalin_data, nrhs, b, ldb );
break;
}
clockStop(timer);
......@@ -228,6 +252,8 @@ pastix_task_solve( pastix_data_t *pastix_data,
solverBackupRestore( pastix_data->solvmatr, sbackup );
solverBackupExit( sbackup );
dump_rhs( "Final", csc->gN, b );
/* Invalidate following steps, and add factorization step to the ones performed */
pastix_data->steps &= ~( STEP_SOLVE |
STEP_REFINE );
......
......@@ -20,7 +20,7 @@
#include "bcsc.h"
#include "sopalin_data.h"
static void (*sopalinFacto[4][4])(sopalin_data_t*) =
static void (*sopalinFacto[4][4])(pastix_data_t *, sopalin_data_t*) =
{
{ sopalin_spotrf, sopalin_dpotrf, sopalin_cpotrf, sopalin_zpotrf },
{ sopalin_ssytrf, sopalin_dsytrf, sopalin_csytrf, sopalin_zsytrf },
......@@ -233,19 +233,18 @@ pastix_task_sopalin( pastix_data_t *pastix_data,
sopalin_data.diagthreshold = pastix_data->dparm[ DPARM_EPSILON_MAGN_CTRL ] * pastix_data->dparm[DPARM_A_NORM];
}
}
sopalin_data.sched = pastix_data->isched;
sbackup = solverBackupInit( pastix_data->solvmatr );
pastix_data->solvmatr->restore = 2;
{
void (*factofct)( sopalin_data_t *);
void (*factofct)( pastix_data_t *, sopalin_data_t *);
double timer;
factofct = sopalinFacto[ pastix_data->iparm[IPARM_FACTORIZATION] ][csc->flttype-2];
assert(sopalinFacto);
clockStart(timer);
factofct( &sopalin_data );
factofct( pastix_data, &sopalin_data );
clockStop(timer);
pastix_print( 0, 0, OUT_TIME_FACT, clockVal(timer) );
......
......@@ -22,12 +22,13 @@
#include "pastix_zcores.h"
void
sequential_zdiag( sopalin_data_t *sopalin_data,
sequential_zdiag( pastix_data_t *pastix_data, sopalin_data_t *sopalin_data,
int nrhs, pastix_complex64_t *b, int ldb )
{
SolverMatrix *datacode = sopalin_data->solvmtx;
SolverCblk *cblk;
pastix_int_t i, j, k;
(void)pastix_data;
cblk = datacode->cblktab;
for (i=0; i<datacode->cblknbr; i++, cblk++){
......
......@@ -16,7 +16,6 @@
* @precisions normal z -> s d c
*
**/
#define _GNU_SOURCE
#include "common.h"
#include "isched.h"
#include "sopalin_data.h"
......@@ -30,16 +29,19 @@
#endif
void
sequential_zgetrf( sopalin_data_t *sopalin_data )
sequential_zgetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
SolverMatrix *datacode = sopalin_data->solvmtx;
SolverMatrix *datacode = pastix_data->solvmatr;
SolverCblk *cblk;
double threshold = sopalin_data->diagthreshold;
pastix_int_t i;
(void)pastix_data;
cblk = datacode->cblktab;
for (i=0; i<datacode->cblknbr; i++, cblk++){
/* Compute */
core_zgetrfsp1d( datacode, cblk, sopalin_data->diagthreshold );
core_zgetrfsp1d( datacode, cblk, threshold );
}
#if defined(PASTIX_DEBUG_FACTO)
......@@ -78,32 +80,43 @@ thread_pzgetrf( int rank, void *args )
#endif
}
void
thread_zgetrf( sopalin_data_t *sopalin_data )
thread_zgetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
isched_parallel_call( sopalin_data->sched, thread_pzgetrf, sopalin_data );
isched_parallel_call( pastix_data->isched, thread_pzgetrf, sopalin_data );
}
#if defined(PASTIX_WITH_PARSEC)
void
parsec_zgetrf( sopalin_data_t *sopalin_data )
parsec_zgetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
sparse_matrix_desc_t desc;
dague_context_t *ctx;
int argc = 0;
ctx = dague_init( -1, &argc, NULL );
/* Start PaRSEC */
if (pastix_data->parsec == NULL) {
pastix_data->parsec = dague_init( -1, &argc, NULL );
}
ctx = pastix_data->parsec;
/* Create the matrix descriptor */
sparse_matrix_init( &desc, sopalin_data->solvmtx,
pastix_size_of( PastixComplex64 ), 1, 0 );
/* Run the facto */
dsparse_zgetrf_sp( ctx, &desc, sopalin_data );
/* Destroy the decriptor */
sparse_matrix_destroy( &desc );
dague_fini( &ctx );
dague_fini( &(pastix_data->parsec) );
}
#endif
static void (*zgetrf_table[4])(sopalin_data_t *) = {
static void (*zgetrf_table[4])(pastix_data_t *, sopalin_data_t *) = {
sequential_zgetrf,
thread_zgetrf,
#if defined(PASTIX_WITH_PARSEC)
......@@ -115,13 +128,14 @@ static void (*zgetrf_table[4])(sopalin_data_t *) = {
};
void
sopalin_zgetrf( sopalin_data_t *sopalin_data )
sopalin_zgetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
int sched = 0;
void (*zgetrf)(sopalin_data_t *) = zgetrf_table[ sched ];
int sched = pastix_data->iparm[IPARM_SCHEDULER];
void (*zgetrf)(pastix_data_t *, sopalin_data_t *) = zgetrf_table[ sched ];
if (zgetrf == NULL) {
zgetrf = thread_zgetrf;
}
zgetrf( sopalin_data );
zgetrf( pastix_data, sopalin_data );
}
......@@ -21,17 +21,27 @@
#include "sopalin_data.h"
#include "pastix_zcores.h"
#if defined(PASTIX_WITH_PARSEC)
#include <dague.h>
#include <dague/data.h>
#include <dague/data_distribution.h>
#include "parsec/sparse-matrix.h"
#endif
void
sequential_zhetrf( sopalin_data_t *sopalin_data )
sequential_zhetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
SolverMatrix *datacode = sopalin_data->solvmtx;
SolverMatrix *datacode = pastix_data->solvmatr;
SolverCblk *cblk;
double threshold = sopalin_data->diagthreshold;
pastix_int_t i;
(void)pastix_data;
cblk = datacode->cblktab;
for (i=0; i<datacode->cblknbr; i++, cblk++){
/* Compute */
core_zhetrfsp1d( datacode, cblk, sopalin_data->diagthreshold );
core_zhetrfsp1d( datacode, cblk, threshold );
}
#if defined(PASTIX_DEBUG_FACTO)
......@@ -39,7 +49,6 @@ sequential_zhetrf( sopalin_data_t *sopalin_data )
#endif
}
void
thread_pzhetrf( int rank, void *args )
{
......@@ -73,26 +82,61 @@ thread_pzhetrf( int rank, void *args )
void
thread_zhetrf( sopalin_data_t *sopalin_data )
thread_zhetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
isched_parallel_call( sopalin_data->sched, thread_pzhetrf, sopalin_data );
isched_parallel_call( pastix_data->isched, thread_pzhetrf, sopalin_data );
}
static void (*zhetrf_table[4])(sopalin_data_t *) = {
#if defined(PASTIX_WITH_PARSEC)
void
parsec_zhetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
sparse_matrix_desc_t desc;
dague_context_t *ctx;
int argc = 0;
/* Start PaRSEC */
if (pastix_data->parsec == NULL) {
pastix_data->parsec = dague_init( -1, &argc, NULL );
}
ctx = pastix_data->parsec;
/* Create the matrix descriptor */
sparse_matrix_init( &desc, sopalin_data->solvmtx,
pastix_size_of( PastixComplex64 ), 1, 0 );
/* Run the facto */
dsparse_zhetrf_sp( ctx, &desc, sopalin_data );
/* Destroy the decriptor */
sparse_matrix_destroy( &desc );
dague_fini( &(pastix_data->parsec) );
}
#endif
static void (*zhetrf_table[4])(pastix_data_t *, sopalin_data_t *) = {
sequential_zhetrf,
thread_zhetrf,
#if defined(PASTIX_WITH_PARSEC)
parsec_zhetrf,
#else
NULL,
#endif
NULL
};
void
sopalin_zhetrf( sopalin_data_t *sopalin_data )
sopalin_zhetrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
int sched = 0;
void (*zhetrf)(sopalin_data_t *) = zhetrf_table[ sched ];
int sched = pastix_data->iparm[IPARM_SCHEDULER];
void (*zhetrf)(pastix_data_t *, sopalin_data_t *) = zhetrf_table[ sched ];
if (zhetrf == NULL) {
zhetrf = thread_zhetrf;
}
zhetrf( sopalin_data );
zhetrf( pastix_data, sopalin_data );
}
......@@ -21,17 +21,27 @@
#include "sopalin_data.h"
#include "pastix_zcores.h"
#if defined(PASTIX_WITH_PARSEC)
#include <dague.h>
#include <dague/data.h>
#include <dague/data_distribution.h>
#include "parsec/sparse-matrix.h"
#endif
void
sequential_zpotrf( sopalin_data_t *sopalin_data )
sequential_zpotrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
SolverMatrix *datacode = sopalin_data->solvmtx;
SolverMatrix *datacode = pastix_data->solvmatr;
SolverCblk *cblk;
double threshold = sopalin_data->diagthreshold;
pastix_int_t i;
(void)pastix_data;
cblk = datacode->cblktab;
for (i=0; i<datacode->cblknbr; i++, cblk++){
/* Compute */
core_zpotrfsp1d( datacode, cblk, sopalin_data->diagthreshold );
core_zpotrfsp1d( datacode, cblk, threshold );
}
#if defined(PASTIX_DEBUG_FACTO)
......@@ -70,28 +80,62 @@ thread_pzpotrf( int rank, void *args )
#endif
}
void
thread_zpotrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
isched_parallel_call( pastix_data->isched, thread_pzpotrf, sopalin_data );
}
#if defined(PASTIX_WITH_PARSEC)
void
thread_zpotrf( sopalin_data_t *sopalin_data )
parsec_zpotrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
isched_parallel_call( sopalin_data->sched, thread_pzpotrf, sopalin_data );
sparse_matrix_desc_t desc;
dague_context_t *ctx;
int argc = 0;
/* Start PaRSEC */
if (pastix_data->parsec == NULL) {
pastix_data->parsec = dague_init( -1, &argc, NULL );
}
ctx = pastix_data->parsec;
/* Create the matrix descriptor */
sparse_matrix_init( &desc, sopalin_data->solvmtx,
pastix_size_of( PastixComplex64 ), 1, 0 );
/* Run the facto */
dsparse_zpotrf_sp( ctx, &desc, sopalin_data );
/* Destroy the decriptor */
sparse_matrix_destroy( &desc );
dague_fini( &(pastix_data->parsec) );
}
#endif
static void (*zpotrf_table[4])(sopalin_data_t *) = {
static void (*zpotrf_table[4])(pastix_data_t *, sopalin_data_t *) = {
sequential_zpotrf,
thread_zpotrf,
#if defined(PASTIX_WITH_PARSEC)
parsec_zpotrf,
#else
NULL,
#endif
NULL
};
void
sopalin_zpotrf( sopalin_data_t *sopalin_data )
sopalin_zpotrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
int sched = 0;
void (*zpotrf)(sopalin_data_t *) = zpotrf_table[ sched ];
int sched = pastix_data->iparm[IPARM_SCHEDULER];
void (*zpotrf)(pastix_data_t *, sopalin_data_t *) = zpotrf_table[ sched ];
if (zpotrf == NULL) {
zpotrf = thread_zpotrf;
}
zpotrf( sopalin_data );
zpotrf( pastix_data, sopalin_data );
}
......@@ -21,17 +21,27 @@
#include "sopalin_data.h"
#include "pastix_zcores.h"
#if defined(PASTIX_WITH_PARSEC)
#include <dague.h>
#include <dague/data.h>
#include <dague/data_distribution.h>
#include "parsec/sparse-matrix.h"
#endif
void
sequential_zsytrf( sopalin_data_t *sopalin_data )
sequential_zsytrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
SolverMatrix *datacode = sopalin_data->solvmtx;
SolverMatrix *datacode = pastix_data->solvmatr;
SolverCblk *cblk;
double threshold = sopalin_data->diagthreshold;
pastix_int_t i;
(void)pastix_data;
cblk = datacode->cblktab;
for (i=0; i<datacode->cblknbr; i++, cblk++){
/* Compute */
core_zsytrfsp1d( datacode, cblk, sopalin_data->diagthreshold );
core_zsytrfsp1d( datacode, cblk, threshold );
}
#if defined(PASTIX_DEBUG_FACTO)
......@@ -72,26 +82,61 @@ thread_pzsytrf( int rank, void *args )
void
thread_zsytrf( sopalin_data_t *sopalin_data )
thread_zsytrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
isched_parallel_call( pastix_data->isched, thread_pzsytrf, sopalin_data );
}
#if defined(PASTIX_WITH_PARSEC)
void
parsec_zsytrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
isched_parallel_call( sopalin_data->sched, thread_pzsytrf, sopalin_data );
sparse_matrix_desc_t desc;
dague_context_t *ctx;
int argc = 0;
/* Start PaRSEC */
if (pastix_data->parsec == NULL) {
pastix_data->parsec = dague_init( -1, &argc, NULL );
}
ctx = pastix_data->parsec;
/* Create the matrix descriptor */
sparse_matrix_init( &desc, sopalin_data->solvmtx,
pastix_size_of( PastixComplex64 ), 1, 0 );
/* Run the facto */
dsparse_zsytrf_sp( ctx, &desc, sopalin_data );
/* Destroy the decriptor */
sparse_matrix_destroy( &desc );
dague_fini( &(pastix_data->parsec) );
}
#endif
static void (*zsytrf_table[4])(sopalin_data_t *) = {
static void (*zsytrf_table[4])(pastix_data_t *, sopalin_data_t *) = {
sequential_zsytrf,
thread_zsytrf,
#if defined(PASTIX_WITH_PARSEC)
parsec_zsytrf,
#else
NULL,
#endif
NULL
};
void
sopalin_zsytrf( sopalin_data_t *sopalin_data )
sopalin_zsytrf( pastix_data_t *pastix_data,
sopalin_data_t *sopalin_data )
{
int sched = 0;
void (*zsytrf)(sopalin_data_t *) = zsytrf_table[ sched ];
int sched = pastix_data->iparm[IPARM_SCHEDULER];
void (*zsytrf)(pastix_data_t *, sopalin_data_t *) = zsytrf_table[ sched ];
if (zsytrf == NULL) {
zsytrf = thread_zsytrf;
}
zsytrf( sopalin_data );
zsytrf( pastix_data, sopalin_data );
}
......@@ -26,7 +26,7 @@ static pastix_complex64_t zone = 1.0;
static pastix_complex64_t mzone = -1.0;
void
sequential_ztrsm( int side, int uplo, int trans, int diag,
sequential_ztrsm( pastix_data_t *pastix_data, int side, int uplo, int trans, int diag,
sopalin_data_t *sopalin_data,
int nrhs, pastix_complex64_t *b, int ldb )
{
......@@ -35,6 +35,7 @@ sequential_ztrsm( int side, int uplo, int trans, int diag,
SolverBlok *blok;
pastix_complex64_t *coeftab;
pastix_int_t i, j, tempm, tempn;
(void)pastix_data;
/*
* Left / Upper / NoTrans
......