petscdmlibmeshimpl.C File Reference

Go to the source code of this file.

Classes
struct	DM_libMesh
struct	DMVec_libMesh
struct	token

Functions
PetscErrorCode	DMlibMeshGetVec_Private (DM, const char *, Vec *)
PetscErrorCode	DMlibMeshSetUpName_Private (DM dm)
PetscErrorCode	DMlibMeshSetSystem_libMesh (DM dm, NonlinearImplicitSystem &sys)
PetscErrorCode	DMlibMeshGetSystem_libMesh (DM dm, NonlinearImplicitSystem *&sys)
PetscErrorCode	DMlibMeshGetBlocks (DM dm, PetscInt *n, char ***blocknames)
PetscErrorCode	DMlibMeshGetVariables (DM dm, PetscInt *n, char ***varnames)
static PetscErrorCode	DMCreateFieldDecomposition_libMesh (DM dm, PetscInt *len, char ***namelist, IS **islist, DM **dmlist)
static PetscErrorCode	DMCreateDomainDecomposition_libMesh (DM dm, PetscInt *len, char ***namelist, IS **innerislist, IS **outerislist, DM **dmlist)
PetscErrorCode	DMlibMeshCreateFieldDecompositionDM (DM dm, PetscInt dnumber, PetscInt *dsizes, char ***dvarlists, DM *ddm)
PetscErrorCode	DMlibMeshCreateDomainDecompositionDM (DM dm, PetscInt dnumber, PetscInt *dsizes, char ***dblocklists, DM *ddm)
static PetscErrorCode	DMlibMeshFunction (DM dm, Vec x, Vec r)
static PetscErrorCode	SNESFunction_DMlibMesh (SNES, Vec x, Vec r, void *ctx)
static PetscErrorCode	DMlibMeshJacobian (DM dm, Vec x, Mat jac, Mat pc)
static PetscErrorCode	SNESJacobian_DMlibMesh (SNES, Vec x, Mat jac, Mat pc, void *ctx)
static PetscErrorCode	DMVariableBounds_libMesh (DM dm, Vec xl, Vec xu)
static PetscErrorCode	DMCreateGlobalVector_libMesh (DM dm, Vec *x)
static PetscErrorCode	DMCreateMatrix_libMesh (DM dm, Mat *A)
static PetscErrorCode	DMView_libMesh (DM dm, PetscViewer viewer)
static PetscErrorCode	DMSetUp_libMesh (DM dm)
static PetscErrorCode	DMDestroy_libMesh (DM dm)
PetscErrorCode	DMCreate_libMesh (DM dm)

Function Documentation

DMCreate_libMesh()

PetscErrorCode DMCreate_libMesh

(

DM

dm

)

Definition at line 944 of file petscdmlibmeshimpl.C.

{
  DM_libMesh     * dlm;
 
  PetscFunctionBegin;
  PetscValidHeaderSpecific(dm,DM_CLASSID,1);
  LibmeshPetscCallQ(PetscNew(&dlm));
  dm->data = dlm;
 
  /* DMlibMesh impl */
  dlm->varids     = new(std::map<std::string, unsigned int>);
  dlm->blockids   = new(std::map<std::string, unsigned int>);
  dlm->varnames   = new(std::map<unsigned int, std::string>);
  dlm->blocknames = new(std::map<unsigned int, std::string>);
  dlm->decomposition   = LIBMESH_PETSC_NULLPTR;
  dlm->decomposition_type  = DMLIBMESH_NO_DECOMPOSITION;
 
  /* DM API */
  dm->ops->createglobalvector = DMCreateGlobalVector_libMesh;
  dm->ops->createlocalvector  = 0; // DMCreateLocalVector_libMesh;
  dm->ops->getcoloring        = 0; // DMGetColoring_libMesh;
  dm->ops->creatematrix       = DMCreateMatrix_libMesh;
  dm->ops->createinterpolation= 0; // DMCreateInterpolation_libMesh;
 
  dm->ops->refine             = 0; // DMRefine_libMesh;
  dm->ops->coarsen            = 0; // DMCoarsen_libMesh;
 
  // * dm->ops->getinjection was renamed to dm->ops->createinjection in PETSc 5a84ad338 (5 Jul 2019)
  // * dm->ops-getaggregates was removed in PETSc 97779f9a (5 Jul 2019)
  // * Both changes were merged into PETSc master in 94aad3ce (7 Jul 2019).
#if PETSC_VERSION_LESS_THAN(3,12,0)
  dm->ops->getinjection       = 0; // DMGetInjection_libMesh;
  dm->ops->getaggregates      = 0; // DMGetAggregates_libMesh;
#else
  dm->ops->createinjection = 0;
#endif
 
 
  dm->ops->createfielddecomposition    = DMCreateFieldDecomposition_libMesh;
  dm->ops->createdomaindecomposition   = DMCreateDomainDecomposition_libMesh;
 
  dm->ops->destroy            = DMDestroy_libMesh;
  dm->ops->view               = DMView_libMesh;
  dm->ops->setfromoptions     = 0; // DMSetFromOptions_libMesh;
  dm->ops->setup              = DMSetUp_libMesh;
 
  /* DMlibMesh API */
  LibmeshPetscCallQ(PetscObjectComposeFunction((PetscObject)dm,"DMlibMeshSetSystem_C",DMlibMeshSetSystem_libMesh));
  LibmeshPetscCallQ(PetscObjectComposeFunction((PetscObject)dm,"DMlibMeshGetSystem_C",DMlibMeshGetSystem_libMesh));
 
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::blockids, DM_libMesh::blocknames, DM_libMesh::decomposition, DM_libMesh::decomposition_type, DMCreateDomainDecomposition_libMesh(), DMCreateFieldDecomposition_libMesh(), DMCreateGlobalVector_libMesh(), DMCreateMatrix_libMesh(), DMDestroy_libMesh(), DMlibMeshGetSystem_libMesh(), DMlibMeshSetSystem_libMesh(), DMSetUp_libMesh(), DMView_libMesh(), libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), DM_libMesh::varids, and DM_libMesh::varnames.

DMCreateDomainDecomposition_libMesh()

static PetscErrorCode DMCreateDomainDecomposition_libMesh ( DM  dm, PetscInt *  len, char ***  namelist, IS **  innerislist, IS **  outerislist, DM **  dmlist  )

static

Definition at line 378 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  DM_libMesh     * dlm = (DM_libMesh *)(dm->data);
  NonlinearImplicitSystem * sys = dlm->sys;
  IS emb;
  if (dlm->decomposition_type != DMLIBMESH_DOMAIN_DECOMPOSITION) PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
  *len = cast_int<unsigned int>(dlm->decomposition->size());
  if (namelist)      {LibmeshPetscCallQ(PetscMalloc(*len*sizeof(char *), namelist));}
  if (innerislist)   {LibmeshPetscCallQ(PetscMalloc(*len*sizeof(IS),    innerislist));}
  if (outerislist)   *outerislist = LIBMESH_PETSC_NULLPTR; /* FIX: allow mesh-based overlap. */
  if (dmlist)        {LibmeshPetscCallQ(PetscMalloc(*len*sizeof(DM),    dmlist));}
  for (auto d : index_range(*dlm->decomposition)) {
    std::set<numeric_index_type>               dindices;
    std::string                          dname;
    std::map<std::string, unsigned int>  dblockids;
    std::map<unsigned int,std::string>   dblocknames;
    unsigned int                         dbcount = 0;
    for (const auto & b : (*dlm->decomposition)[d]) {
      std::string bname = (*dlm->blocknames)[b];
      dblockids.insert(std::pair<std::string, unsigned int>(bname,b));
      dblocknames.insert(std::pair<unsigned int,std::string>(b,bname));
      if (!dbcount) dname = bname;
      else   dname += "_" + bname;
      ++dbcount;
      if (!innerislist) continue;
      const subdomain_id_type b_sbd_id = cast_int<subdomain_id_type>(b);
      MeshBase::const_element_iterator       el     = sys->get_mesh().active_local_subdomain_elements_begin(b_sbd_id);
      const MeshBase::const_element_iterator end_el = sys->get_mesh().active_local_subdomain_elements_end(b_sbd_id);
      for ( ; el != end_el; ++el) {
        const Elem * elem = *el;
        std::vector<numeric_index_type> evindices;
        // Iterate only over this DM's variables.
        for (const auto & pr : *(dlm->varids)) {
          // Get the degree of freedom indices for the given variable off the current element.
          sys->get_dof_map().dof_indices(elem, evindices, pr.second);
          for (const auto & dof : evindices) {
            if (dof >= sys->get_dof_map().first_dof() && dof < sys->get_dof_map().end_dof()) // might want to use variable_first/last_local_dof instead
              dindices.insert(dof);
          }
        }
      }
    }
    if (namelist) {
      LibmeshPetscCallQ(PetscStrallocpy(dname.c_str(),(*namelist)+d));
    }
    if (innerislist) {
      PetscInt * darray;
      IS dis;
      LibmeshPetscCallQ(PetscMalloc(sizeof(PetscInt)*dindices.size(), &darray));
      numeric_index_type i = 0;
      for (const auto & id : dindices) {
        darray[i] = id;
        ++i;
      }
      LibmeshPetscCallQ(ISCreateGeneral(((PetscObject)dm)->comm,
                                        cast_int<PetscInt>(dindices.size()),
                                        darray, PETSC_OWN_POINTER, &dis));
      if (dlm->embedding) {
        /* Create a relative embedding into the parent's index space. */
        LibmeshPetscCallQ(ISEmbed(dis,dlm->embedding, PETSC_TRUE, &emb));
        PetscInt elen, dlen;
        LibmeshPetscCallQ(ISGetLocalSize(emb, &elen));
        LibmeshPetscCallQ(ISGetLocalSize(dis, &dlen));
        if (elen != dlen) LIBMESH_SETERRQ1(((PetscObject)dm)->comm, PETSC_ERR_PLIB, "Failed to embed field %zu" , d);
        LibmeshPetscCallQ(ISDestroy(&dis));
        dis = emb;
      }
      else {
        emb = dis;
      }
      if (innerislist) {
        LibmeshPetscCallQ(PetscObjectReference((PetscObject)dis));
        (*innerislist)[d] = dis;
      }
      LibmeshPetscCallQ(ISDestroy(&dis));
    }
    if (dmlist) {
      DM ddm;
      LibmeshPetscCallQ(DMCreate(((PetscObject)dm)->comm, &ddm));
      LibmeshPetscCallQ(DMSetType(ddm, DMLIBMESH));
      DM_libMesh * ddlm = (DM_libMesh *)(ddm->data);
      ddlm->sys = dlm->sys;
      /* copy over the varids and varnames */
      *ddlm->varids    = *dlm->varids;
      *ddlm->varnames  = *dlm->varnames;
      /* set the blocks from the d-th part of the decomposition. */
      *ddlm->blockids    = dblockids;
      *ddlm->blocknames  = dblocknames;
      LibmeshPetscCallQ(PetscObjectReference((PetscObject)emb));
      ddlm->embedding = emb;
      ddlm->embedding_type = DMLIBMESH_DOMAIN_EMBEDDING;
 
      LibmeshPetscCallQ(DMlibMeshSetUpName_Private(ddm));
      LibmeshPetscCallQ(DMSetFromOptions(ddm));
      (*dmlist)[d] = ddm;
    }
  }
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References b, DM_libMesh::blockids, DM_libMesh::blocknames, DM_libMesh::decomposition, DM_libMesh::decomposition_type, DMlibMeshSetUpName_Private(), libMesh::DofMap::dof_indices(), DM_libMesh::embedding, DM_libMesh::embedding_type, libMesh::DofMapBase::end_dof(), libMesh::DofMapBase::first_dof(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::index_range(), libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), DM_libMesh::sys, DM_libMesh::varids, and DM_libMesh::varnames.

Referenced by DMCreate_libMesh().

DMCreateFieldDecomposition_libMesh()

static PetscErrorCode DMCreateFieldDecomposition_libMesh ( DM  dm, PetscInt *  len, char ***  namelist, IS **  islist, DM **  dmlist  )

static

Definition at line 278 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  DM_libMesh     * dlm = (DM_libMesh *)(dm->data);
  NonlinearImplicitSystem * sys = dlm->sys;
  IS emb;
  if (dlm->decomposition_type != DMLIBMESH_FIELD_DECOMPOSITION) PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
 
  *len = cast_int<unsigned int>(dlm->decomposition->size());
  if (namelist) {LibmeshPetscCallQ(PetscMalloc(*len*sizeof(char *), namelist));}
  if (islist)   {LibmeshPetscCallQ(PetscMalloc(*len*sizeof(IS),    islist));}
  if (dmlist)   {LibmeshPetscCallQ(PetscMalloc(*len*sizeof(DM),    dmlist));}
  DofMap & dofmap = dlm->sys->get_dof_map();
  for (auto d : index_range(*dlm->decomposition)) {
    std::set<numeric_index_type>         dindices;
    std::string                          dname;
    std::map<std::string, unsigned int>  dvarids;
    std::map<unsigned int, std::string>  dvarnames;
    unsigned int                         dvcount = 0;
    for (const auto & v : (*dlm->decomposition)[d]) {
      std::string vname = (*dlm->varnames)[v];
      dvarids.insert(std::pair<std::string, unsigned int>(vname,v));
      dvarnames.insert(std::pair<unsigned int,std::string>(v,vname));
      if (!dvcount) dname = vname;
      else   dname += "_" + vname;
      ++dvcount;
      if (!islist) continue;
      // Iterate only over this DM's blocks.
      for (const auto & pr : *(dlm->blockids)) {
        const subdomain_id_type sbd_id = cast_int<subdomain_id_type>(pr.second);
        for (const auto & elem :
               as_range(sys->get_mesh().active_local_subdomain_elements_begin(sbd_id),
                        sys->get_mesh().active_local_subdomain_elements_end(sbd_id))) {
          //unsigned int e_subdomain = elem->subdomain_id();
          std::vector<numeric_index_type> evindices;
          // Get the degree of freedom indices for the given variable off the current element.
          dofmap.dof_indices(elem, evindices, v);
          for (numeric_index_type dof : evindices) {
            if (dof >= dofmap.first_dof() && dof < dofmap.end_dof()) // might want to use variable_first/last_local_dof instead
              dindices.insert(dof);
          }
        }
      }
    }
    if (namelist) {
      LibmeshPetscCallQ(PetscStrallocpy(dname.c_str(),(*namelist)+d));
    }
    if (islist) {
      IS dis;
      PetscInt * darray;
      LibmeshPetscCallQ(PetscMalloc(sizeof(PetscInt)*dindices.size(), &darray));
      numeric_index_type i = 0;
      for (const auto & id : dindices) {
        darray[i] = id;
        ++i;
      }
      LibmeshPetscCallQ(ISCreateGeneral(((PetscObject)dm)->comm,
                                        cast_int<PetscInt>(dindices.size()),
                                        darray, PETSC_OWN_POINTER, &dis));
      if (dlm->embedding) {
        /* Create a relative embedding into the parent's index space. */
        LibmeshPetscCallQ(ISEmbed(dis,dlm->embedding, PETSC_TRUE, &emb));
        PetscInt elen, dlen;
        LibmeshPetscCallQ(ISGetLocalSize(emb, &elen));
        LibmeshPetscCallQ(ISGetLocalSize(dis, &dlen));
        if (elen != dlen) LIBMESH_SETERRQ1(((PetscObject)dm)->comm, PETSC_ERR_PLIB, "Failed to embed subdomain %zu", d);
        LibmeshPetscCallQ(ISDestroy(&dis));
        dis = emb;
      }
      else {
        emb = dis;
      }
      (*islist)[d] = dis;
    }
    if (dmlist) {
      DM ddm;
      LibmeshPetscCallQ(DMCreate(((PetscObject)dm)->comm, &ddm));
      LibmeshPetscCallQ(DMSetType(ddm, DMLIBMESH));
      DM_libMesh * ddlm = (DM_libMesh *)(ddm->data);
      ddlm->sys = dlm->sys;
      /* copy over the block ids and names */
      *ddlm->blockids = *dlm->blockids;
      *ddlm->blocknames = *dlm->blocknames;
      /* set the vars from the d-th part of the decomposition. */
      *ddlm->varids     = dvarids;
      *ddlm->varnames   = dvarnames;
      LibmeshPetscCallQ(PetscObjectReference((PetscObject)emb));
      ddlm->embedding = emb;
      ddlm->embedding_type = DMLIBMESH_FIELD_EMBEDDING;
 
      LibmeshPetscCallQ(DMlibMeshSetUpName_Private(ddm));
      LibmeshPetscCallQ(DMSetFromOptions(ddm));
      (*dmlist)[d] = ddm;
    }
  }
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::as_range(), DM_libMesh::blockids, DM_libMesh::blocknames, DM_libMesh::decomposition, DM_libMesh::decomposition_type, DMlibMeshSetUpName_Private(), libMesh::DofMap::dof_indices(), DM_libMesh::embedding, DM_libMesh::embedding_type, libMesh::DofMapBase::end_dof(), libMesh::DofMapBase::first_dof(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::index_range(), libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), DM_libMesh::sys, DM_libMesh::varids, and DM_libMesh::varnames.

Referenced by DMCreate_libMesh().

DMCreateGlobalVector_libMesh()

static PetscErrorCode DMCreateGlobalVector_libMesh ( DM  dm, Vec *  x  )

static

Definition at line 762 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  DM_libMesh     * dlm = (DM_libMesh *)(dm->data);
  PetscBool eq;
 
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH, &eq));
 
  if (!eq)
    LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
 
  if (!dlm->sys)
    SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "No libMesh system set for DM_libMesh");
 
  NumericVector<Number> * nv = (dlm->sys->solution).get();
  PetscVector<Number> *   pv = dynamic_cast<PetscVector<Number> *>(nv);
  Vec                    v  = pv->vec();
  /* Unfortunately, currently this does not produce a ghosted vector, so nonlinear subproblem solves aren't going to be easily available.
     Should work fine for getting vectors out for linear subproblem solvers. */
  if (dlm->embedding) {
    PetscInt n;
    LibmeshPetscCallQ(VecCreate(((PetscObject)v)->comm, x));
    LibmeshPetscCallQ(ISGetLocalSize(dlm->embedding, &n));
    LibmeshPetscCallQ(VecSetSizes(*x,n,PETSC_DETERMINE));
    LibmeshPetscCallQ(VecSetType(*x,((PetscObject)v)->type_name));
    LibmeshPetscCallQ(VecSetFromOptions(*x));
    LibmeshPetscCallQ(VecSetUp(*x));
  }
  else {
    LibmeshPetscCallQ(VecDuplicate(v,x));
  }
 
#if PETSC_VERSION_LESS_THAN(3,13,0)
  LibmeshPetscCallQ(PetscObjectCompose((PetscObject)*x,"DM",(PetscObject)dm));
#else
  LibmeshPetscCallQ(VecSetDM(*x, dm));
#endif
 
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::embedding, libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), libMesh::System::solution, DM_libMesh::sys, and libMesh::PetscVector< T >::vec().

Referenced by DMCreate_libMesh().

DMCreateMatrix_libMesh()

static PetscErrorCode DMCreateMatrix_libMesh ( DM  dm, Mat *  A  )

static

Definition at line 808 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  DM_libMesh     * dlm = (DM_libMesh *)(dm->data);
  PetscBool eq;
 
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH, &eq));
 
  if (!eq)
    LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
 
  if (!dlm->sys)
    SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "No libMesh system set for DM_libMesh");
 
  *A = (dynamic_cast<PetscMatrixBase<Number> *>(dlm->sys->matrix))->mat();
  LibmeshPetscCallQ(PetscObjectReference((PetscObject)(*A)));
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::LibmeshPetscCallQ(), libMesh::ImplicitSystem::matrix, libMesh::PetscFunctionReturn(), and DM_libMesh::sys.

Referenced by DMCreate_libMesh().

DMDestroy_libMesh()

static PetscErrorCode DMDestroy_libMesh ( DM  dm )

static

Definition at line 926 of file petscdmlibmeshimpl.C.

{
  DM_libMesh * dlm = (DM_libMesh *)(dm->data);
  PetscFunctionBegin;
  delete dlm->varids;
  delete dlm->varnames;
  delete dlm->blockids;
  delete dlm->blocknames;
  delete dlm->decomposition;
  LibmeshPetscCallQ(ISDestroy(&dlm->embedding));
  LibmeshPetscCallQ(PetscFree(dm->data));
 
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::blockids, DM_libMesh::blocknames, DM_libMesh::decomposition, DM_libMesh::embedding, libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), DM_libMesh::varids, and DM_libMesh::varnames.

Referenced by DMCreate_libMesh().

DMlibMeshCreateDomainDecompositionDM()

PetscErrorCode DMlibMeshCreateDomainDecompositionDM	(	DM	dm,
		PetscInt	dnumber,
		PetscInt *	dsizes,
		char ***	dblocklists,
		DM *	ddm
	)

Definition at line 537 of file petscdmlibmeshimpl.C.

{
  PetscBool islibmesh;
  PetscFunctionBegin;
  PetscValidHeaderSpecific(dm,DM_CLASSID,1);
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH,&islibmesh));
  if (!islibmesh) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Got DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
  if (dnumber < 0) LIBMESH_SETERRQ1(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Negative number %" LIBMESH_PETSCINT_FMT " of decomposition parts", dnumber);
#if PETSC_RELEASE_GREATER_EQUALS(3,20,0)
  PetscAssertPointer(ddm,5);
#else
  PetscValidPointer(ddm,5);
#endif
  DM_libMesh * dlm = (DM_libMesh *)(dm->data);
  LibmeshPetscCallQ(DMCreate(((PetscObject)dm)->comm, ddm));
  LibmeshPetscCallQ(DMSetType(*ddm, DMLIBMESH));
  DM_libMesh * ddlm = (DM_libMesh *)((*ddm)->data);
  ddlm->sys = dlm->sys;
  ddlm->varids   = dlm->varids;
  ddlm->varnames = dlm->varnames;
  ddlm->blockids   = dlm->blockids;
  ddlm->blocknames = dlm->blocknames;
  ddlm->decomposition = new(std::vector<std::set<unsigned int>>);
  ddlm->decomposition_type = DMLIBMESH_DOMAIN_DECOMPOSITION;
  if (dnumber) {
    for (PetscInt d = 0; d < dnumber; ++d) {
      if (dsizes[d] < 0) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Negative size %" LIBMESH_PETSCINT_FMT " of decomposition part %" LIBMESH_PETSCINT_FMT, dsizes[d],d);
      ddlm->decomposition->push_back(std::set<unsigned int>());
      for (PetscInt b = 0; b < dsizes[d]; ++b) {
        std::string bname(dblocklists[d][b]);
        std::map<std::string, unsigned int>::const_iterator bit = dlm->blockids->find(bname);
        if (bit == dlm->blockids->end())
          LIBMESH_SETERRQ3(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Block %" LIBMESH_PETSCINT_FMT " on the %" LIBMESH_PETSCINT_FMT "-th list with name %s is not owned by this DM", b, d, dblocklists[d][b]);
        unsigned int bid = bit->second;
        (*ddlm->decomposition)[d].insert(bid);
      }
    }
  }
  else { // Empty splits indicate default: split all blocks with one per split.
    PetscInt d = 0;
    for (const auto & pr : (*ddlm->blockids)) {
      ddlm->decomposition->push_back(std::set<unsigned int>());
      unsigned int bid = pr.second;
      (*ddlm->decomposition)[d].insert(bid);
      ++d;
    }
  }
  LibmeshPetscCallQ(DMlibMeshSetUpName_Private(*ddm));
  LibmeshPetscCallQ(DMSetFromOptions(*ddm));
  LibmeshPetscCallQ(DMSetUp(*ddm));
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References b, DM_libMesh::blockids, DM_libMesh::blocknames, DM_libMesh::decomposition, DM_libMesh::decomposition_type, DMlibMeshSetUpName_Private(), libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), DM_libMesh::sys, DM_libMesh::varids, and DM_libMesh::varnames.

DMlibMeshCreateFieldDecompositionDM()

PetscErrorCode DMlibMeshCreateFieldDecompositionDM	(	DM	dm,
		PetscInt	dnumber,
		PetscInt *	dsizes,
		char ***	dvarlists,
		DM *	ddm
	)

Definition at line 482 of file petscdmlibmeshimpl.C.

{
  PetscBool islibmesh;
  PetscFunctionBegin;
  PetscValidHeaderSpecific(dm,DM_CLASSID,1);
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH,&islibmesh));
  if (!islibmesh) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Got DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
  if (dnumber < 0) LIBMESH_SETERRQ1(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Negative number %" LIBMESH_PETSCINT_FMT " of decomposition parts", dnumber);
#if PETSC_RELEASE_GREATER_EQUALS(3,20,0)
  PetscAssertPointer(ddm,5);
#else
  PetscValidPointer(ddm,5);
#endif
  DM_libMesh * dlm = (DM_libMesh *)(dm->data);
  LibmeshPetscCallQ(DMCreate(((PetscObject)dm)->comm, ddm));
  LibmeshPetscCallQ(DMSetType(*ddm, DMLIBMESH));
  DM_libMesh * ddlm = (DM_libMesh *)((*ddm)->data);
  ddlm->sys = dlm->sys;
  ddlm->varids = dlm->varids;
  ddlm->varnames = dlm->varnames;
  ddlm->blockids = dlm->blockids;
  ddlm->blocknames = dlm->blocknames;
  ddlm->decomposition = new(std::vector<std::set<unsigned int>>);
  ddlm->decomposition_type = DMLIBMESH_FIELD_DECOMPOSITION;
  if (dnumber) {
    for (PetscInt d = 0; d < dnumber; ++d) {
      if (dsizes[d] < 0) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Negative size %" LIBMESH_PETSCINT_FMT " of decomposition part %" LIBMESH_PETSCINT_FMT, dsizes[d],d);
      ddlm->decomposition->push_back(std::set<unsigned int>());
      for (PetscInt v = 0; v < dsizes[d]; ++v) {
        std::string vname(dvarlists[d][v]);
        std::map<std::string, unsigned int>::const_iterator vit = dlm->varids->find(vname);
        if (vit == dlm->varids->end())
          LIBMESH_SETERRQ3(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Variable %" LIBMESH_PETSCINT_FMT " on the %" LIBMESH_PETSCINT_FMT "-th list with name %s is not owned by this DM", v, d, dvarlists[d][v]);
        unsigned int vid = vit->second;
        (*ddlm->decomposition)[d].insert(vid);
      }
    }
  }
  else { // Empty splits indicate default: split all variables with one per split.
    PetscInt d = 0;
    for (const auto & pr : (*ddlm->varids)) {
      ddlm->decomposition->push_back(std::set<unsigned int>());
      unsigned int vid = pr.second;
      (*ddlm->decomposition)[d].insert(vid);
      ++d;
    }
  }
  LibmeshPetscCallQ(DMlibMeshSetUpName_Private(*ddm));
  LibmeshPetscCallQ(DMSetFromOptions(*ddm));
  LibmeshPetscCallQ(DMSetUp(*ddm));
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::blockids, DM_libMesh::blocknames, DM_libMesh::decomposition, DM_libMesh::decomposition_type, DMlibMeshSetUpName_Private(), libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), DM_libMesh::sys, DM_libMesh::varids, and DM_libMesh::varnames.

DMlibMeshFunction()

static PetscErrorCode DMlibMeshFunction ( DM  dm, Vec  x, Vec  r  )

static

Definition at line 599 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  libmesh_assert(x);
  libmesh_assert(r);
 
  NonlinearImplicitSystem * _sys;
  LibmeshPetscCallQ(DMlibMeshGetSystem(dm, _sys));
  NonlinearImplicitSystem & sys = *_sys;
  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
  PetscVector<Number> & R_sys = *cast_ptr<PetscVector<Number> *>(sys.rhs);
  PetscVector<Number> X_global(x, _sys->comm()), R(r, _sys->comm());
 
  // Use the systems update() to get a good local version of the parallel solution
  X_global.swap(X_sys);
  R.swap(R_sys);
 
  _sys->get_dof_map().enforce_constraints_exactly(*_sys);
  _sys->update();
 
  // Swap back
  X_global.swap(X_sys);
  R.swap(R_sys);
  R.zero();
 
  // if the user has provided both function pointers and objects only the pointer
  // will be used, so catch that as an error
  libmesh_error_msg_if(_sys->nonlinear_solver->residual && _sys->nonlinear_solver->residual_object,
                       "ERROR: cannot specify both a function and object to compute the Residual!");
 
  libmesh_error_msg_if(_sys->nonlinear_solver->matvec && _sys->nonlinear_solver->residual_and_jacobian_object,
                       "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
 
  if (_sys->nonlinear_solver->residual != nullptr)
    _sys->nonlinear_solver->residual(*(_sys->current_local_solution.get()), R, *_sys);
 
  else if (_sys->nonlinear_solver->residual_object != nullptr)
    _sys->nonlinear_solver->residual_object->residual(*(_sys->current_local_solution.get()), R, *_sys);
 
  else if (_sys->nonlinear_solver->matvec   != nullptr)
    _sys->nonlinear_solver->matvec(*(_sys->current_local_solution.get()), &R, nullptr, *_sys);
 
  else if (_sys->nonlinear_solver->residual_and_jacobian_object != nullptr)
    _sys->nonlinear_solver->residual_and_jacobian_object->residual_and_jacobian(*(_sys->current_local_solution.get()), &R, nullptr, *_sys);
 
  else
    libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
 
  R.close();
  X_global.close();
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, DMlibMeshGetSystem(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::System::get_dof_map(), libMesh::libmesh_assert(), libMesh::LibmeshPetscCallQ(), libMesh::NonlinearImplicitSystem::nonlinear_solver, libMesh::PetscFunctionReturn(), libMesh::ExplicitSystem::rhs, libMesh::System::solution, and libMesh::System::update().

Referenced by SNESFunction_DMlibMesh().

DMlibMeshGetBlocks()

PetscErrorCode DMlibMeshGetBlocks	(	DM	dm,
		PetscInt *	n,
		char ***	blocknames
	)

Definition at line 168 of file petscdmlibmeshimpl.C.

{
  PetscInt i;
  PetscFunctionBegin;
  PetscValidHeaderSpecific(dm,DM_CLASSID,1);
  PetscBool islibmesh;
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH,&islibmesh));
  if (!islibmesh) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Got DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
  DM_libMesh * dlm = (DM_libMesh *)(dm->data);
#if PETSC_RELEASE_GREATER_EQUALS(3,20,0)
  PetscAssertPointer(n,2);
#else
  PetscValidPointer(n,2);
#endif
  *n = cast_int<unsigned int>(dlm->blockids->size());
  if (!blocknames) PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
  LibmeshPetscCallQ(PetscMalloc(*n*sizeof(char *), blocknames));
  i = 0;
  for (const auto & pr : *(dlm->blockids))
    {
      LibmeshPetscCallQ(PetscStrallocpy(pr.first.c_str(), *blocknames+i));
      ++i;
    }
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::blockids, libMesh::LibmeshPetscCallQ(), and libMesh::PetscFunctionReturn().

DMlibMeshGetSystem_libMesh()

PetscErrorCode DMlibMeshGetSystem_libMesh	(	DM	dm,
		NonlinearImplicitSystem *&	sys
	)

Definition at line 153 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  PetscValidHeaderSpecific(dm,DM_CLASSID,1);
  PetscBool islibmesh;
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH,&islibmesh));
  if (!islibmesh) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Got DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
  DM_libMesh * dlm = (DM_libMesh *)(dm->data);
  sys = dlm->sys;
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), and DM_libMesh::sys.

Referenced by DMCreate_libMesh().

DMlibMeshGetVariables()

PetscErrorCode DMlibMeshGetVariables	(	DM	dm,
		PetscInt *	n,
		char ***	varnames
	)

Definition at line 196 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  PetscValidHeaderSpecific(dm,DM_CLASSID,1);
  PetscBool islibmesh;
  PetscInt i;
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH,&islibmesh));
  if (!islibmesh) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Got DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
  DM_libMesh * dlm = (DM_libMesh *)(dm->data);
#if PETSC_RELEASE_GREATER_EQUALS(3,20,0)
  PetscAssertPointer(n,2);
#else
  PetscValidPointer(n,2);
#endif
  *n = cast_int<unsigned int>(dlm->varids->size());
  if (!varnames) PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
  LibmeshPetscCallQ(PetscMalloc(*n*sizeof(char *), varnames));
  i = 0;
  for (const auto & pr : *(dlm->varids))
    {
      LibmeshPetscCallQ(PetscStrallocpy(pr.first.c_str(), *varnames+i));
      ++i;
    }
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), and DM_libMesh::varids.

DMlibMeshGetVec_Private()

PetscErrorCode DMlibMeshGetVec_Private	(	DM	,
		const char *	,
		Vec *
	)

Definition at line 80 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
 
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::PetscFunctionReturn().

DMlibMeshJacobian()

static PetscErrorCode DMlibMeshJacobian ( DM  dm, Vec  x, Mat  jac, Mat  pc  )

static

Definition at line 665 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  NonlinearImplicitSystem * _sys;
  LibmeshPetscCallQ(DMlibMeshGetSystem(dm, _sys));
  NonlinearImplicitSystem & sys = *_sys;
 
  libmesh_assert(pc);
  libmesh_assert(jac);
  PetscMatrixBase<Number> & the_pc = *PetscMatrixBase<Number>::get_context(pc, sys.comm());
  PetscMatrixBase<Number> & Jac = *PetscMatrixBase<Number>::get_context(jac, sys.comm());
  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
  PetscMatrixBase<Number> & Jac_sys = *cast_ptr<PetscMatrixBase<Number> *>(sys.matrix);
  PetscVector<Number> X_global(x, sys.comm());
 
  // Set the dof maps
  the_pc.attach_dof_map(sys.get_dof_map());
  Jac.attach_dof_map(sys.get_dof_map());
 
  // Use the systems update() to get a good local version of the parallel solution
  X_global.swap(X_sys);
  Jac.swap(Jac_sys);
 
  sys.get_dof_map().enforce_constraints_exactly(sys);
  sys.update();
 
  X_global.swap(X_sys);
  Jac.swap(Jac_sys);
 
  the_pc.zero();
 
  // if the user has provided both function pointers and objects only the pointer
  // will be used, so catch that as an error
  libmesh_error_msg_if(sys.nonlinear_solver->jacobian && sys.nonlinear_solver->jacobian_object,
                       "ERROR: cannot specify both a function and object to compute the Jacobian!");
 
  libmesh_error_msg_if(sys.nonlinear_solver->matvec && sys.nonlinear_solver->residual_and_jacobian_object,
                       "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
 
  if (sys.nonlinear_solver->jacobian != nullptr)
    sys.nonlinear_solver->jacobian(*(sys.current_local_solution.get()), the_pc, sys);
 
  else if (sys.nonlinear_solver->jacobian_object != nullptr)
    sys.nonlinear_solver->jacobian_object->jacobian(*(sys.current_local_solution.get()), the_pc, sys);
 
  else if (sys.nonlinear_solver->matvec != nullptr)
    sys.nonlinear_solver->matvec(*(sys.current_local_solution.get()), nullptr, &the_pc, sys);
 
  else if (sys.nonlinear_solver->residual_and_jacobian_object != nullptr)
    sys.nonlinear_solver->residual_and_jacobian_object->residual_and_jacobian(*(sys.current_local_solution.get()), nullptr, &the_pc, sys);
 
  else
    libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
 
  the_pc.close();
  Jac.close();
  X_global.close();
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::SparseMatrix< T >::attach_dof_map(), libMesh::PetscMatrixBase< T >::close(), libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, DMlibMeshGetSystem(), libMesh::DofMap::enforce_constraints_exactly(), libMesh::PetscMatrixBase< T >::get_context(), libMesh::System::get_dof_map(), libMesh::libmesh_assert(), libMesh::LibmeshPetscCallQ(), libMesh::ImplicitSystem::matrix, libMesh::NonlinearImplicitSystem::nonlinear_solver, libMesh::PetscFunctionReturn(), libMesh::System::solution, libMesh::PetscMatrixBase< T >::swap(), libMesh::System::update(), and libMesh::SparseMatrix< T >::zero().

Referenced by SNESJacobian_DMlibMesh().

DMlibMeshSetSystem_libMesh()

PetscErrorCode DMlibMeshSetSystem_libMesh	(	DM	dm,
		NonlinearImplicitSystem &	sys
	)

Definition at line 93 of file petscdmlibmeshimpl.C.

{
  const Parallel::Communicator & comm = sys.comm();
 
  PetscFunctionBegin;
  PetscValidHeaderSpecific(dm,DM_CLASSID,1);
  PetscBool islibmesh;
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH,&islibmesh));
  if (!islibmesh) LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "Got DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
 
  if (dm->setupcalled) SETERRQ(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot reset the libMesh system after DM has been set up.");
  DM_libMesh * dlm = (DM_libMesh *)(dm->data);
  dlm->sys =&sys;
  /* Initially populate the sets of active blockids and varids using all of the
     existing blocks/variables (only variables are supported at the moment). */
  DofMap & dofmap = dlm->sys->get_dof_map();
  dlm->varids->clear();
  dlm->varnames->clear();
  for (auto v : make_range(dofmap.n_variables())) {
    std::string vname = dofmap.variable(v).name();
    dlm->varids->insert(std::pair<std::string,unsigned int>(vname,v));
    dlm->varnames->insert(std::pair<unsigned int,std::string>(v,vname));
  }
  const MeshBase & mesh = dlm->sys->get_mesh();
  dlm->blockids->clear();
  dlm->blocknames->clear();
  std::set<subdomain_id_type> blocks;
  /* The following effectively is a verbatim copy of MeshBase::n_subdomains(). */
  // This requires an inspection on every processor
  libmesh_parallel_only(mesh.comm());
  for (const auto & elem : mesh.active_element_ptr_range())
    blocks.insert(elem->subdomain_id());
  // Some subdomains may only live on other processors
  comm.set_union(blocks);
 
  std::set<subdomain_id_type>::iterator bit = blocks.begin();
  std::set<subdomain_id_type>::iterator bend = blocks.end();
  if (bit == bend) SETERRQ(((PetscObject)dm)->comm, PETSC_ERR_PLIB, "No mesh blocks found.");
 
  for (; bit != bend; ++bit) {
    subdomain_id_type bid = *bit;
    std::string bname = mesh.subdomain_name(bid);
    if (!bname.length()) {
      /* Block names are currently implemented for Exodus II meshes
         only, so we might have to make up our own block names and
         maintain our own mapping of block ids to names.
      */
      std::ostringstream ss;
      ss << "dm" << bid;
      bname = ss.str();
    }
    dlm->blockids->insert(std::pair<std::string,unsigned int>(bname,bid));
    dlm->blocknames->insert(std::pair<unsigned int,std::string>(bid,bname));
  }
  LibmeshPetscCallQ(DMlibMeshSetUpName_Private(dm));
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::blockids, DM_libMesh::blocknames, libMesh::ParallelObject::comm(), DMlibMeshSetUpName_Private(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::LibmeshPetscCallQ(), libMesh::make_range(), mesh, libMesh::DofMap::n_variables(), libMesh::Variable::name(), libMesh::PetscFunctionReturn(), libMesh::Parallel::Communicator::set_union(), libMesh::MeshBase::subdomain_name(), DM_libMesh::sys, libMesh::DofMap::variable(), DM_libMesh::varids, and DM_libMesh::varnames.

Referenced by DMCreate_libMesh().

DMlibMeshSetUpName_Private()

PetscErrorCode DMlibMeshSetUpName_Private

(

DM

dm

)

Definition at line 224 of file petscdmlibmeshimpl.C.

{
  DM_libMesh * dlm = (DM_libMesh *)dm->data;
  PetscFunctionBegin;
  std::string name = dlm->sys->name();
  std::map<unsigned int, std::string> * dnames = LIBMESH_PETSC_NULLPTR,
                                      * enames = LIBMESH_PETSC_NULLPTR;
  if (dlm->decomposition_type == DMLIBMESH_FIELD_DECOMPOSITION) {
    name += ":dec:var:";
    dnames = dlm->varnames;
  }
  if (dlm->decomposition_type == DMLIBMESH_DOMAIN_DECOMPOSITION) {
    name += ":dec:block:";
    dnames = dlm->blocknames;
  }
  if (dnames) {
    for (auto decomp : *dlm->decomposition) {
      for (std::set<unsigned int>::iterator dit_begin = decomp.begin(),
                                            dit = dit_begin,
                                            dit_end   = decomp.end();
           dit != dit_end; ++dit) {
        unsigned int id = *dit;
        if (dit != dit_begin)
          name += ",";
        name += (*dnames)[id];
      }
      name += ";";
    }
  }
  if (dlm->embedding_type == DMLIBMESH_FIELD_EMBEDDING) {
    name += ":emb:var:";
    enames = dlm->varnames;
  }
  if (dlm->embedding_type == DMLIBMESH_DOMAIN_EMBEDDING) {
    name += ":emb:block:";
    enames = dlm->blocknames;
  }
  if (enames) {
    for (auto eit     = enames->begin(),
              eit_end = enames->end(); eit != eit_end; ++eit) {
      std::string & ename = eit->second;
      if (eit != enames->begin())
        name += ",";
      name += ename;
    }
    name += ";";
  }
  LibmeshPetscCallQ(PetscObjectSetName((PetscObject)dm, name.c_str()));
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::blocknames, DM_libMesh::decomposition, DM_libMesh::decomposition_type, DM_libMesh::embedding_type, libMesh::LibmeshPetscCallQ(), libMesh::System::name(), libMesh::PetscFunctionReturn(), DM_libMesh::sys, and DM_libMesh::varnames.

Referenced by DMCreateDomainDecomposition_libMesh(), DMCreateFieldDecomposition_libMesh(), DMlibMeshCreateDomainDecompositionDM(), DMlibMeshCreateFieldDecompositionDM(), and DMlibMeshSetSystem_libMesh().

DMSetUp_libMesh()

static PetscErrorCode DMSetUp_libMesh ( DM  dm )

static

Definition at line 888 of file petscdmlibmeshimpl.C.

{
  PetscFunctionBegin;
  DM_libMesh     * dlm = (DM_libMesh *)(dm->data);
  PetscBool eq;
 
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)dm, DMLIBMESH, &eq));
 
  if (!eq)
    LIBMESH_SETERRQ2(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "DM of type %s, not of type %s", ((PetscObject)dm)->type_name, DMLIBMESH);
 
  if (!dlm->sys)
    SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "No libMesh system set for DM_libMesh");
  /*
    Do not evaluate function, Jacobian or bounds for an embedded DM -- the subproblem might not have enough information for that.
  */
  if (!dlm->embedding) {
    LibmeshPetscCallQ(DMSNESSetFunction(dm, SNESFunction_DMlibMesh, (void *)dm));
    LibmeshPetscCallQ(DMSNESSetJacobian(dm, SNESJacobian_DMlibMesh, (void *)dm));
    if (dlm->sys->nonlinear_solver->bounds || dlm->sys->nonlinear_solver->bounds_object)
      LibmeshPetscCallQ(DMSetVariableBounds(dm, DMVariableBounds_libMesh));
  }
  else {
    /*
      Fow now we don't implement even these, although a linear "Dirichlet" subproblem is well-defined.
      Creating the submatrix, however, might require extracting the submatrix preallocation from an unassembled matrix.
    */
    dm->ops->createglobalvector = 0;
    dm->ops->creatematrix = 0;
  }
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DMVariableBounds_libMesh(), DM_libMesh::embedding, libMesh::LibmeshPetscCallQ(), libMesh::NonlinearImplicitSystem::nonlinear_solver, libMesh::PetscFunctionReturn(), SNESFunction_DMlibMesh(), SNESJacobian_DMlibMesh(), and DM_libMesh::sys.

Referenced by DMCreate_libMesh().

DMVariableBounds_libMesh()

static PetscErrorCode DMVariableBounds_libMesh ( DM  dm, Vec  xl, Vec  xu  )

static

Definition at line 737 of file petscdmlibmeshimpl.C.

{
  NonlinearImplicitSystem * _sys;
  LibmeshPetscCallQ(DMlibMeshGetSystem(dm, _sys));
  NonlinearImplicitSystem & sys = *_sys;
  PetscVector<Number> XL(xl, sys.comm());
  PetscVector<Number> XU(xu, sys.comm());
  PetscFunctionBegin;
  // Workaround for nonstandard Q suffix warning with quad precision
  const PetscReal petsc_inf = std::numeric_limits<PetscReal>::max() / 4;
  LibmeshPetscCallQ(VecSet(xl, -petsc_inf));
  LibmeshPetscCallQ(VecSet(xu, petsc_inf));
  if (sys.nonlinear_solver->bounds != nullptr)
    sys.nonlinear_solver->bounds(XL,XU,sys);
  else if (sys.nonlinear_solver->bounds_object != nullptr)
    sys.nonlinear_solver->bounds_object->bounds(XL,XU, sys);
  else
    SETERRQ(((PetscObject)dm)->comm, PETSC_ERR_ARG_WRONG, "No bounds calculation in this libMesh object");
 
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::ParallelObject::comm(), DMlibMeshGetSystem(), libMesh::LibmeshPetscCallQ(), libMesh::NonlinearImplicitSystem::nonlinear_solver, and libMesh::PetscFunctionReturn().

Referenced by DMSetUp_libMesh().

DMView_libMesh()

static PetscErrorCode DMView_libMesh ( DM  dm, PetscViewer  viewer  )

static

Definition at line 830 of file petscdmlibmeshimpl.C.

{
  PetscBool isascii;
  const char * name, * prefix;
  DM_libMesh * dlm = (DM_libMesh *)dm->data;
  PetscFunctionBegin;
  LibmeshPetscCallQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
  if (isascii) {
    LibmeshPetscCallQ(PetscObjectGetName((PetscObject)dm, &name));
    LibmeshPetscCallQ(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix));
    LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "DM libMesh with name %s and prefix %s\n", name, prefix));
    LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "blocks:"));
    std::map<std::string,unsigned int>::iterator bit = dlm->blockids->begin();
    std::map<std::string,unsigned int>::const_iterator bend = dlm->blockids->end();
    for (; bit != bend; ++bit) {
      LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "(%s,%d) ", bit->first.c_str(), bit->second));
    }
    LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "\n"));
    LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "variables:"));
    std::map<std::string,unsigned int>::iterator vit = dlm->varids->begin();
    std::map<std::string,unsigned int>::const_iterator vend = dlm->varids->end();
    for (; vit != vend; ++vit) {
      LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "(%s,%d) ", vit->first.c_str(), vit->second));
    }
    LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "\n"));
    if (dlm->decomposition_type == DMLIBMESH_NO_DECOMPOSITION) {
      LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "No decomposition\n"));
    }
    else {
      if (dlm->decomposition_type == DMLIBMESH_FIELD_DECOMPOSITION) {
        LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "Field decomposition by variable: "));
      }
      else if (dlm->decomposition_type == DMLIBMESH_DOMAIN_DECOMPOSITION) {
        LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "Domain decomposition by block: "));
      }
      else LIBMESH_SETERRQ1(((PetscObject)dm)->comm, PETSC_ERR_PLIB, "Unexpected decomposition type: %d", dlm->decomposition_type);
      /* FIX: decompositions might have different sizes and components on different ranks. */
      for (auto d : index_range(*dlm->decomposition)) {
        std::set<unsigned int>::iterator dbegin  = (*dlm->decomposition)[d].begin();
        std::set<unsigned int>::iterator dit     = (*dlm->decomposition)[d].begin();
        std::set<unsigned int>::iterator dend    = (*dlm->decomposition)[d].end();
        for (; dit != dend; ++dit) {
          if (dit != dbegin) {
            LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, ","));
          }
          LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "%u", *dit));
        }
        LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, ";"));
      }
      LibmeshPetscCallQ(PetscViewerASCIIPrintf(viewer, "\n"));
    }
  }
 
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References DM_libMesh::blockids, DM_libMesh::decomposition, DM_libMesh::decomposition_type, libMesh::index_range(), libMesh::LibmeshPetscCallQ(), libMesh::PetscFunctionReturn(), and DM_libMesh::varids.

Referenced by DMCreate_libMesh().

SNESFunction_DMlibMesh()

static PetscErrorCode SNESFunction_DMlibMesh ( SNES  , Vec  x, Vec  r, void *  ctx  )

static

Definition at line 654 of file petscdmlibmeshimpl.C.

{
  DM dm = (DM)ctx;
  PetscFunctionBegin;
  LibmeshPetscCallQ(DMlibMeshFunction(dm,x,r));
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::ctx, DMlibMeshFunction(), libMesh::LibmeshPetscCallQ(), and libMesh::PetscFunctionReturn().

Referenced by DMSetUp_libMesh().

SNESJacobian_DMlibMesh()

static PetscErrorCode SNESJacobian_DMlibMesh ( SNES  , Vec  x, Mat  jac, Mat  pc, void *  ctx  )

static

Definition at line 727 of file petscdmlibmeshimpl.C.

{
  DM dm = (DM)ctx;
  PetscFunctionBegin;
  LibmeshPetscCallQ(DMlibMeshJacobian(dm,x,jac,pc));
  PetscFunctionReturn(LIBMESH_PETSC_SUCCESS);
}

References libMesh::ctx, DMlibMeshJacobian(), libMesh::LibmeshPetscCallQ(), and libMesh::PetscFunctionReturn().

Referenced by DMSetUp_libMesh().