Actual source code: mpiaijspooles.c
2: /*
3: Provides an interface to the Spooles parallel sparse solver (MPI SPOOLES)
4: */
7: #include <../src/mat/impls/aij/mpi/mpiaij.h>
8: #include <../src/mat/impls/aij/seq/spooles/spooles.h>
10: /* Note the Petsc r and c permutations are ignored */
13: PetscErrorCode MatLUFactorSymbolic_MPIAIJSpooles(Mat F,Mat A,IS r,IS c,const MatFactorInfo *info)
14: {
16: if (!info->dtcol) {
17: Mat_Spooles *lu = (Mat_Spooles*) F->spptr;
18: lu->options.pivotingflag = SPOOLES_NO_PIVOTING;
19: }
20: F->ops->lufactornumeric = MatFactorNumeric_MPISpooles;
21: return(0);
22: }
27: PetscErrorCode MatFactorGetSolverPackage_spooles(Mat A,const MatSolverPackage *type)
28: {
30: *type = MATSOLVERSPOOLES;
31: return(0);
32: }
38: PetscErrorCode MatGetFactor_mpiaij_spooles(Mat A,MatFactorType ftype,Mat *F)
39: {
40: Mat_Spooles *lu;
41: Mat B;
45: /* Create the factorization matrix F */
46: MatCreate(((PetscObject)A)->comm,&B);
47: MatSetSizes(B,A->rmap->n,A->cmap->n,A->rmap->N,A->cmap->N);
48: MatSetType(B,((PetscObject)A)->type_name);
49: MatMPIAIJSetPreallocation(B,0,PETSC_NULL,0,PETSC_NULL);
51: PetscNewLog(B,Mat_Spooles,&lu);
52: B->spptr = lu;
53: lu->flg = DIFFERENT_NONZERO_PATTERN;
54: lu->options.useQR = PETSC_FALSE;
56: if (ftype == MAT_FACTOR_LU) {
57: B->ops->lufactorsymbolic = MatLUFactorSymbolic_MPIAIJSpooles;
58: B->ops->view = MatView_Spooles;
59: B->ops->destroy = MatDestroy_MPIAIJSpooles;
60: PetscObjectComposeFunctionDynamic((PetscObject)B,"MatFactorGetSolverPackage_C","MatFactorGetSolverPackage_spooles",MatFactorGetSolverPackage_spooles);
61:
63: lu->options.symflag = SPOOLES_NONSYMMETRIC;
64: lu->options.pivotingflag = SPOOLES_PIVOTING;
65: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Only LU for AIJ matrices, use SBAIJ for Cholesky");
67: B->factortype = ftype;
68: MPI_Comm_dup(((PetscObject)A)->comm,&(lu->comm_spooles));
69: *F = B;
70: return(0);
71: }