Actual source code: factorschur.c

petsc-3.11.4 2019-09-28
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:   MatFactorInfo    info;
  8:   PetscErrorCode   ierr;

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

 26:   MatDestroy(&St);
 27:   return(0);
 28: }

 30: PETSC_INTERN PetscErrorCode MatFactorUpdateSchurStatus_Private(Mat F)
 31: {
 32:   Mat            S = F->schur;

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

 59: /* Schur status updated in the interface */
 60: PETSC_INTERN PetscErrorCode MatFactorFactorizeSchurComplement_Private(Mat F)
 61: {
 62:   MatFactorInfo  info;

 66:   if (F->factortype == MAT_FACTOR_CHOLESKY) { /* LDL^t regarded as Cholesky */
 67:     MatCholeskyFactor(F->schur,NULL,&info);
 68:   } else {
 69:     MatLUFactor(F->schur,NULL,NULL,&info);
 70:   }
 71:   return(0);
 72: }

 74: /* Schur status updated in the interface */
 75: PETSC_INTERN PetscErrorCode MatFactorInvertSchurComplement_Private(Mat F)
 76: {
 77:   Mat S = F->schur;

 80:   if (S) {
 81:     PetscMPIInt    size;
 82:     PetscBool      isdense;

 85:     MPI_Comm_size(PetscObjectComm((PetscObject)S),&size);
 86:     if (size > 1) SETERRQ(PetscObjectComm((PetscObject)S),PETSC_ERR_SUP,"Not yet implemented");
 87:     PetscObjectTypeCompare((PetscObject)S,MATSEQDENSE,&isdense);
 88:     if (!isdense) SETERRQ1(PetscObjectComm((PetscObject)S),PETSC_ERR_SUP,"Not implemented for type %s",((PetscObject)S)->type_name);
 89:     MatSeqDenseInvertFactors_Private(S);
 90:   }
 91:   return(0);
 92: }