Actual source code: aijspooles.c

  1: #define PETSCMAT_DLL

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

 7:  #include ../src/mat/impls/aij/seq/spooles/spooles.h
 8:  #include ../src/mat/impls/aij/seq/aij.h

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

 19:   MatView_SeqAIJ(A,viewer);

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

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

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

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

 62: PetscErrorCode MatFactorGetSolverPackage_seqaij_spooles(Mat A,const MatSolverPackage *type)
 63: {
 65:   *type = MAT_SOLVER_SPOOLES;
 66:   return(0);
 67: }
 69: 
 73: PetscErrorCode MatGetFactor_seqaij_spooles(Mat A,MatFactorType ftype,Mat *F)
 74: {
 75:   Mat            B;
 76:   Mat_Spooles    *lu;
 78:   int            m=A->rmap->n,n=A->cmap->n;

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

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

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

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