31 typedef struct solve_param {
33 double dparm[DPARM_SIZE];
50 static void *solve_smp(
void *arg)
55 void *x, *b, *x0 = NULL;
61 solve_param_t param = *(solve_param_t *)arg;
67 scatter = param.scatter;
81 param.iparm, param.dparm,
86 #if defined(PASTIX_WITH_MPI)
89 MPI_Comm_size( MPI_COMM_WORLD, &size );
90 MPI_Comm_rank( MPI_COMM_WORLD, &rank );
93 fprintf( stderr,
"\nWarning: Reentrant doesn't work with multiple MPI instances\n" );
94 fprintf( stderr,
"Quitting now\n" );
105 spm = malloc(
sizeof( spmatrix_t ) );
107 rc = spmReadDriverDist( param.driver, param.filename, spm, MPI_COMM_WORLD );
110 rc = spmReadDriver( param.driver, param.filename, spm );
112 if ( rc != SPM_SUCCESS ) {
117 spmPrintInfo( spm, stdout );
119 rc = spmCheckAndCorrect( spm, &spm2 );
129 if ( spm->flttype == SpmPattern ) {
130 spmGenFakeValues( spm );
141 double normA = spmNorm( SpmFrobeniusNorm, spm );
142 spmScal( 1./normA, spm );
153 size = pastix_size_of( spm->flttype ) * spm->nexp * nrhs;
162 spmGenRHS( SpmRhsRndX, nrhs, spm, x0, spm->nexp, b, spm->nexp );
163 memcpy( x, b, size );
166 spmGenRHS( SpmRhsRndB, nrhs, spm, NULL, spm->nexp, x, spm->nexp );
169 spmScalMat( 1./normA, spm, nrhs, b, spm->nexp );
172 memcpy( b, x, size );
184 spm, x0, spm->nexp, b, spm->nexp, x, spm->nexp );
200 int main (
int argc,
char **argv)
203 double dparm[DPARM_SIZE];
205 char *filename = NULL;
206 int nbcallingthreads = 2;
207 solve_param_t *solve_param;
224 &check, &scatter, &driver, &filename );
229 solve_param = (solve_param_t*) malloc(nbcallingthreads *
sizeof(solve_param_t));
230 threads = (pthread_t*) malloc(nbcallingthreads *
sizeof(pthread_t));
232 for (i = 0; i < nbcallingthreads; i++)
234 memcpy(solve_param[i].iparm, iparm,
sizeof(solve_param[i].iparm));
235 memcpy(solve_param[i].dparm, dparm,
sizeof(solve_param[i].dparm));
236 solve_param[i].check = check;
237 solve_param[i].scatter = scatter;
238 solve_param[i].id = i;
239 solve_param[i].driver = driver;
240 solve_param[i].filename = filename;
241 solve_param[i].rc = 0;
246 pthread_create(&threads[i], NULL, solve_smp, (
void *)&solve_param[i]);
252 for (i = 0; i < nbcallingthreads; i++) {
253 pthread_join(threads[i],(
void**)NULL);
254 rc += solve_param[i].rc;
BEGIN_C_DECLS typedef int pastix_int_t
void pastixFinalize(pastix_data_t **pastix_data)
Finalize the solver instance.
void pastixInitParam(pastix_int_t *iparm, double *dparm)
Initialize the iparm and dparm arrays to their default values.
void pastixInitWithAffinity(pastix_data_t **pastix_data, PASTIX_Comm pastix_comm, pastix_int_t *iparm, double *dparm, const int *bindtab)
Initialize the solver instance with a bintab array to specify the thread binding.
@ 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.
int pastix_task_refine(pastix_data_t *pastix_data, pastix_int_t n, pastix_int_t nrhs, void *B, pastix_int_t ldb, void *X, pastix_int_t ldx)
Perform iterative refinement.
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_solve(pastix_data_t *pastix_data, pastix_int_t m, pastix_int_t nrhs, void *B, pastix_int_t ldb)
Solve the given problem.
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.