1 #include <petsc/private/petscfvimpl.h> /*I "petscfv.h" I*/
2
3 #include <petscfvceed.h>
4
5 /*@C
6 PetscFVSetCeed - Set the `Ceed` object to a `PetscFV`
7
8 Not Collective
9
10 Input Parameters:
11 + fv - The `PetscFV`
12 - ceed - The `Ceed` object
13
14 Level: intermediate
15
16 .seealso: `PetscFV`, `PetscFVGetCeedBasis()`, `DMGetCeed()`
17 @*/
PetscFVSetCeed(PetscFV fv,Ceed ceed)18 PetscErrorCode PetscFVSetCeed(PetscFV fv, Ceed ceed)
19 {
20 PetscFunctionBegin;
21 PetscValidHeaderSpecific(fv, PETSCFV_CLASSID, 1);
22 if (fv->ceed == ceed) PetscFunctionReturn(PETSC_SUCCESS);
23 PetscCallCEED(CeedReferenceCopy(ceed, &fv->ceed));
24 PetscFunctionReturn(PETSC_SUCCESS);
25 }
26
27 /*@C
28 PetscFVGetCeedBasis - Get the `Ceed` object mirroring this `PetscFV`
29
30 Not Collective
31
32 Input Parameter:
33 . fv - The `PetscFV`
34
35 Output Parameter:
36 . basis - The `CeedBasis`
37
38 Level: intermediate
39
40 Note:
41 This is a borrowed reference, so it is not freed.
42
43 .seealso: `PetscFV`, `PetscFVSetCeed()`, `DMGetCeed()`
44 @*/
PetscFVGetCeedBasis(PetscFV fv,CeedBasis * basis)45 PetscErrorCode PetscFVGetCeedBasis(PetscFV fv, CeedBasis *basis)
46 {
47 PetscQuadrature q;
48 PetscInt dim, Nc, ord;
49
50 PetscFunctionBegin;
51 PetscValidHeaderSpecific(fv, PETSCFV_CLASSID, 1);
52 PetscAssertPointer(basis, 2);
53 if (!fv->ceedBasis && fv->ceed) {
54 PetscCall(PetscFVGetSpatialDimension(fv, &dim));
55 PetscCall(PetscFVGetNumComponents(fv, &Nc));
56 PetscCall(PetscFVGetQuadrature(fv, &q));
57 PetscCall(PetscQuadratureGetOrder(q, &ord));
58 PetscCallCEED(CeedBasisCreateTensorH1Lagrange(fv->ceed, dim, Nc, 1, (ord + 1) / 2, CEED_GAUSS, &fv->ceedBasis));
59 }
60 *basis = fv->ceedBasis;
61 PetscFunctionReturn(PETSC_SUCCESS);
62 }
63