Actual source code: factorschur.c

petsc-3.8.4 2018-03-24
Report Typos and Errors
  1:  #include <petsc/private/matimpl.h>
  2:  #include <../src/mat/impls/dense/seq/dense.h>

  4: PETSC_INTERN PetscErrorCode MatFactorSetUpInPlaceSchur_Private(Mat F)
  5: {
  6:   Mat              St, S = F->schur;
  7:   MatSolverPackage solvertype;
  8:   MatFactorInfo    info;
  9:   PetscErrorCode   ierr;

 12:   if (S->solvertype) {
 13:     PetscStrallocpy(S->solvertype,&solvertype);
 14:   } else {
 15:     PetscStrallocpy(MATSOLVERPETSC,&solvertype);
 16:   }
 17:   MatSetUnfactored(S);
 18:   MatGetFactor(S,solvertype,F->factortype,&St);
 19:   if (St->factortype == MAT_FACTOR_CHOLESKY) { /* LDL^t regarded as Cholesky */
 20:     MatCholeskyFactorSymbolic(St,S,NULL,&info);
 21:   } else {
 22:     MatLUFactorSymbolic(St,S,NULL,NULL,&info);
 23:   }
 24:   S->ops->solve             = St->ops->solve;
 25:   S->ops->matsolve          = St->ops->matsolve;
 26:   S->ops->solvetranspose    = St->ops->solvetranspose;
 27:   S->ops->matsolvetranspose = St->ops->matsolvetranspose;
 28:   S->ops->solveadd          = St->ops->solveadd;
 29:   S->ops->solvetransposeadd = St->ops->solvetransposeadd;
 30:   S->factortype             = St->factortype;

 32:   MatDestroy(&St);
 33:   PetscFree(solvertype);
 34:   return(0);
 35: }

 37: PETSC_INTERN PetscErrorCode MatFactorUpdateSchurStatus_Private(Mat F)
 38: {
 39:   Mat            S = F->schur;

 43:   switch(F->schur_status) {
 44:   case MAT_FACTOR_SCHUR_UNFACTORED:
 45:   case MAT_FACTOR_SCHUR_INVERTED:
 46:     if (S) {
 47:       S->ops->solve             = NULL;
 48:       S->ops->matsolve          = NULL;
 49:       S->ops->solvetranspose    = NULL;
 50:       S->ops->matsolvetranspose = NULL;
 51:       S->ops->solveadd          = NULL;
 52:       S->ops->solvetransposeadd = NULL;
 53:       S->factortype             = MAT_FACTOR_NONE;
 54:       PetscFree(S->solvertype);
 55:     }
 56:     break;
 57:   case MAT_FACTOR_SCHUR_FACTORED:
 58:     MatFactorSetUpInPlaceSchur_Private(F);
 59:     break;
 60:   default:
 61:     SETERRQ1(PetscObjectComm((PetscObject)F),PETSC_ERR_SUP,"Unhandled MatFactorSchurStatus %D",F->schur_status);
 62:   }
 63:   return(0);
 64: }

 66: /* Schur status updated in the interface */
 67: PETSC_INTERN PetscErrorCode MatFactorFactorizeSchurComplement_Private(Mat F)
 68: {
 69:   MatFactorInfo  info;

 73:   if (F->factortype == MAT_FACTOR_CHOLESKY) { /* LDL^t regarded as Cholesky */
 74:     MatCholeskyFactor(F->schur,NULL,&info);
 75:   } else {
 76:     MatLUFactor(F->schur,NULL,NULL,&info);
 77:   }
 78:   return(0);
 79: }

 81: /* Schur status updated in the interface */
 82: PETSC_INTERN PetscErrorCode MatFactorInvertSchurComplement_Private(Mat F)
 83: {
 84:   Mat S = F->schur;

 87:   if (S) {
 88:     PetscMPIInt    size;
 89:     PetscBool      isdense;

 92:     MPI_Comm_size(PetscObjectComm((PetscObject)S),&size);
 93:     if (size > 1) SETERRQ(PetscObjectComm((PetscObject)S),PETSC_ERR_SUP,"Not yet implemented");
 94:     PetscObjectTypeCompare((PetscObject)S,MATSEQDENSE,&isdense);
 95:     if (!isdense) SETERRQ1(PetscObjectComm((PetscObject)S),PETSC_ERR_SUP,"Not implemented for type %s",((PetscObject)S)->type_name);
 96:     MatSeqDenseInvertFactors_Private(S);
 97:   }
 98:   return(0);
 99: }