1f312b944SMark Adams const char help[] = "Simple example to get equally space points in high-order elements (and XGC mirror)";
2f312b944SMark Adams
3f312b944SMark Adams #include <petscfe.h>
4f312b944SMark Adams #include <petscdmplex.h>
x(PetscInt dim,PetscReal time,const PetscReal x[],PetscInt Nf_unused,PetscScalar * u,void * actx)5*2a8381b2SBarry Smith static PetscErrorCode x(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf_unused, PetscScalar *u, void *actx)
6d71ae5a4SJacob Faibussowitsch {
7f312b944SMark Adams PetscFunctionBegin;
8f312b944SMark Adams u[0] = x[0];
93ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
10f312b944SMark Adams }
11f312b944SMark Adams
main(int argc,char ** argv)12d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv)
13d71ae5a4SJacob Faibussowitsch {
14f312b944SMark Adams PetscInt dim = 2, cells[] = {1, 1, 1};
15f312b944SMark Adams DM K;
16f312b944SMark Adams PetscReal radius = 2, lo[] = {-radius, -radius, -radius}, hi[] = {radius, radius, radius};
17f312b944SMark Adams DMBoundaryType periodicity[3] = {DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE};
18f312b944SMark Adams PetscFE fe;
19f312b944SMark Adams PetscErrorCode (*initu[1])(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *);
20f312b944SMark Adams Vec X;
21f312b944SMark Adams
22f312b944SMark Adams PetscFunctionBeginUser;
23f312b944SMark Adams PetscCall(PetscInitialize(&argc, &argv, NULL, help));
24f312b944SMark Adams PetscOptionsBegin(PETSC_COMM_WORLD, "", "Options for PETSCDUALSPACELAGRANGE test", "none");
25f312b944SMark Adams PetscCall(PetscOptionsRangeInt("-dim", "The spatial dimension", "ex1.c", dim, &dim, NULL, 0, 3));
26f312b944SMark Adams PetscOptionsEnd();
27f312b944SMark Adams
2842108689Sksagiyam PetscCall(DMPlexCreateBoxMesh(PETSC_COMM_SELF, dim, PETSC_FALSE, cells, lo, hi, periodicity, PETSC_TRUE, 0, PETSC_TRUE, &K));
29f312b944SMark Adams PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, PETSC_FALSE, NULL, PETSC_DECIDE, &fe));
30f312b944SMark Adams PetscCall(DMSetField(K, 0, NULL, (PetscObject)fe));
31f312b944SMark Adams PetscCall(PetscFEDestroy(&fe));
32f312b944SMark Adams PetscCall(DMCreateDS(K));
33f312b944SMark Adams
34f312b944SMark Adams initu[0] = x;
35f312b944SMark Adams PetscCall(DMCreateGlobalVector(K, &X));
36f312b944SMark Adams PetscCall(DMProjectFunction(K, 0.0, initu, NULL, INSERT_ALL_VALUES, X));
37f312b944SMark Adams PetscCall(DMViewFromOptions(K, NULL, "-dual_dm_view"));
38f312b944SMark Adams PetscCall(VecViewFromOptions(X, NULL, "-dual_vec_view"));
39f312b944SMark Adams PetscCall(DMDestroy(&K));
40f312b944SMark Adams PetscCall(VecDestroy(&X));
41f312b944SMark Adams PetscCall(PetscFinalize());
42f312b944SMark Adams return 0;
43f312b944SMark Adams }
44f312b944SMark Adams
45f312b944SMark Adams /*TEST
46f312b944SMark Adams
47f312b944SMark Adams testset:
48f312b944SMark Adams filter: grep -v DM_
49f312b944SMark Adams diff_args: -j
50f312b944SMark Adams args: -petscspace_degree 4 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1 -dual_dm_view -dual_vec_view
51f312b944SMark Adams test:
52f312b944SMark Adams requires: !single
53f312b944SMark Adams suffix: 0
54f312b944SMark Adams test:
55f312b944SMark Adams requires: !single
56f312b944SMark Adams suffix: 3d
57f312b944SMark Adams args: -dim 3
58f312b944SMark Adams
59f312b944SMark Adams TEST*/
60