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