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/aij.h
8: #if defined(PETSC_HAVE_SPOOLES) && !defined(PETSC_USE_SINGLE) && !defined(PETSC_USE_COMPLEX)
9: #include src/mat/impls/aij/seq/spooles.h
11: /* Note the Petsc r and c permutations are ignored */
12: #undef __FUNCT__
14: int MatLUFactorSymbolic_SeqAIJ_Spooles(Mat A,IS r,IS c,MatLUInfo *info,Mat *F)
15: {
16: Mat_Spooles *lu;
17: int ierr,m=A->m,n=A->n;
20: /* Create the factorization matrix F */
21: MatCreateSeqAIJ(A->comm,m,n,PETSC_NULL,PETSC_NULL,F);
23: (*F)->ops->lufactornumeric = MatFactorNumeric_SeqAIJ_Spooles;
24: (*F)->factor = FACTOR_LU;
26: ierr = PetscNew(Mat_Spooles,&lu);
27: (*F)->spptr = (void*)lu;
28: lu->options.symflag = SPOOLES_NONSYMMETRIC;
29: lu->options.pivotingflag = SPOOLES_PIVOTING;
30: lu->flg = DIFFERENT_NONZERO_PATTERN;
31: lu->options.useQR = PETSC_FALSE;
33: if (info->dtcol == 0.0) {
34: lu->options.pivotingflag = SPOOLES_NO_PIVOTING;
35: }
36:
37: return(0);
38: }
40: /* Note the Petsc r and c permutations are ignored */
41: #undef __FUNCT__
43: int MatQRFactorSymbolic_SeqAIJ_Spooles(Mat A,IS r,IS c,MatLUInfo *info,Mat *F)
44: {
45: Mat_Spooles *lu;
46: int ierr,m=A->m,n=A->n;
49: /* Create the factorization matrix F */
50: MatCreateSeqAIJ(A->comm,m,n,PETSC_NULL,PETSC_NULL,F);
52: (*F)->ops->lufactornumeric = MatFactorNumeric_SeqAIJ_Spooles;
53: (*F)->factor = FACTOR_LU;
55: ierr = PetscNew(Mat_Spooles,&lu);
56: (*F)->spptr = (void*)lu;
57: lu->options.symflag = SPOOLES_NONSYMMETRIC;
58: lu->options.pivotingflag = SPOOLES_NO_PIVOTING;
59: lu->flg = DIFFERENT_NONZERO_PATTERN;
60: lu->options.useQR = PETSC_TRUE;
62: return(0);
63: }
65: #undef __FUNCT__
67: int MatUseSpooles_SeqAIJ(Mat A)
68: {
69: int ierr;
70: PetscTruth useQR=PETSC_FALSE;
73: PetscOptionsHasName(A->prefix,"-mat_aij_spooles_qr",&useQR);
74: if (useQR){
75: A->ops->lufactorsymbolic = MatQRFactorSymbolic_SeqAIJ_Spooles;
76: } else {
77: A->ops->lufactorsymbolic = MatLUFactorSymbolic_SeqAIJ_Spooles;
78: }
79: return(0);
80: }
82: #else
84: #undef __FUNCT__
86: int MatUseSpooles_SeqAIJ(Mat A)
87: {
89: return(0);
90: }
92: #endif