Actual source code: factorschur.c

petsc-3.13.6 2020-09-29
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;

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

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

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

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

 64:   PetscLogEventBegin(MAT_FactorFactS,F,0,0,0);
 65:   if (F->factortype == MAT_FACTOR_CHOLESKY) { /* LDL^t regarded as Cholesky */
 66:     MatCholeskyFactor(F->schur,NULL,&info);
 67:   } else {
 68:     MatLUFactor(F->schur,NULL,NULL,&info);
 69:   }
 70:   PetscLogEventEnd(MAT_FactorFactS,F,0,0,0);
 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,isdensecuda;

 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:     PetscObjectTypeCompare((PetscObject)S,MATSEQDENSECUDA,&isdensecuda);
 89:     PetscLogEventBegin(MAT_FactorInvS,F,0,0,0);
 90:     if (isdense) {
 91:       MatSeqDenseInvertFactors_Private(S);
 92: #if defined(PETSC_HAVE_CUDA)
 93:     } else if (isdensecuda) {
 94:       MatSeqDenseCUDAInvertFactors_Private(S);
 95: #endif
 96:     } else SETERRQ1(PetscObjectComm((PetscObject)S),PETSC_ERR_SUP,"Not implemented for type %s",((PetscObject)S)->type_name);
 97:     PetscLogEventEnd(MAT_FactorInvS,F,0,0,0);
 98:   }
 99:   return(0);
100: }