Skip to content

Commit 5f0717d

Browse files
jeremyltjeremylt
authored
Merge pull request #1820 from CEED/jeremy/bp-multi-eval-mode
ex - add PETSc BP1+3 and BP2+4 examples
2 parents 9f6b563 + 3c11f1f commit 5f0717d

File tree

9 files changed

+347
-121
lines changed

9 files changed

+347
-121
lines changed

‎examples/petsc/bps.c

Lines changed: 8 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -26,6 +26,8 @@
2626
//
2727
//TESTARGS(name="BP3, tet elements") -ceed {ceed_resource} -test -problem bp3 -degree 3 -ksp_max_it_clip 50,50 -simplex
2828
//TESTARGS(name="BP5, hex elements") -ceed {ceed_resource} -test -problem bp5 -degree 3 -ksp_max_it_clip 18,18
29+
//TESTARGS(name="BP1+3, hex elements") -ceed {ceed_resource} -test -problem bp1_3 -degree 3 -ksp_max_it_clip 18,18
30+
//TESTARGS(name="BP2+4, hex elements") -ceed {ceed_resource} -test -problem bp2_4 -degree 3 -ksp_max_it_clip 18,18
2931

3032
/// @file
3133
/// CEED BPs example using PETSc with DMPlex
@@ -183,9 +185,9 @@ static PetscErrorCode RunWithDM(RunParams rp, DM dm, const char *ceed_resource)
183185
{
184186
PC pc;
185187
PetscCall(KSPGetPC(ksp, &pc));
186-
if (rp->bp_choice == CEED_BP1 || rp->bp_choice == CEED_BP2) {
188+
if (rp->bp_choice == CEED_BP1 || rp->bp_choice == CEED_BP2 || rp->bp_choice == CEED_BP13 || rp->bp_choice == CEED_BP24) {
187189
PetscCall(PCSetType(pc, PCJACOBI));
188-
if (rp->simplex) {
190+
if (rp->simplex || rp->bp_choice == CEED_BP13 || rp->bp_choice == CEED_BP24) {
189191
PetscCall(PCJacobiSetType(pc, PC_JACOBI_DIAGONAL));
190192
} else {
191193
PetscCall(PCJacobiSetType(pc, PC_JACOBI_ROWSUM));
@@ -255,7 +257,10 @@ static PetscErrorCode RunWithDM(RunParams rp, DM dm, const char *ceed_resource)
255257
PetscCall(SetupErrorOperatorCtx(rp->comm, dm, ceed, ceed_data, X_loc, op_error, op_error_ctx));
256258
PetscScalar l2_error;
257259
PetscCall(ComputeL2Error(X, &l2_error, op_error_ctx));
258-
PetscReal tol = 5e-2;
260+
// Tighter tol for BP1, BP2
261+
// Looser tol for BP3, BP4, BP5, and BP6 with extra for vector valued problems
262+
// BP1+3 and BP2+4 follow the pattern for BP3 and BP4
263+
PetscReal tol = rp->bp_choice < CEED_BP3 ? 5e-4 : (5e-2 + (rp->bp_choice % 2 == 1 ? 5e-3 : 0));
259264
if (!rp->test_mode || l2_error > tol) {
260265
PetscCall(MPI_Allreduce(&my_rt, &rt_min, 1, MPI_DOUBLE, MPI_MIN, rp->comm));
261266
PetscCall(MPI_Allreduce(&my_rt, &rt_max, 1, MPI_DOUBLE, MPI_MAX, rp->comm));

‎examples/petsc/bps.h

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -17,4 +17,4 @@ static const char *const mem_types[] = {"host", "device", "memType", "CEED_MEM_"
1717
typedef enum { COARSEN_UNIFORM = 0, COARSEN_LOGARITHMIC = 1 } CoarsenType;
1818
static const char *const coarsen_types[] = {"uniform", "logarithmic", "CoarsenType", "COARSEN", 0};
1919

20-
static const char *const bp_types[] = {"bp1", "bp2", "bp3", "bp4", "bp5", "bp6", "BPType", "CEED_BP", 0};
20+
static const char *const bp_types[] = {"bp1", "bp2", "bp3", "bp4", "bp5", "bp6", "bp1_3", "bp2_4", "BPType", "CEED_BP", 0};

‎examples/petsc/include/bpsproblemdata.h

Lines changed: 139 additions & 103 deletions
Original file line numberDiff line numberDiff line change
@@ -14,7 +14,9 @@
1414

1515
#include "../include/structs.h"
1616
#include "../qfunctions/bps/bp1.h"
17+
#include "../qfunctions/bps/bp13.h"
1718
#include "../qfunctions/bps/bp2.h"
19+
#include "../qfunctions/bps/bp24.h"
1820
#include "../qfunctions/bps/bp3.h"
1921
#include "../qfunctions/bps/bp4.h"
2022
#include "../qfunctions/bps/common.h"
@@ -23,107 +25,141 @@
2325
// BP Option Data
2426
// -----------------------------------------------------------------------------
2527

26-
BPData bp_options[6] = {
27-
[CEED_BP1] = {.num_comp_u = 1,
28-
.num_comp_x = 3,
29-
.topo_dim = 3,
30-
.q_data_size = 1,
31-
.q_extra = 1,
32-
.setup_geo = SetupMassGeo,
33-
.setup_rhs = SetupMassRhs,
34-
.apply = Mass,
35-
.error = Error,
36-
.setup_geo_loc = SetupMassGeo_loc,
37-
.setup_rhs_loc = SetupMassRhs_loc,
38-
.apply_loc = Mass_loc,
39-
.error_loc = Error_loc,
40-
.in_mode = CEED_EVAL_INTERP,
41-
.out_mode = CEED_EVAL_INTERP,
42-
.q_mode = CEED_GAUSS,
43-
.enforce_bc = PETSC_FALSE},
44-
[CEED_BP2] = {.num_comp_u = 3,
45-
.num_comp_x = 3,
46-
.topo_dim = 3,
47-
.q_data_size = 1,
48-
.q_extra = 1,
49-
.setup_geo = SetupMassGeo,
50-
.setup_rhs = SetupMassRhs3,
51-
.apply = Mass3,
52-
.error = Error3,
53-
.setup_geo_loc = SetupMassGeo_loc,
54-
.setup_rhs_loc = SetupMassRhs3_loc,
55-
.apply_loc = Mass3_loc,
56-
.error_loc = Error3_loc,
57-
.in_mode = CEED_EVAL_INTERP,
58-
.out_mode = CEED_EVAL_INTERP,
59-
.q_mode = CEED_GAUSS,
60-
.enforce_bc = PETSC_FALSE},
61-
[CEED_BP3] = {.num_comp_u = 1,
62-
.num_comp_x = 3,
63-
.topo_dim = 3,
64-
.q_data_size = 7,
65-
.q_extra = 1,
66-
.setup_geo = SetupDiffGeo,
67-
.setup_rhs = SetupDiffRhs,
68-
.apply = Diff,
69-
.error = Error,
70-
.setup_geo_loc = SetupDiffGeo_loc,
71-
.setup_rhs_loc = SetupDiffRhs_loc,
72-
.apply_loc = Diff_loc,
73-
.error_loc = Error_loc,
74-
.in_mode = CEED_EVAL_GRAD,
75-
.out_mode = CEED_EVAL_GRAD,
76-
.q_mode = CEED_GAUSS,
77-
.enforce_bc = PETSC_TRUE },
78-
[CEED_BP4] = {.num_comp_u = 3,
79-
.num_comp_x = 3,
80-
.topo_dim = 3,
81-
.q_data_size = 7,
82-
.q_extra = 1,
83-
.setup_geo = SetupDiffGeo,
84-
.setup_rhs = SetupDiffRhs3,
85-
.apply = Diff3,
86-
.error = Error3,
87-
.setup_geo_loc = SetupDiffGeo_loc,
88-
.setup_rhs_loc = SetupDiffRhs3_loc,
89-
.apply_loc = Diff3_loc,
90-
.error_loc = Error3_loc,
91-
.in_mode = CEED_EVAL_GRAD,
92-
.out_mode = CEED_EVAL_GRAD,
93-
.q_mode = CEED_GAUSS,
94-
.enforce_bc = PETSC_TRUE },
95-
[CEED_BP5] = {.num_comp_u = 1,
96-
.num_comp_x = 3,
97-
.topo_dim = 3,
98-
.q_data_size = 7,
99-
.q_extra = 0,
100-
.setup_geo = SetupDiffGeo,
101-
.setup_rhs = SetupDiffRhs,
102-
.apply = Diff,
103-
.error = Error,
104-
.setup_geo_loc = SetupDiffGeo_loc,
105-
.setup_rhs_loc = SetupDiffRhs_loc,
106-
.apply_loc = Diff_loc,
107-
.error_loc = Error_loc,
108-
.in_mode = CEED_EVAL_GRAD,
109-
.out_mode = CEED_EVAL_GRAD,
110-
.q_mode = CEED_GAUSS_LOBATTO,
111-
.enforce_bc = PETSC_TRUE },
112-
[CEED_BP6] = {.num_comp_u = 3,
113-
.num_comp_x = 3,
114-
.topo_dim = 3,
115-
.q_data_size = 7,
116-
.q_extra = 0,
117-
.setup_geo = SetupDiffGeo,
118-
.setup_rhs = SetupDiffRhs3,
119-
.apply = Diff3,
120-
.error = Error3,
121-
.setup_geo_loc = SetupDiffGeo_loc,
122-
.setup_rhs_loc = SetupDiffRhs3_loc,
123-
.apply_loc = Diff3_loc,
124-
.error_loc = Error3_loc,
125-
.in_mode = CEED_EVAL_GRAD,
126-
.out_mode = CEED_EVAL_GRAD,
127-
.q_mode = CEED_GAUSS_LOBATTO,
128-
.enforce_bc = PETSC_TRUE }
28+
BPData bp_options[8] = {
29+
[CEED_BP1] = {.num_comp_u = 1,
30+
.num_comp_x = 3,
31+
.topo_dim = 3,
32+
.q_data_size = 1,
33+
.q_extra = 1,
34+
.setup_geo = SetupMassGeo,
35+
.setup_rhs = SetupMassRhs,
36+
.apply = Mass,
37+
.error = Error,
38+
.setup_geo_loc = SetupMassGeo_loc,
39+
.setup_rhs_loc = SetupMassRhs_loc,
40+
.apply_loc = Mass_loc,
41+
.error_loc = Error_loc,
42+
.in_mode = CEED_EVAL_INTERP,
43+
.out_mode = CEED_EVAL_INTERP,
44+
.q_mode = CEED_GAUSS,
45+
.enforce_bc = PETSC_FALSE},
46+
[CEED_BP2] = {.num_comp_u = 3,
47+
.num_comp_x = 3,
48+
.topo_dim = 3,
49+
.q_data_size = 1,
50+
.q_extra = 1,
51+
.setup_geo = SetupMassGeo,
52+
.setup_rhs = SetupMassRhs3,
53+
.apply = Mass3,
54+
.error = Error3,
55+
.setup_geo_loc = SetupMassGeo_loc,
56+
.setup_rhs_loc = SetupMassRhs3_loc,
57+
.apply_loc = Mass3_loc,
58+
.error_loc = Error3_loc,
59+
.in_mode = CEED_EVAL_INTERP,
60+
.out_mode = CEED_EVAL_INTERP,
61+
.q_mode = CEED_GAUSS,
62+
.enforce_bc = PETSC_FALSE},
63+
[CEED_BP3] = {.num_comp_u = 1,
64+
.num_comp_x = 3,
65+
.topo_dim = 3,
66+
.q_data_size = 7,
67+
.q_extra = 1,
68+
.setup_geo = SetupDiffGeo,
69+
.setup_rhs = SetupDiffRhs,
70+
.apply = Diff,
71+
.error = Error,
72+
.setup_geo_loc = SetupDiffGeo_loc,
73+
.setup_rhs_loc = SetupDiffRhs_loc,
74+
.apply_loc = Diff_loc,
75+
.error_loc = Error_loc,
76+
.in_mode = CEED_EVAL_GRAD,
77+
.out_mode = CEED_EVAL_GRAD,
78+
.q_mode = CEED_GAUSS,
79+
.enforce_bc = PETSC_TRUE },
80+
[CEED_BP4] = {.num_comp_u = 3,
81+
.num_comp_x = 3,
82+
.topo_dim = 3,
83+
.q_data_size = 7,
84+
.q_extra = 1,
85+
.setup_geo = SetupDiffGeo,
86+
.setup_rhs = SetupDiffRhs3,
87+
.apply = Diff3,
88+
.error = Error3,
89+
.setup_geo_loc = SetupDiffGeo_loc,
90+
.setup_rhs_loc = SetupDiffRhs3_loc,
91+
.apply_loc = Diff3_loc,
92+
.error_loc = Error3_loc,
93+
.in_mode = CEED_EVAL_GRAD,
94+
.out_mode = CEED_EVAL_GRAD,
95+
.q_mode = CEED_GAUSS,
96+
.enforce_bc = PETSC_TRUE },
97+
[CEED_BP5] = {.num_comp_u = 1,
98+
.num_comp_x = 3,
99+
.topo_dim = 3,
100+
.q_data_size = 7,
101+
.q_extra = 0,
102+
.setup_geo = SetupDiffGeo,
103+
.setup_rhs = SetupDiffRhs,
104+
.apply = Diff,
105+
.error = Error,
106+
.setup_geo_loc = SetupDiffGeo_loc,
107+
.setup_rhs_loc = SetupDiffRhs_loc,
108+
.apply_loc = Diff_loc,
109+
.error_loc = Error_loc,
110+
.in_mode = CEED_EVAL_GRAD,
111+
.out_mode = CEED_EVAL_GRAD,
112+
.q_mode = CEED_GAUSS_LOBATTO,
113+
.enforce_bc = PETSC_TRUE },
114+
[CEED_BP6] = {.num_comp_u = 3,
115+
.num_comp_x = 3,
116+
.topo_dim = 3,
117+
.q_data_size = 7,
118+
.q_extra = 0,
119+
.setup_geo = SetupDiffGeo,
120+
.setup_rhs = SetupDiffRhs3,
121+
.apply = Diff3,
122+
.error = Error3,
123+
.setup_geo_loc = SetupDiffGeo_loc,
124+
.setup_rhs_loc = SetupDiffRhs3_loc,
125+
.apply_loc = Diff3_loc,
126+
.error_loc = Error3_loc,
127+
.in_mode = CEED_EVAL_GRAD,
128+
.out_mode = CEED_EVAL_GRAD,
129+
.q_mode = CEED_GAUSS_LOBATTO,
130+
.enforce_bc = PETSC_TRUE },
131+
[CEED_BP13] = {.num_comp_u = 1,
132+
.num_comp_x = 3,
133+
.topo_dim = 3,
134+
.q_data_size = 7,
135+
.q_extra = 1,
136+
.setup_geo = SetupDiffGeo,
137+
.setup_rhs = SetupMassDiffRhs,
138+
.apply = MassDiff,
139+
.error = Error,
140+
.setup_geo_loc = SetupDiffGeo_loc,
141+
.setup_rhs_loc = SetupMassDiffRhs_loc,
142+
.apply_loc = MassDiff_loc,
143+
.error_loc = Error_loc,
144+
.in_mode = CEED_EVAL_INTERP + CEED_EVAL_GRAD,
145+
.out_mode = CEED_EVAL_INTERP + CEED_EVAL_GRAD,
146+
.q_mode = CEED_GAUSS,
147+
.enforce_bc = PETSC_TRUE },
148+
[CEED_BP24] = {.num_comp_u = 3,
149+
.num_comp_x = 3,
150+
.topo_dim = 3,
151+
.q_data_size = 7,
152+
.q_extra = 1,
153+
.setup_geo = SetupDiffGeo,
154+
.setup_rhs = SetupMassDiffRhs3,
155+
.apply = MassDiff3,
156+
.error = Error3,
157+
.setup_geo_loc = SetupDiffGeo_loc,
158+
.setup_rhs_loc = SetupMassDiffRhs3_loc,
159+
.apply_loc = MassDiff3_loc,
160+
.error_loc = Error3_loc,
161+
.in_mode = CEED_EVAL_INTERP + CEED_EVAL_GRAD,
162+
.out_mode = CEED_EVAL_INTERP + CEED_EVAL_GRAD,
163+
.q_mode = CEED_GAUSS,
164+
.enforce_bc = PETSC_TRUE },
129165
};

‎examples/petsc/include/structs.h

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -65,7 +65,7 @@ typedef struct {
6565
} BPData;
6666

6767
// BP options
68-
typedef enum { CEED_BP1 = 0, CEED_BP2 = 1, CEED_BP3 = 2, CEED_BP4 = 3, CEED_BP5 = 4, CEED_BP6 = 5 } BPType;
68+
typedef enum { CEED_BP1 = 0, CEED_BP2 = 1, CEED_BP3 = 2, CEED_BP4 = 3, CEED_BP5 = 4, CEED_BP6 = 5, CEED_BP13 = 6, CEED_BP24 = 7 } BPType;
6969

7070
// -----------------------------------------------------------------------------
7171
// Parameter structure for running problems

‎examples/petsc/qfunctions/bps/bp13.h

Lines changed: 74 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,74 @@
1+
// Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and other CEED contributors.
2+
// All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
3+
//
4+
// SPDX-License-Identifier: BSD-2-Clause
5+
//
6+
// This file is part of CEED: http://github.com/ceed
7+
8+
/// @file
9+
/// libCEED QFunctions for diffusion operator example using PETSc
10+
11+
#include <ceed/types.h>
12+
#ifndef CEED_RUNNING_JIT_PASS
13+
#include <math.h>
14+
#endif
15+
16+
// -----------------------------------------------------------------------------
17+
// This QFunction sets up the rhs and true solution for the problem
18+
// -----------------------------------------------------------------------------
19+
CEED_QFUNCTION(SetupMassDiffRhs)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
20+
#ifndef M_PI
21+
#define M_PI 3.14159265358979323846
22+
#endif
23+
const CeedScalar *x = in[0], *w = in[1];
24+
CeedScalar *true_soln = out[0], *rhs = out[1];
25+
26+
// Quadrature Point Loop
27+
CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
28+
const CeedScalar c[3] = {0, 1., 2.};
29+
const CeedScalar k[3] = {1., 2., 3.};
30+
31+
true_soln[i] = sin(M_PI * (c[0] + k[0] * x[i + Q * 0])) * sin(M_PI * (c[1] + k[1] * x[i + Q * 1])) * sin(M_PI * (c[2] + k[2] * x[i + Q * 2]));
32+
33+
rhs[i] = w[i + Q * 0] * (M_PI * M_PI * (k[0] * k[0] + k[1] * k[1] + k[2] * k[2]) + 1.0) * true_soln[i];
34+
} // End of Quadrature Point Loop
35+
return 0;
36+
}
37+
38+
// -----------------------------------------------------------------------------
39+
// This QFunction applies the mass + diffusion operator for a scalar field.
40+
//
41+
// Inputs:
42+
// u - Input vector at quadrature points
43+
// ug - Input vector gradient at quadrature points
44+
// q_data - Geometric factors
45+
//
46+
// Output:
47+
// v - Output vector (test functions) at quadrature points
48+
// vg - Output vector (test functions) gradient at quadrature points
49+
// -----------------------------------------------------------------------------
50+
CEED_QFUNCTION(MassDiff)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
51+
const CeedScalar *u = in[0], *ug = in[1], *q_data = in[2];
52+
CeedScalar *v = out[0], *vg = out[1];
53+
54+
// Quadrature Point Loop
55+
CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
56+
// Read spatial derivatives of u
57+
const CeedScalar du[3] = {ug[i + Q * 0], ug[i + Q * 1], ug[i + Q * 2]};
58+
// Read q_data (dXdxdXdx_T symmetric matrix)
59+
const CeedScalar dXdxdXdx_T[3][3] = {
60+
{q_data[i + 1 * Q], q_data[i + 2 * Q], q_data[i + 3 * Q]},
61+
{q_data[i + 2 * Q], q_data[i + 4 * Q], q_data[i + 5 * Q]},
62+
{q_data[i + 3 * Q], q_data[i + 5 * Q], q_data[i + 6 * Q]}
63+
};
64+
65+
// Mass
66+
v[i] = q_data[i + 0 * Q] * u[i];
67+
// Diff
68+
for (int j = 0; j < 3; j++) { // j = direction of vg
69+
vg[i + j * Q] = (du[0] * dXdxdXdx_T[0][j] + du[1] * dXdxdXdx_T[1][j] + du[2] * dXdxdXdx_T[2][j]);
70+
}
71+
} // End of Quadrature Point Loop
72+
return 0;
73+
}
74+
// -----------------------------------------------------------------------------

0 commit comments

Comments
 (0)