Actual source code: aijspooles.c
1: /*
2: Provides an interface to the Spooles serial sparse solver
3: */
5: #include src/mat/impls/aij/seq/spooles/spooles.h
9: PetscErrorCode MatView_SeqAIJSpooles(Mat A,PetscViewer viewer)
10: {
12: PetscTruth iascii;
13: PetscViewerFormat format;
14: Mat_Spooles *lu=(Mat_Spooles*)(A->spptr);
17: (*lu->MatView)(A,viewer);
19: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
20: if (iascii) {
21: PetscViewerGetFormat(viewer,&format);
22: if (format == PETSC_VIEWER_ASCII_FACTOR_INFO) {
23: MatFactorInfo_Spooles(A,viewer);
24: }
25: }
26: return(0);
27: }
31: PetscErrorCode MatAssemblyEnd_SeqAIJSpooles(Mat A,MatAssemblyType mode) {
33: Mat_Spooles *lu=(Mat_Spooles *)(A->spptr);
36: (*lu->MatAssemblyEnd)(A,mode);
38: lu->MatLUFactorSymbolic = A->ops->lufactorsymbolic;
39: lu->MatCholeskyFactorSymbolic = A->ops->choleskyfactorsymbolic;
40: if (lu->useQR){
41: A->ops->lufactorsymbolic = MatQRFactorSymbolic_SeqAIJSpooles;
42: } else {
43: A->ops->choleskyfactorsymbolic = MatCholeskyFactorSymbolic_SeqAIJSpooles;
44: A->ops->lufactorsymbolic = MatLUFactorSymbolic_SeqAIJSpooles;
45: }
46: return(0);
47: }
49: /* Note the Petsc r and c permutations are ignored */
52: PetscErrorCode MatLUFactorSymbolic_SeqAIJSpooles(Mat A,IS r,IS c,MatFactorInfo *info,Mat *F)
53: {
54: Mat B;
55: Mat_Spooles *lu;
57: int 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) {
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: PetscErrorCode MatQRFactorSymbolic_SeqAIJSpooles(Mat A,IS r,IS c,MatFactorInfo *info,Mat *F)
85: {
86: Mat B;
87: Mat_Spooles *lu;
89: int m=A->m,n=A->n;
92: SETERRQ(PETSC_ERR_SUP,"QR Factorization is unsupported as the Spooles implementation of QR is invalid.");
93: /* Create the factorization matrix */
94: MatCreate(A->comm,m,n,PETSC_NULL,PETSC_NULL,&B);
95: MatSetType(B,A->type_name);
96: MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);
98: B->ops->lufactornumeric = MatFactorNumeric_SeqAIJSpooles;
99: B->factor = FACTOR_LU;
101: lu = (Mat_Spooles*)(B->spptr);
102: lu->options.symflag = SPOOLES_NONSYMMETRIC;
103: lu->options.pivotingflag = SPOOLES_NO_PIVOTING;
104: lu->flg = DIFFERENT_NONZERO_PATTERN;
105: lu->options.useQR = PETSC_TRUE;
107: *F = B;
108: return(0);
109: }
111: /* Note the Petsc r permutation is ignored */
114: PetscErrorCode MatCholeskyFactorSymbolic_SeqAIJSpooles(Mat A,IS r,MatFactorInfo *info,Mat *F)
115: {
116: Mat B;
117: Mat_Spooles *lu;
119: int m=A->m,n=A->n;
122: /* Create the factorization matrix */
123: MatCreate(A->comm,m,n,PETSC_NULL,PETSC_NULL,&B);
124: MatSetType(B,A->type_name);
125: MatSeqAIJSetPreallocation(B,PETSC_NULL,PETSC_NULL);
127: B->ops->choleskyfactornumeric = MatFactorNumeric_SeqAIJSpooles;
128: B->ops->getinertia = MatGetInertia_SeqSBAIJSpooles;
129: B->factor = FACTOR_CHOLESKY;
131: lu = (Mat_Spooles*)(B->spptr);
132: lu->options.pivotingflag = SPOOLES_NO_PIVOTING;
133: lu->options.symflag = SPOOLES_SYMMETRIC; /* default */
134: lu->flg = DIFFERENT_NONZERO_PATTERN;
135: lu->options.useQR = PETSC_FALSE;
137: *F = B;
138: return(0);
139: }