Actual source code: aijspooles.c

  2: /* 
  3:    Provides an interface to the Spooles serial sparse solver
  4: */

  6: #include <../src/mat/impls/aij/seq/spooles/spooles.h>
  7: #include <../src/mat/impls/aij/seq/aij.h>

 11: PetscErrorCode MatView_Spooles(Mat A,PetscViewer viewer)
 12: {
 13:   PetscErrorCode    ierr;
 14:   PetscBool         iascii;
 15:   PetscViewerFormat format;

 18:   MatView_SeqAIJ(A,viewer);

 20:   PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
 21:   if (iascii) {
 22:     PetscViewerGetFormat(viewer,&format);
 23:     if (format == PETSC_VIEWER_ASCII_INFO) {
 24:       MatFactorInfo_Spooles(A,viewer);
 25:     }
 26:   }
 27:   return(0);
 28: }

 30: /* Note the Petsc r and c permutations are ignored */
 33: PetscErrorCode MatLUFactorSymbolic_SeqAIJSpooles(Mat F,Mat A,IS r,IS c,const MatFactorInfo *info)
 34: {
 35:   Mat_Spooles    *lu = (Mat_Spooles*)(F->spptr);;

 38:   F->ops->lufactornumeric =  MatFactorNumeric_SeqSpooles;
 39:   if (!info->dtcol) {
 40:     lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
 41:   }
 42:   return(0);
 43: }

 45: /* Note the Petsc r permutation is ignored */
 48: PetscErrorCode MatCholeskyFactorSymbolic_SeqAIJSpooles(Mat F,Mat A,IS r,const MatFactorInfo *info)
 49: {
 51:   F->ops->choleskyfactornumeric  = MatFactorNumeric_SeqSpooles;
 52: #if !defined(PETSC_USE_COMPLEX)
 53:   F->ops->getinertia             = MatGetInertia_SeqSBAIJSpooles;
 54: #endif
 55:   return(0);
 56: }

 61: PetscErrorCode MatFactorGetSolverPackage_seqaij_spooles(Mat A,const MatSolverPackage *type)
 62: {
 64:   *type = MATSOLVERSPOOLES;
 65:   return(0);
 66: }
 68: 
 72: PetscErrorCode MatGetFactor_seqaij_spooles(Mat A,MatFactorType ftype,Mat *F)
 73: {
 74:   Mat            B;
 75:   Mat_Spooles    *lu;
 77:   int            m=A->rmap->n,n=A->cmap->n;

 80:   /* Create the factorization matrix */
 81:   MatCreate(((PetscObject)A)->comm,&B);
 82:   MatSetSizes(B,m,n,PETSC_DETERMINE,PETSC_DETERMINE);
 83:   MatSetType(B,((PetscObject)A)->type_name);
 84:   MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);
 85: 
 86:   PetscNewLog(B,Mat_Spooles,&lu);
 87:   B->spptr = lu;
 88:   lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
 89:   lu->flg                   = DIFFERENT_NONZERO_PATTERN;
 90:   lu->options.useQR         = PETSC_FALSE;

 92:   if (ftype == MAT_FACTOR_LU) {
 93:     B->ops->lufactorsymbolic = MatLUFactorSymbolic_SeqAIJSpooles;

 95:     lu->options.symflag      = SPOOLES_NONSYMMETRIC;
 96:     lu->options.pivotingflag = SPOOLES_PIVOTING;
 97:   } else if (ftype == MAT_FACTOR_CHOLESKY) {
 98:     B->ops->choleskyfactorsymbolic = MatCholeskyFactorSymbolic_SeqAIJSpooles;
 99:     lu->options.symflag            = SPOOLES_SYMMETRIC;   /* default */
100:   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Spooles only supports LU and Cholesky factorizations");
101:   B->ops->view    = MatView_Spooles;
102:   B->ops->destroy = MatDestroy_SeqAIJSpooles;
103:   PetscObjectComposeFunctionDynamic((PetscObject)B,"MatFactorGetSolverPackage_C","MatFactorGetSolverPackage_seqaij_spooles",MatFactorGetSolverPackage_seqaij_spooles);
104:   B->factortype   = ftype;

106:   *F = B;
107:   return(0);
108: }
109: