pastix.h

Summary
pastix.h
PaStiX interfaceParallel Sparse matriX solver interface
Main PaStiX functions
pastixComputes steps of the resolution of Ax=b linear system, using direct methods.
Examplefrom file simple.c :
dpastixComputes steps of the resolution of Ax=b linear system, using direct methods.
Examplefrom file simple_dist.c :
Thread functions
pastix_bindThreadsSet bindtab in pastix_data, it gives for each thread the CPU to bind in to.
Checking the matrix.
pastix_checkMatrixCheck the matrix :
Getting solver distribution.
pastix_getLocalNodeNbrReturn the node number in the new distribution computed by analyze step.
pastix_getLocalNodeLstFill in nodelst with the list of local nodes/columns.
pastix_getLocalUnknownNbrReturn the unknown number in the new distribution computed by analyze step.
pastix_getLocalUnknownLstFill in unknownlst with the list of local unknowns/clumns.
About the Schur complement.
pastix_setSchurUnknownListSet the list of unknowns to isolate at the end of the matrix via permutations.
pastix_getSchurLocalNodeNbrCompute the number of nodes in the local part of the Schur.
pastix_getSchurLocalUnkownNbrCompute the number of unknowns in the local part of the Schur.
pastix_getSchurLocalNodeListCompute the list of nodes in the local part of the Schur.
pastix_getSchurLocalUnkownListCompute the list of unknowns in the local part of the Schur.
pastix_setSchurArrayGive user memory area to store Schur in PaStiX.
pastix_getSchurGet the Schur complement from PaStiX.
About parameters.
pastix_initParamSets default parameters for iparm and dparm
About Scaling.Working on it...

PaStiX interface

ComplexFloat_

Parallel Sparse matriX solver interface

Authors

Astrid CASADEIastrid..nosp@m.casadei@inri.nosp@m.a.fr
Mathieu FAVERGEfav.nosp@m.erge@labr.nosp@m.i.fr
Xavier LACOSTExavier..nosp@m.lacoste@inri.nosp@m.a.fr
Pierre RAMETra.nosp@m.met@labr.nosp@m.i.fr

Blabla

Main PaStiX functions

pastix

void pastix(pastix_data_t **pastix_data,
MPI_Comm pastix_comm,
INT n,
INT *colptr,
INT *row,
FLOAT *avals,
INT *perm,
INT *invp,
FLOAT *b,
INT rhs,
INT *iparm,
double *dparm)

Computes steps of the resolution of Ax=b linear system, using direct methods.

The matrix is given in CSC format.

Parameters

pastix_dataData used for a step by step execution.
pastix_commMPI communicator which compute the resolution.
nSize of the system.
colptrTabular containing the start of each column in row and avals tabulars.
rowTabular containing the row number for each element sorted by column.
avalsTabular containing the values of each element sorted by column.
permPermutation tabular for the renumerotation of the unknowns.
invpReverse permutation tabular for the renumerotation of the unknowns.
bRight hand side vector(s).
rhsNumber of right hand side vector(s).
iparmInteger parameters given to pastix.
dparmDouble parameters given to pâstix.
Summary
Examplefrom file simple.c :

Example

from file simple.c :

/\*******************************************\/
/\*    Check Matrix format                  *\/
/\*******************************************\/
/\*
 * Matrix needs :
 *    - to be in fortran numbering
 *    - to have only the lower triangular part in symmetric case
 *    - to have a graph with a symmetric structure in unsymmetric case
 *\/
pastix_checkMatrix( MPI_COMM_WORLD, verbosemode,
                    (MTX_ISSYM(type) ? API_SYM_YES : API_SYM_NO),
                    API_YES,
                    ncol, &colptr, &rows, &values, NULL);

/\*******************************************\/
/\* Initialize parameters to default values *\/
/\*******************************************\/
iparm[IPARM_MODIFY_PARAMETER] = API_NO;
pastix(&pastix_data, MPI_COMM_WORLD,
       ncol, colptr, rows, values,
       perm, invp, rhs, 1, iparm, dparm);

/\*******************************************\/
/\*       Customize some parameters         *\/
/\*******************************************\/
iparm[IPARM_THREAD_NBR] = nbthread;
if (MTX_ISSYM(type)) {
  :wiparm[IPARM_SYM]           = API_SYM_YES;
  iparm[IPARM_FACTORIZATION] = API_FACT_LDLT;
}
else{
  iparm[IPARM_SYM]           = API_SYM_NO;
  iparm[IPARM_FACTORIZATION] = API_FACT_LU;
}
iparm[IPARM_START_TASK]          = API_TASK_ORDERING;
iparm[IPARM_END_TASK]            = API_TASK_CLEAN;

/\*******************************************\/
/\*           Save the rhs                  *\/
/\*    (it will be replaced by solution)    *\/
/\*******************************************\/
rhssaved = malloc(ncol*sizeof(pastix_float_t));
memcpy(rhssaved, rhs, ncol*sizeof(pastix_float_t));

/\*******************************************\/
/\*           Call pastix                   *\/
/\*******************************************\/
perm = malloc(ncol*sizeof(pastix_int_t));
invp = malloc(ncol*sizeof(pastix_int_t));

pastix(&pastix_data, MPI_COMM_WORLD,
 ncol, colptr, rows, values,
 perm, invp, rhs, 1, iparm, dparm);

dpastix

Computes steps of the resolution of Ax=b linear system, using direct methods.  Here the matrix is given distributed.

The matrix is given in CSCD format.

Parameters

pastix_dataData used for a step by step execution.
pastix_commMPI communicator which compute the resolution.
nSize of the system.
colptrTabular containing the start of each column in row and avals tabulars.
rowTabular containing the row number for each element sorted by column.
avalsTabular containing the values of each element sorted by column.
loc2globGlobal column number of the local columns.
permPermutation tabular for the renumerotation of the unknowns.
invpReverse permutation tabular for the renumerotation of the unknowns.
bRight hand side vector(s).
rhsNumber of right hand side vector(s).
iparmInteger parameters given to pastix.
dparmDouble parameters given to pâstix.

Example

from file simple_dist.c :

/\*******************************************\/
/\*    Check Matrix format                  *\/
/\*******************************************\/
/\*
 * Matrix needs :
 *    - to be in fortran numbering
 *    - to have only the lower triangular part in symmetric case
 *    - to have a graph with a symmetric structure in unsymmetric case
 *\/
pastix_checkMatrix( MPI_COMM_WORLD, verbosemode,
                    (MTX_ISSYM(type) ? API_SYM_YES : API_SYM_NO),
                    API_YES,
                    ncol, &colptr, &rows, &values, &loc2glob, 1);

/\*******************************************\/
/\* Initialize parameters to default values *\/
/\*******************************************\/
iparm[IPARM_MODIFY_PARAMETER] = API_NO;
dpastix(&pastix_data, MPI_COMM_WORLD,
        ncol, colptr, rows, values, loc2glob,
        perm, invp, rhs, 1, iparm, dparm);

/\*******************************************\/
/\*       Customize some parameters         *\/
/\*******************************************\/
iparm[IPARM_THREAD_NBR] = nbthread;
if (MTX_ISSYM(type))
  {
    iparm[IPARM_SYM]           = API_SYM_YES;
    iparm[IPARM_FACTORIZATION] = API_FACT_LDLT;
  }
else
  {
    iparm[IPARM_SYM]           = API_SYM_NO;
    iparm[IPARM_FACTORIZATION] = API_FACT_LU;
  }
iparm[IPARM_MATRIX_VERIFICATION] = API_NO;

iparm[IPARM_START_TASK]          = API_TASK_ORDERING;
iparm[IPARM_END_TASK]            = API_TASK_BLEND;

/\*******************************************\/
/\*           Call pastix                   *\/
/\*******************************************\/
perm = malloc(ncol*sizeof(pastix_int_t));
/\* No need to allocate invp in dpastix *\/

dpastix(&pastix_data, MPI_COMM_WORLD,
        ncol, colptr, rows, NULL, loc2glob,
        perm, NULL, NULL, 1, iparm, dparm);

/\* Redistributing the matrix *\/

ncol2 = pastix_getLocalNodeNbr(&pastix_data);

if (NULL == (loc2glob2 = malloc(ncol2 * sizeof(pastix_int_t))))
  {
    fprintf(stderr, "Malloc error\n");
    return EXIT_FAILURE;
  }

pastix_getLocalNodeLst(&pastix_data, loc2glob2);

if (EXIT_SUCCESS != cscd_redispatch(ncol,   colptr,   rows,
                                    values,    rhs,  loc2glob,
                                    ncol2, &colptr2, &rows2,
                                    &values2, &rhs2, loc2glob2,
                                    MPI_COMM_WORLD))
  return EXIT_FAILURE;

free(colptr);
free(rows);
free(values);
free(rhs);
free(loc2glob);
free(perm);

iparm[IPARM_START_TASK]          = API_TASK_NUMFACT;
iparm[IPARM_END_TASK]            = API_TASK_CLEAN;


dpastix(&pastix_data, MPI_COMM_WORLD,
        ncol2, colptr2, rows2, values2, loc2glob2,
        perm, invp, rhs2, 1, iparm, dparm);

Thread functions

pastix_bindThreads

Set bindtab in pastix_data, it gives for each thread the CPU to bind in to. bindtab follows this organisation :

bindtab[threadnum] = cpu to set thread threadnum.

Parameters

pastix_dataStructure de donnée pour l’utilisation step by step
thrdnbrNombre de threads / Taille du tableau
bindtabTableau de correspondance entre chaque thread et coeur de la machine

Checking the matrix.

pastix_checkMatrix

Check the matrix :

  • Renumbers in Fortran numerotation (base 1) if needed (base 0)
  • Can scale the matrix if compiled with -DMC64 -DSCALING (untested)
  • Checks the symetry of the graph in non symmetric mode.  With non distributed matrices, with flagcor == API_YES, tries to correct the matrix.
  • sort the CSC.

Parameters

pastix_commPaStiX MPI communicator
verbLevel of prints (API_VERBOSE_[NOT|NO|YES])
flagsymIndicate if the given matrix is symetric (API_SYM_YES or API_SYM_NO)
flagcorIndicate if we permit the function to reallocate the matrix.
nSize of the matrix.
colptrFirst element of each row in row and avals.
rowRow of each element of the matrix.
avalsValue of each element of the matrix.
loc2globGlobal column number of local columns.
dofNumber of degrees of freedom.

Getting solver distribution.

pastix_getLocalNodeNbr

Return the node number in the new distribution computed by analyze step.  Needs analyze step to be runned with pastix_data before.

Parameters

pastix_dataData used for a step by step execution.

Returns

Number of local nodes/columns in new distribution.

pastix_getLocalNodeLst

Fill in nodelst with the list of local nodes/columns.  Needs nodelst to be allocated with nodenbr*sizeof(pastix_int_t), where nodenbr has been computed by pastix_getLocalNodeNbr.

Parameters

pastix_dataData used for a step by step execution.
nodelstAn array where to write the list of local nodes/columns.

pastix_getLocalUnknownNbr

Return the unknown number in the new distribution computed by analyze step.  Needs analyze step to be runned with pastix_data before.

Parameters

pastix_dataData used for a step by step execution.

Returns

Number of local unknowns/columns in new distribution.

pastix_getLocalUnknownLst

Fill in unknownlst with the list of local unknowns/clumns.  Needs unknownlst to be allocated with unknownnbr*sizeof(pastix_int_t), where unknownnbr has been computed by pastix_getLocalUnknownNbr.

Parameters

pastix_dataData used for a step by step execution.
unknownlstAn array where to write the list of local unknowns/columns.

About the Schur complement.

pastix_setSchurUnknownList

Set the list of unknowns to isolate at the end of the matrix via permutations.

Has to be called if using IPARM_SCHUR = API_YES.

Parameters

pastix_dataData used for a step by step execution.
nNumber of unknowns.
listList of unknowns.

pastix_getSchurLocalNodeNbr

Compute the number of nodes in the local part of the Schur.

Parameters

pastix_dataCommon data structure for PaStiX calls.
nodeNbr(out) Number of nodes in Schur (local).

Returns

NO_ERRFor the moment

TODO: Error management.

pastix_getSchurLocalUnkownNbr

Compute the number of unknowns in the local part of the Schur.

Parameters

pastix_dataCommon data structure for PaStiX calls.
unknownNbr(out) Number of unknowns in Schur (local).

Returns

NO_ERRFor the moment

TODO: Error management.

pastix_getSchurLocalNodeList

Compute the list of nodes in the local part of the Schur.

Parameters

pastix_dataCommon data structure for PaStiX calls.
nodes(out) Nodes in Schur (local).

Returns

NO_ERRFor the moment

TODO: Error management.

pastix_getSchurLocalUnkownList

Compute the list of unknowns in the local part of the Schur.

Parameters

pastix_dataCommon data structure for PaStiX calls.
unknowns(out) Unknowns in Schur (local).

Returns

NO_ERRFor the moment

TODO: Error management.

pastix_setSchurArray

Give user memory area to store Schur in PaStiX.

Parameters

pastix_dataCommon data structure for PaStiX calls.
arrayMemory area to store the Schur.

Returns

NO_ERRFor the moment

TODO: Error management.

pastix_getSchur

Get the Schur complement from PaStiX.

Schur complement is a dense block in a column scheme.

Parameters

pastix_dataData used for a step by step execution.
schurArray to fill-in with Schur complement.

About parameters.

pastix_initParam

Sets default parameters for iparm and dparm

Parameters

iparmtabular of IPARM_SIZE integer parameters.
dparmtabular of DPARM_SIZE double parameters.

About Scaling.

Working on it...

ComplexFloat_
Parallel Sparse matriX solver interface
void pastix(pastix_data_t **pastix_data,
MPI_Comm pastix_comm,
INT n,
INT *colptr,
INT *row,
FLOAT *avals,
INT *perm,
INT *invp,
FLOAT *b,
INT rhs,
INT *iparm,
double *dparm)
Computes steps of the resolution of Ax=b linear system, using direct methods.
A simple example : read the matrix, check it is correct and correct it if needed, then run pastix in one call.
A simple example with a distributed matrix : read the matrix, check it is correct and correct it if needed, distribute it then run pastix in one call.
Return the node number in the new distribution computed by analyze step.
Return the unknown number in the new distribution computed by analyze step.
Close