Actual source code: ex1.c
1: const char help[] = "Test basic creation and evaluation of PETSCSPACEPTRIMMED";
3: #include <petscfe.h>
5: static PetscErrorCode test(PetscInt dim, PetscInt formDegree, PetscInt degree, PetscInt nCopies)
6: {
7: MPI_Comm comm = PETSC_COMM_SELF;
8: PetscSpace sp;
9: PetscInt Nf, Nb;
10: PetscInt maxDexp, maxD, d;
11: PetscInt Nbexp, Bsize, Dsize, Hsize;
12: PetscReal *B, *D, *H;
13: PetscQuadrature quad;
14: PetscInt npoints;
15: const PetscReal *points;
17: PetscSpaceCreate(comm, &sp);
18: PetscObjectSetName((PetscObject)sp, "ptrimmed");
19: PetscSpaceSetType(sp, PETSCSPACEPTRIMMED);
20: PetscSpaceSetNumVariables(sp, dim);
21: PetscDTBinomialInt(dim, PetscAbsInt(formDegree), &Nf);
22: PetscSpaceSetNumComponents(sp, Nf * nCopies);
23: PetscSpaceSetDegree(sp, degree, PETSC_DETERMINE);
24: PetscSpacePTrimmedSetFormDegree(sp, formDegree);
25: PetscSpaceSetUp(sp);
26: PetscSpaceView(sp, NULL);
28: PetscDTPTrimmedSize(dim, formDegree == 0 ? degree : degree + 1, PetscAbsInt(formDegree), &Nbexp);
29: Nbexp *= nCopies;
30: PetscSpaceGetDimension(sp, &Nb);
33: maxDexp = (PetscAbsInt(formDegree) == dim || formDegree == 0) ? degree : degree + 1;
34: PetscSpaceGetDegree(sp, &d, &maxD);
38: PetscDTStroudConicalQuadrature(dim, 1, maxD + 1, -1., 1., &quad);
39: PetscQuadratureGetData(quad, NULL, NULL, &npoints, &points, NULL);
41: Bsize = npoints * Nb * Nf * nCopies;
42: Dsize = dim * Bsize;
43: Hsize = dim * Dsize;
44: PetscMalloc3(Bsize, &B, Dsize, &D, Hsize, &H);
45: PetscSpaceEvaluate(sp, npoints, points, B, D, H);
49: PetscFree3(B, D, H);
50: PetscQuadratureDestroy(&quad);
51: PetscSpaceDestroy(&sp);
52: return 0;
53: }
55: int main(int argc, char *argv[])
56: {
58: PetscInitialize(&argc, &argv, NULL, help);
59: for (PetscInt dim = 0; dim <= 3; dim++) {
60: for (PetscInt formDegree = -dim; formDegree <= dim; formDegree++) {
61: for (PetscInt degree = 0; degree <= 4; degree++) {
62: if (formDegree == 0 && degree == 0) continue;
63: for (PetscInt nCopies = 1; nCopies <= PetscMax(2, dim); nCopies++) test(dim, formDegree, degree, nCopies);
64: }
65: }
66: }
67: PetscFinalize();
68: return 0;
69: }
71: /*TEST
73: test:
75: TEST*/