-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathperiodicMeshOutput.c
103 lines (88 loc) · 3.89 KB
/
periodicMeshOutput.c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
static char help[] = "Tests for periodic mesh output\n\n";
#include <petscdmplex.h>
PetscErrorCode CheckMesh(DM dm) {
PetscReal detJ, J[9];
PetscReal vol;
PetscInt dim, depth, cStart, cEnd, c;
PetscFunctionBegin;
PetscCall(DMGetDimension(dm, &dim));
PetscCall(DMPlexGetDepth(dm, &depth));
PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
for (c = cStart; c < cEnd; ++c) {
PetscCall(DMPlexComputeCellGeometryFEM(dm, c, NULL, NULL, J, NULL, &detJ));
PetscCheck(detJ > 0.0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG,
"Mesh cell %" PetscInt_FMT " is inverted, |J| = %g", c, (double) detJ);
if (depth > 1) {
PetscCall(DMPlexComputeCellGeometryFVM(dm, c, &vol, NULL, NULL));
PetscCheck(vol > 0.0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG,
"Mesh cell %" PetscInt_FMT " is inverted, vol = %g", c, (double) vol);
}
}
PetscFunctionReturn(0);
}
PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm) {
PetscFunctionBegin;
PetscCall(DMCreate(comm, dm));
PetscCall(DMSetType(*dm, DMPLEX));
PetscCall(DMSetFromOptions(*dm));
PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
PetscFunctionReturn(0);
// PetscFunctionBegin;
// PetscInt faces[2] = {3, 3};
// PetscReal lower[2] = {0.0, 0.0};
// PetscReal upper[2] = {1.0, 1.0};
// DMBoundaryType boundaryTypes[2] = {DM_BOUNDARY_PERIODIC, DM_BOUNDARY_PERIODIC};
//
// PetscCall(DMPlexCreateBoxMesh(comm, 2, PETSC_FALSE, faces, lower, upper, boundaryTypes, PETSC_TRUE, dm));
//// PetscCall(DMSetBasicAdjacency(*dm, PETSC_TRUE, PETSC_FALSE));
// PetscCall(DMSetFromOptions(*dm));
// PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
// PetscFunctionReturn(0);
}
int main(int argc, char **argv) {
DM dm;
PetscFunctionBeginUser;
PetscCall(PetscInitialize(&argc, &argv, NULL, help));
PetscCall(CreateMesh(PETSC_COMM_WORLD, &dm));
PetscCall(CheckMesh(dm));
PetscCall(DMDestroy(&dm));
PetscCall(PetscFinalize());
return 0;
}
/*TEST
test:
suffix: 0
args: -dm_plex_simplex 0 -dm_plex_box_faces 3,1,0 -dm_plex_box_bd periodic,none -dm_view ::ascii_info_detail
test:
suffix: 1
nsize: 2
args: -dm_plex_simplex 0 -dm_plex_box_faces 3,1,0 -dm_plex_box_bd periodic,none -petscpartitioner_type simple -dm_view ::ascii_info_detail
test:
suffix: 2
nsize: 2
args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -petscpartitioner_type simple -dm_view ::ascii_info_detail
test:
suffix: 3
nsize: 4
args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -petscpartitioner_type simple -dm_view ::ascii_info_detail
test:
suffix: 4
nsize: 2
args: -dm_plex_simplex 0 -dm_plex_box_faces 3,1,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view ::ascii_info_detail
test:
suffix: 5
nsize: 2
args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view ::ascii_info_detail
# This checks that the SF with extra root for periodic cut still checks
test:
suffix: 5_hdf5
requires: hdf5
nsize: 2
args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view hdf5:mesh.h5
test:
suffix: 6
nsize: 4
args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view ::ascii_info_detail
-dm_plex_simplex 0 -dm_plex_box_faces 3,3,0 -dm_plex_box_bd periodic,periodic -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view hdf5:mesh.h5
-dm_plex_simplex 0 -dm_plex_box_faces 3,3,0 -dm_plex_box_bd periodic,periodic -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view :mesh.tex:ascii_latex
TEST*/