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

 11: PetscErrorCode MatView_SeqAIJSpooles(Mat A,PetscViewer viewer)
 12: {
 14:   PetscTruth        iascii;
 15:   PetscViewerFormat format;
 16:   Mat_Spooles       *lu=(Mat_Spooles*)(A->spptr);

 19:   (*lu->MatView)(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: }

 33: PetscErrorCode MatAssemblyEnd_SeqAIJSpooles(Mat A,MatAssemblyType mode) {
 35:   Mat_Spooles *lu=(Mat_Spooles *)(A->spptr);

 38:   (*lu->MatAssemblyEnd)(A,mode);

 40:   lu->MatLUFactorSymbolic          = A->ops->lufactorsymbolic;
 41:   lu->MatCholeskyFactorSymbolic    = A->ops->choleskyfactorsymbolic;
 42:   if (lu->useQR){
 43:     A->ops->lufactorsymbolic       = MatQRFactorSymbolic_SeqAIJSpooles;
 44:   } else {
 45:     A->ops->choleskyfactorsymbolic = MatCholeskyFactorSymbolic_SeqAIJSpooles;
 46:     A->ops->lufactorsymbolic       = MatLUFactorSymbolic_SeqAIJSpooles;
 47:   }
 48:   return(0);
 49: }

 51: /* Note the Petsc r and c permutations are ignored */
 54: PetscErrorCode MatLUFactorSymbolic_SeqAIJSpooles(Mat A,IS r,IS c,MatFactorInfo *info,Mat *F)
 55: {
 56:   Mat          B;
 57:   Mat_Spooles  *lu;
 59:   int          m=A->rmap.n,n=A->cmap.n;

 62:   /* Create the factorization matrix */
 63:   MatCreate(A->comm,&B);
 64:   MatSetSizes(B,m,n,PETSC_DETERMINE,PETSC_DETERMINE);
 65:   MatSetType(B,A->type_name);
 66:   MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);

 68:   B->ops->lufactornumeric  = MatFactorNumeric_SeqAIJSpooles;
 69:   B->factor                = FACTOR_LU;

 71:   lu                        = (Mat_Spooles*)(B->spptr);
 72:   lu->options.symflag       = SPOOLES_NONSYMMETRIC;
 73:   lu->options.pivotingflag  = SPOOLES_PIVOTING;
 74:   lu->flg                   = DIFFERENT_NONZERO_PATTERN;
 75:   lu->options.useQR         = PETSC_FALSE;

 77:   if (!info->dtcol) {
 78:     lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
 79:   }
 80:   *F = B;
 81:   return(0);
 82: }

 84: /* Note the Petsc r and c permutations are ignored */
 87: PetscErrorCode MatQRFactorSymbolic_SeqAIJSpooles(Mat A,IS r,IS c,MatFactorInfo *info,Mat *F)
 88: {
 89:   Mat          B;
 90:   Mat_Spooles  *lu;
 92:   int          m=A->rmap.n,n=A->cmap.n;

 95:   SETERRQ(PETSC_ERR_SUP,"QR Factorization is unsupported as the Spooles implementation of QR is invalid.");
 96:   /* Create the factorization matrix */
 97:   MatCreate(A->comm,&B);
 98:   MatSetSizes(B,m,n,PETSC_DETERMINE,PETSC_DETERMINE);
 99:   MatSetType(B,A->type_name);
100:   MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);

102:   B->ops->lufactornumeric  = MatFactorNumeric_SeqAIJSpooles;
103:   B->factor                = FACTOR_LU;

105:   lu                        = (Mat_Spooles*)(B->spptr);
106:   lu->options.symflag       = SPOOLES_NONSYMMETRIC;
107:   lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
108:   lu->flg                   = DIFFERENT_NONZERO_PATTERN;
109:   lu->options.useQR         = PETSC_TRUE;

111:   *F = B;
112:   return(0);
113: }

115: /* Note the Petsc r permutation is ignored */
118: PetscErrorCode MatCholeskyFactorSymbolic_SeqAIJSpooles(Mat A,IS r,MatFactorInfo *info,Mat *F)
119: {
120:   Mat         B;
121:   Mat_Spooles *lu;
123:   int         m=A->rmap.n,n=A->cmap.n;

126:   /* Create the factorization matrix */
127:   MatCreate(A->comm,&B);
128:   MatSetSizes(B,m,n,PETSC_DETERMINE,PETSC_DETERMINE);
129:   MatSetType(B,A->type_name);
130:   MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);

132:   B->ops->choleskyfactornumeric  = MatFactorNumeric_SeqAIJSpooles;
133:   B->ops->getinertia             = MatGetInertia_SeqSBAIJSpooles;
134:   B->factor                      = FACTOR_CHOLESKY;

136:   lu                        = (Mat_Spooles*)(B->spptr);
137:   lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
138:   lu->options.symflag       = SPOOLES_SYMMETRIC;   /* default */
139:   lu->flg                   = DIFFERENT_NONZERO_PATTERN;
140:   lu->options.useQR         = PETSC_FALSE;

142:   *F = B;
143:   return(0);
144: }