Actual source code: feceed.c

  1: #include <petsc/private/petscfeimpl.h>

  3: #ifdef PETSC_HAVE_LIBCEED
  4: #include <petscfeceed.h>

  6: /*@C
  7:   PetscFESetCeed - Set the Ceed object

  9:   Not Collective

 11:   Input Parameters:
 12: + fe   - The PetscFE
 13: - ceed - The Ceed object

 15:   Level: intermediate

 17: .seealso: PetscFEGetCeedBasis(), DMGetCeed()
 18: @*/
 19: PetscErrorCode PetscFESetCeed(PetscFE fe, Ceed ceed)
 20: {

 25:   if (fe->ceed == ceed) return(0);
 26:   CeedReferenceCopy(ceed, &fe->ceed);
 27:   return(0);
 28: }

 30: /*@C
 31:   PetscFEGetCeedBasis - Get the Ceed object mirroring this FE

 33:   Not Collective

 35:   Input Parameter:
 36: . fe - The PetscFE

 38:   Output Parameter:
 39: . basis - The CeedBasis

 41:   Note: This is a borrowed reference, so it is not freed.

 43:   Level: intermediate

 45: .seealso: PetscFESetCeed(), DMGetCeed()
 46: @*/
 47: PetscErrorCode PetscFEGetCeedBasis(PetscFE fe, CeedBasis *basis)
 48: {
 49:   PetscSpace      sp;
 50:   PetscQuadrature q;
 51:   PetscInt        dim, Nc, deg, ord;
 52:   PetscErrorCode  ierr;

 57:   if (!fe->ceedBasis && fe->ceed) {
 58:     PetscFEGetSpatialDimension(fe, &dim);
 59:     PetscFEGetNumComponents(fe, &Nc);
 60:     PetscFEGetBasisSpace(fe, &sp);
 61:     PetscSpaceGetDegree(sp, &deg, NULL);
 62:     PetscFEGetQuadrature(fe, &q);
 63:     PetscQuadratureGetOrder(q, &ord);
 64:     CeedBasisCreateTensorH1Lagrange(fe->ceed, dim, Nc, deg+1, (ord+1)/2, CEED_GAUSS, &fe->ceedBasis);
 65:   }
 66:   *basis = fe->ceedBasis;
 67:   return(0);
 68: }

 70: #endif