Actual source code: sbaijspooles.c
1: /*$Id: sbaijspooles.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/sbaij/seq/sbaij.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: /*
12: input:
13: F: numeric factor
14: output:
15: nneg, nzero, npos: matrix inertia
16: */
18: #undef __FUNCT__
20: int MatGetInertia_SeqSBAIJ_Spooles(Mat F,int *nneg,int *nzero,int *npos)
21: {
22: Mat_Spooles *lu = (Mat_Spooles*)F->spptr;
23: int ierr,neg,zero,pos;
26: FrontMtx_inertia(lu->frontmtx, &neg, &zero, &pos) ;
27: if(nneg) *nneg = neg;
28: if(nzero) *nzero = zero;
29: if(npos) *npos = pos;
30: return(0);
31: }
33: /* Note the Petsc r permutation is ignored */
34: #undef __FUNCT__
36: int MatCholeskyFactorSymbolic_SeqSBAIJ_Spooles(Mat A,IS r,PetscReal f,Mat *F)
37: {
38: Mat_Spooles *lu;
39: int ierr,m=A->m,n=A->n;
42: /* Create the factorization matrix F */
43: MatCreateSeqAIJ(A->comm,m,n,PETSC_NULL,PETSC_NULL,F);
45: (*F)->ops->choleskyfactornumeric = MatFactorNumeric_SeqAIJ_Spooles;
46: (*F)->ops->getinertia = MatGetInertia_SeqSBAIJ_Spooles;
47: (*F)->factor = FACTOR_CHOLESKY;
49: ierr = PetscNew(Mat_Spooles,&lu);
50: (*F)->spptr = (void*)lu;
51: lu->options.symflag = SPOOLES_SYMMETRIC;
52: lu->options.pivotingflag = SPOOLES_NO_PIVOTING;
53: lu->flg = DIFFERENT_NONZERO_PATTERN;
54: lu->options.useQR = PETSC_FALSE;
56: return(0);
57: }
59: #undef __FUNCT__
61: int MatUseSpooles_SeqSBAIJ(Mat A)
62: {
63: Mat_SeqSBAIJ *sbaij = (Mat_SeqSBAIJ*)A->data;
64: int bs = sbaij->bs;
67: if (bs > 1) SETERRQ1(1,"Block size %d not supported by Spooles",bs);
68: A->ops->choleskyfactorsymbolic = MatCholeskyFactorSymbolic_SeqSBAIJ_Spooles;
69: return(0);
70: }
72: #else
74: #undef __FUNCT__
76: int MatUseSpooles_SeqSBAIJ(Mat A)
77: {
79: return(0);
80: }
82: #endif