35 assert( ipiv != NULL );
37 *ipiv = malloc( N *
sizeof(
int) );
44 info = LAPACKE_spotrf_work( LAPACK_COL_MAJOR,
'L', N, S, lds );
47 info = LAPACKE_sgetrf_work( LAPACK_COL_MAJOR, N, N, S, lds, *ipiv );
50 fprintf(stderr,
"Factorization type not handled by Schur example\n");
56 info = LAPACKE_cpotrf_work( LAPACK_COL_MAJOR,
'L', N, S, lds );
59 info = LAPACKE_cgetrf_work( LAPACK_COL_MAJOR, N, N, S, lds, *ipiv );
62 fprintf(stderr,
"Factorization type not handled by Schur example\n");
68 info = LAPACKE_zpotrf_work( LAPACK_COL_MAJOR,
'L', N, S, lds );
71 info = LAPACKE_zgetrf_work( LAPACK_COL_MAJOR, N, N, S, lds, *ipiv );
74 fprintf(stderr,
"Factorization type not handled by Schur example\n");
80 info = LAPACKE_dpotrf_work( LAPACK_COL_MAJOR,
'L', N, S, lds );
83 info = LAPACKE_dgetrf_work( LAPACK_COL_MAJOR, N, N, S, lds, *ipiv );
86 fprintf(stderr,
"Factorization type not handled by Schur example\n");
90 fprintf(stderr,
"Incorrect arithmetic type\n");
93 fprintf(stderr,
"Error in schurFactorize with info =%d\n", info );
111 assert(ipiv != NULL);
116 float *b = (
float *)bptr;
120 info = LAPACKE_spotrs_work( LAPACK_COL_MAJOR,
'L', Nschur, NRHS, S, lds, b, ldb );
123 info = LAPACKE_sgetrs_work( LAPACK_COL_MAJOR,
'N', Nschur, NRHS, S, lds, *ipiv, b, ldb );
126 fprintf(stderr,
"Factorization type not handled by Schur example\n");
130 case PastixComplex32:
136 info = LAPACKE_cpotrs_work( LAPACK_COL_MAJOR,
'L', Nschur, NRHS, S, lds, b, ldb );
139 info = LAPACKE_cgetrs_work( LAPACK_COL_MAJOR,
'N', Nschur, NRHS, S, lds, *ipiv, b, ldb );
142 fprintf(stderr,
"Factorization type not handled by Schur example\n");
146 case PastixComplex64:
148 pastix_complex64_t *b = (pastix_complex64_t *)bptr;
152 info = LAPACKE_zpotrs_work( LAPACK_COL_MAJOR,
'L', Nschur, NRHS, S, lds, b, ldb );
155 info = LAPACKE_zgetrs_work( LAPACK_COL_MAJOR,
'N', Nschur, NRHS, S, lds, *ipiv, b, ldb );
158 fprintf(stderr,
"Factorization type not handled by Schur example\n");
164 double *b = (
double *)bptr;
168 info = LAPACKE_dpotrs_work( LAPACK_COL_MAJOR,
'L', Nschur, NRHS, S, lds, b, ldb );
171 info = LAPACKE_dgetrs_work( LAPACK_COL_MAJOR,
'N', Nschur, NRHS, S, lds, *ipiv, b, ldb );
174 fprintf(stderr,
"Factorization type not handled by Schur example\n");
179 fprintf(stderr,
"Incorrect arithmetic type\n");
188 fprintf(stderr,
"Error in schurSolve with info =%d\n", info );
194 int main (
int argc,
char **argv)
198 double dparm[DPARM_SIZE];
200 char *filename = NULL;
201 spmatrix_t *spm, spm2;
202 void *x, *b, *S, *x0 = NULL;
223 &check, &scatter, &driver, &filename );
229 fprintf(stderr,
"This types of factorization (LDL^t and LDL^h) are not supported by this example.\n");
236 pastixInit( &pastix_data, MPI_COMM_WORLD, iparm, dparm );
241 spm = malloc(
sizeof( spmatrix_t ) );
243 rc = spmReadDriverDist( driver, filename, spm, MPI_COMM_WORLD );
246 rc = spmReadDriver( driver, filename, spm );
249 if ( rc != SPM_SUCCESS ) {
254 spmPrintInfo( spm, stdout );
256 rc = spmCheckAndCorrect( spm, &spm2 );
265 if ( spm->flttype == SpmPattern ) {
266 spmGenFakeValues( spm );
273 nschur = spm->gN / 3;
275 nschur = (nschur > 5000) ? 5000 : nschur;
282 for (i=0; i<nschur; i++) {
299 double normA = spmNorm( SpmFrobeniusNorm, spm );
300 spmScal( 1./normA, spm );
311 S = malloc( pastix_size_of( spm->flttype ) * nschur * lds );
327 nschur, S, lds, &ipiv );
333 size = pastix_size_of( spm->flttype ) * spm->nexp * nrhs;
343 spmGenRHS( SpmRhsRndX, nrhs, spm, x0, spm->nexp, b, spm->nexp );
344 memcpy( x, b, size );
347 spmGenRHS( SpmRhsRndB, nrhs, spm, NULL, spm->nexp, x, spm->nexp );
350 spmScalMat( 1./normA, spm, nrhs, b, spm->nexp );
359 spm->nexp, nrhs, x, ldb, Xp );
375 size = pastix_size_of( spm->flttype ) * nschur * nrhs;
376 schur_x = malloc( size );
381 nschur, nrhs, S, lds, schur_x, nschur, &ipiv );
401 spm->nexp, nrhs, x, ldb, Xp );
BEGIN_C_DECLS typedef int pastix_int_t
float _Complex pastix_complex32_t
int pastixRhsInit(pastix_rhs_t *rhs)
Initialize an RHS data structure.
int pastixRhsFinalize(pastix_rhs_t rhs)
Cleanup an RHS data structure.
spm_coeftype_t pastix_coeftype_t
Arithmetic types.
void pastixFinalize(pastix_data_t **pastix_data)
Finalize the solver instance.
enum pastix_diag_e pastix_diag_t
Diagonal.
void pastixInitParam(pastix_int_t *iparm, double *dparm)
Initialize the iparm and dparm arrays to their default values.
void pastixInit(pastix_data_t **pastix_data, PASTIX_Comm pastix_comm, pastix_int_t *iparm, double *dparm)
Initialize the solver instance.
enum pastix_factotype_e pastix_factotype_t
Factorization algorithms available for IPARM_FACTORIZATION parameter.
@ DPARM_EPSILON_REFINEMENT
void pastixGetOptions(int argc, char **argv, pastix_int_t *iparm, double *dparm, int *check, int *scatter, spm_driver_t *driver, char **filename)
PaStiX helper function to read command line options in examples.
void pastixSetSchurUnknownList(pastix_data_t *pastix_data, pastix_int_t n, const pastix_int_t *list)
Set the list of unknowns that belongs to the schur complement.
int pastixGetSchur(const pastix_data_t *pastix_data, void *S, pastix_int_t lds)
Return the Schur complement.
int pastixRhsSchurSet(const pastix_data_t *pastix_data, pastix_int_t m, pastix_int_t n, void *B, pastix_int_t ldb, pastix_rhs_t rhsB)
Set the vector in an RHS data structure.
int pastix_subtask_applyorder(pastix_data_t *pastix_data, pastix_dir_t dir, pastix_int_t m, pastix_int_t n, void *B, pastix_int_t ldb, pastix_rhs_t Bp)
Apply a permutation on the right-and-side vector before the solve step.
int pastix_subtask_trsm(pastix_data_t *pastix_data, pastix_side_t side, pastix_uplo_t uplo, pastix_trans_t trans, pastix_diag_t diag, pastix_rhs_t b)
Apply a triangular solve on the right-and-side vectors.
int pastixRhsSchurGet(const pastix_data_t *pastix_data, pastix_int_t m, pastix_int_t n, pastix_rhs_t rhsB, void *B, pastix_int_t ldb)
Get the vector in an RHS data structure.
int pastix_task_analyze(pastix_data_t *pastix_data, const spmatrix_t *spm)
Perform all the preprocessing steps: ordering, symbolic factorization, reordering,...
int pastix_task_numfact(pastix_data_t *pastix_data, spmatrix_t *spm)
Perform all the numerical factorization steps: fill the internal block CSC and the solver matrix stru...
Main PaStiX data structure.
Main PaStiX RHS structure.