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: {
22: if (fe->ceed == ceed) return 0;
23: CeedReferenceCopy(ceed, &fe->ceed);
24: return 0;
25: }
27: /*@C
28: PetscFEGetCeedBasis - Get the Ceed object mirroring this FE
30: Not Collective
32: Input Parameter:
33: . fe - The PetscFE
35: Output Parameter:
36: . basis - The CeedBasis
38: Note: This is a borrowed reference, so it is not freed.
40: Level: intermediate
42: .seealso: PetscFESetCeed(), DMGetCeed()
43: @*/
44: PetscErrorCode PetscFEGetCeedBasis(PetscFE fe, CeedBasis *basis)
45: {
46: PetscSpace sp;
47: PetscQuadrature q;
48: PetscInt dim, Nc, deg, ord;
52: if (!fe->ceedBasis && fe->ceed) {
53: PetscFEGetSpatialDimension(fe, &dim);
54: PetscFEGetNumComponents(fe, &Nc);
55: PetscFEGetBasisSpace(fe, &sp);
56: PetscSpaceGetDegree(sp, °, NULL);
57: PetscFEGetQuadrature(fe, &q);
58: PetscQuadratureGetOrder(q, &ord);
59: CeedBasisCreateTensorH1Lagrange(fe->ceed, dim, Nc, deg+1, (ord+1)/2, CEED_GAUSS, &fe->ceedBasis);
60: }
61: *basis = fe->ceedBasis;
62: return 0;
63: }
65: #endif