Actual source code: mpiaijspooles.c

  1: /*$Id: mpiaijspooles.c,v 1.10 2001/08/15 15:56:50 bsmith Exp $*/
  2: /* 
  3:    Provides an interface to the Spooles parallel sparse solver (MPI SPOOLES)
  4: */


 7:  #include src/mat/impls/aij/mpi/mpiaij.h
 8:  #include src/mat/impls/aij/seq/spooles/spooles.h

 12: int MatAssemblyEnd_MPIAIJSpooles(Mat A,MatAssemblyType mode) {
 13:   int         ierr;
 14:   Mat_Spooles *lu=(Mat_Spooles *)(A->spptr);

 17:   (*lu->MatAssemblyEnd)(A,mode);
 18:   lu->MatLUFactorSymbolic  = A->ops->lufactorsymbolic;
 19:   A->ops->lufactorsymbolic = MatLUFactorSymbolic_MPIAIJSpooles;
 20:   return(0);
 21: }

 23: /* Note the Petsc r and c permutations are ignored */
 26: int MatLUFactorSymbolic_MPIAIJSpooles(Mat A,IS r,IS c,MatFactorInfo *info,Mat *F)
 27: {
 28:   Mat_Spooles *lu;
 29:   Mat         B;
 30:   int         ierr;


 34:   /* Create the factorization matrix F */
 35:   MatCreate(A->comm,A->m,A->n,A->M,A->N,&B);
 36:   MatSetType(B,A->type_name);
 37:   MatMPIAIJSetPreallocation(B,0,PETSC_NULL,0,PETSC_NULL);

 39:   B->ops->lufactornumeric = MatFactorNumeric_MPIAIJSpooles;
 40:   B->factor               = FACTOR_LU;

 42:   lu                       = (Mat_Spooles *)(B->spptr);
 43:   lu->options.symflag      = SPOOLES_NONSYMMETRIC;
 44:   lu->options.pivotingflag = SPOOLES_PIVOTING;
 45:   lu->flg                  = DIFFERENT_NONZERO_PATTERN;
 46:   lu->options.useQR        = PETSC_FALSE;

 48:   MPI_Comm_dup(A->comm,&(lu->comm_spooles));

 50:   if (info->dtcol == 0.0) {
 51:     lu->options.pivotingflag  = SPOOLES_NO_PIVOTING;
 52:   }
 53:   *F = B;
 54:   return(0);
 55: }