Actual source code: aijspooles.c

  1: /*$Id: aijspooles.c,v 1.10 2001/08/15 15:56:50 bsmith Exp $*/
  2: /* 
  3:    Provides an interface to the Spooles serial sparse solver
  4: */

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

 10: int MatView_SeqAIJSpooles(Mat A,PetscViewer viewer)
 11: {
 12:   int               ierr;
 13:   PetscTruth        isascii;
 14:   PetscViewerFormat format;
 15:   Mat_Spooles       *lu=(Mat_Spooles*)(A->spptr);

 18:   (*lu->MatView)(A,viewer);

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

 32: int MatAssemblyEnd_SeqAIJSpooles(Mat A,MatAssemblyType mode) {
 33:   int         ierr;
 34:   Mat_Spooles *lu=(Mat_Spooles *)(A->spptr);

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

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

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

 60:   /* Create the factorization matrix */
 61:   MatCreate(A->comm,m,n,PETSC_NULL,PETSC_NULL,&B);
 62:   MatSetType(B,A->type_name);
 63:   MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);

 65:   B->ops->lufactornumeric  = MatFactorNumeric_SeqAIJSpooles;
 66:   B->factor                = FACTOR_LU;

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

 74:   if (info->dtcol == 0.0) {
 75:     lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
 76:   }
 77:   *F = B;
 78:   return(0);
 79: }

 81: /* Note the Petsc r and c permutations are ignored */
 84: int MatQRFactorSymbolic_SeqAIJSpooles(Mat A,IS r,IS c,MatFactorInfo *info,Mat *F)
 85: {
 86:   Mat          B;
 87:   Mat_Spooles  *lu;
 88:   int          ierr,m=A->m,n=A->n;

 91:   SETERRQ(PETSC_ERR_SUP,"QR Factorization is unsupported as the Spooles implementation of QR is invalid.");
 92:   /* Create the factorization matrix */
 93:   MatCreate(A->comm,m,n,PETSC_NULL,PETSC_NULL,&B);
 94:   MatSetType(B,A->type_name);
 95:   MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);

 97:   B->ops->lufactornumeric  = MatFactorNumeric_SeqAIJSpooles;
 98:   B->factor                = FACTOR_LU;

100:   lu                        = (Mat_Spooles*)(B->spptr);
101:   lu->options.symflag       = SPOOLES_NONSYMMETRIC;
102:   lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
103:   lu->flg                   = DIFFERENT_NONZERO_PATTERN;
104:   lu->options.useQR         = PETSC_TRUE;

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

110: /* Note the Petsc r permutation is ignored */
113: int MatCholeskyFactorSymbolic_SeqAIJSpooles(Mat A,IS r,MatFactorInfo *info,Mat *F)
114: {
115:   Mat         B;
116:   Mat_Spooles *lu;
117:   int         ierr,m=A->m,n=A->n;

120:   /* Create the factorization matrix */
121:   MatCreate(A->comm,m,n,PETSC_NULL,PETSC_NULL,&B);
122:   MatSetType(B,A->type_name);
123:   MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);

125:   B->ops->choleskyfactornumeric  = MatFactorNumeric_SeqAIJSpooles;
126:   B->ops->getinertia             = MatGetInertia_SeqSBAIJSpooles;
127:   B->factor                      = FACTOR_CHOLESKY;

129:   lu                        = (Mat_Spooles*)(B->spptr);
130:   lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
131:   lu->options.symflag       = SPOOLES_SYMMETRIC;   /* default */
132:   lu->flg                   = DIFFERENT_NONZERO_PATTERN;
133:   lu->options.useQR         = PETSC_FALSE;

135:   *F = B;
136:   return(0);
137: }