PICurv 0.1.0
A Parallel Particle-In-Cell Solver for Curvilinear LES
Loading...
Searching...
No Matches
test_poisson_rhs.c
Go to the documentation of this file.
1/**
2 * @file test_poisson_rhs.c
3 * @brief C unit tests for Poisson, RHS, body-force, and diffusivity helpers.
4 */
5
6#include "test_support.h"
7
8#include "poisson.h"
9#include "rhs.h"
10#include "verification_sources.h"
11/**
12 * @brief Allocates Poisson/RHS support vectors required by the tests.
13 */
14
15static PetscErrorCode EnsurePoissonAndRhsVectors(UserCtx *user)
16{
17 PetscFunctionBeginUser;
18 if (!user->Phi) PetscCall(DMCreateGlobalVector(user->da, &user->Phi));
19 if (!user->lPhi) PetscCall(DMCreateLocalVector(user->da, &user->lPhi));
20 if (!user->Aj) PetscCall(DMCreateGlobalVector(user->da, &user->Aj));
21 if (!user->lAj) PetscCall(DMCreateLocalVector(user->da, &user->lAj));
22 if (!user->Diffusivity) PetscCall(DMCreateGlobalVector(user->da, &user->Diffusivity));
23 if (!user->lDiffusivity) PetscCall(DMCreateLocalVector(user->da, &user->lDiffusivity));
24 PetscFunctionReturn(0);
25}
26/**
27 * @brief Tests that pressure updates add the correction potential.
28 */
29
30static PetscErrorCode TestUpdatePressureAddsPhi(void)
31{
32 SimCtx *simCtx = NULL;
33 UserCtx *user = NULL;
34
35 PetscFunctionBeginUser;
36 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
37 PetscCall(EnsurePoissonAndRhsVectors(user));
38 PetscCall(VecSet(user->P, 1.0));
39 PetscCall(VecSet(user->Phi, 0.25));
40
41 PetscCall(UpdatePressure(user));
42 PetscCall(PicurvAssertVecConstant(user->P, 1.25, 1.0e-12, "UpdatePressure should add Phi into P"));
43
44 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
45 PetscFunctionReturn(0);
46}
47/**
48 * @brief Tests that the Poisson RHS is zero for zero-divergence velocity.
49 */
50
51static PetscErrorCode TestPoissonRHSZeroDivergence(void)
52{
53 SimCtx *simCtx = NULL;
54 UserCtx *user = NULL;
55 Vec B = NULL;
56
57 PetscFunctionBeginUser;
58 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
59 PetscCall(EnsurePoissonAndRhsVectors(user));
60
61 simCtx->dt = 0.5;
62 PetscCall(VecSet(user->Ucont, 0.0));
63 PetscCall(VecSet(user->Nvert, 0.0));
64 PetscCall(VecSet(user->Aj, 1.0));
65 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
66 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
67 PetscCall(DMGlobalToLocalBegin(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
68 PetscCall(DMGlobalToLocalEnd(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
69 PetscCall(DMGlobalToLocalBegin(user->da, user->Aj, INSERT_VALUES, user->lAj));
70 PetscCall(DMGlobalToLocalEnd(user->da, user->Aj, INSERT_VALUES, user->lAj));
71
72 PetscCall(VecDuplicate(user->P, &B));
73 PetscCall(PoissonRHS(user, B));
74 PetscCall(PicurvAssertVecConstant(B, 0.0, 1.0e-12, "zero velocity divergence should produce zero Poisson RHS"));
75 PetscCall(PicurvAssertRealNear(0.0, simCtx->summationRHS, 1.0e-12, "global Poisson RHS sum should be zero"));
76
77 PetscCall(VecDestroy(&B));
78 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
79 PetscFunctionReturn(0);
80}
81/**
82 * @brief Tests the body-force dispatcher across supported source terms.
83 */
84
85static PetscErrorCode TestComputeBodyForcesDispatcher(void)
86{
87 SimCtx *simCtx = NULL;
88 UserCtx *user = NULL;
89 Vec rct = NULL;
90 Cmpnts ***rct_arr = NULL;
91 Cmpnts ***l_csi = NULL;
92 Cmpnts ***l_eta = NULL;
93 Cmpnts ***l_zet = NULL;
94
95 PetscFunctionBeginUser;
96 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
97 PetscCall(VecDuplicate(user->Ucont, &rct));
98 PetscCall(VecZeroEntries(rct));
99
100 user->boundary_faces[BC_FACE_NEG_X].face_id = BC_FACE_NEG_X;
101 user->boundary_faces[BC_FACE_NEG_X].mathematical_type = PERIODIC;
102 user->boundary_faces[BC_FACE_NEG_X].handler_type = BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX;
103 user->boundary_faces[BC_FACE_POS_X].face_id = BC_FACE_POS_X;
104 user->boundary_faces[BC_FACE_POS_X].mathematical_type = PERIODIC;
105 user->boundary_faces[BC_FACE_POS_X].handler_type = BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX;
106 simCtx->bulkVelocityCorrection = 2.0;
107 simCtx->dt = 1.0;
108 simCtx->forceScalingFactor = 1.0;
109 simCtx->drivingForceMagnitude = 0.0;
110 PetscCall(VecSet(user->Nvert, 0.0));
111 PetscCall(DMGlobalToLocalBegin(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
112 PetscCall(DMGlobalToLocalEnd(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
113 PetscCall(DMDAVecGetArray(user->fda, user->lCsi, &l_csi));
114 PetscCall(DMDAVecGetArray(user->fda, user->lEta, &l_eta));
115 PetscCall(DMDAVecGetArray(user->fda, user->lZet, &l_zet));
116 l_csi[1][1][1].x = 1.0; l_csi[1][1][1].y = 0.0; l_csi[1][1][1].z = 0.0;
117 l_eta[1][1][1].x = 0.0; l_eta[1][1][1].y = 1.0; l_eta[1][1][1].z = 0.0;
118 l_zet[1][1][1].x = 0.0; l_zet[1][1][1].y = 0.0; l_zet[1][1][1].z = 1.0;
119 PetscCall(DMDAVecRestoreArray(user->fda, user->lZet, &l_zet));
120 PetscCall(DMDAVecRestoreArray(user->fda, user->lEta, &l_eta));
121 PetscCall(DMDAVecRestoreArray(user->fda, user->lCsi, &l_csi));
122
123 PetscCall(ComputeBodyForces(user, rct));
124 PetscCall(DMDAVecGetArrayRead(user->fda, rct, &rct_arr));
125 PetscCall(PicurvAssertRealNear(1.5, simCtx->drivingForceMagnitude, 1.0e-12,
126 "ComputeBodyForces should update the driven-flow controller magnitude"));
127 PetscCall(PicurvAssertRealNear(0.0, rct_arr[1][1][1].y, 1.0e-12, "ComputeBodyForces should keep y unchanged"));
128 PetscCall(PicurvAssertRealNear(0.0, rct_arr[1][1][1].z, 1.0e-12, "ComputeBodyForces should keep z unchanged"));
129 PetscCall(DMDAVecRestoreArrayRead(user->fda, rct, &rct_arr));
130
131 PetscCall(VecDestroy(&rct));
132 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
133 PetscFunctionReturn(0);
134}
135/**
136 * @brief Tests Eulerian diffusivity for the molecular-only configuration.
137 */
138
139static PetscErrorCode TestComputeEulerianDiffusivityMolecularOnly(void)
140{
141 SimCtx *simCtx = NULL;
142 UserCtx *user = NULL;
143 PetscReal ***diff_arr = NULL;
144
145 PetscFunctionBeginUser;
146 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
147 PetscCall(EnsurePoissonAndRhsVectors(user));
148
149 simCtx->ren = 2.0;
150 simCtx->schmidt_number = 4.0;
151 simCtx->les = PETSC_FALSE;
152 simCtx->rans = PETSC_FALSE;
153 PetscCall(VecZeroEntries(user->Diffusivity));
154
155 PetscCall(ComputeEulerianDiffusivity(user));
156 PetscCall(DMDAVecGetArrayRead(user->da, user->Diffusivity, &diff_arr));
157 PetscCall(PicurvAssertRealNear(0.125, diff_arr[1][1][1], 1.0e-12, "molecular diffusivity interior value"));
158 PetscCall(DMDAVecRestoreArrayRead(user->da, user->Diffusivity, &diff_arr));
159
160 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
161 PetscFunctionReturn(0);
162}
163/**
164 * @brief Tests that convection vanishes for a quiescent field.
165 */
166
167static PetscErrorCode TestConvectionZeroField(void)
168{
169 SimCtx *simCtx = NULL;
170 UserCtx *user = NULL;
171 Vec conv = NULL;
172
173 PetscFunctionBeginUser;
174 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
175 PetscCall(VecDuplicate(user->lUcont, &conv));
176 PetscCall(VecSet(user->Ucont, 0.0));
177 PetscCall(VecSet(user->Ucat, 0.0));
178 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
179 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
180 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
181 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
182
183 PetscCall(Convection(user, user->lUcont, user->lUcat, conv));
184 PetscCall(PicurvAssertVecConstant(conv, 0.0, 1.0e-12, "Convection should vanish for a quiescent field"));
185
186 PetscCall(VecDestroy(&conv));
187 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
188 PetscFunctionReturn(0);
189}
190/**
191 * @brief Tests that viscous terms vanish for a uniform field.
192 */
193
194static PetscErrorCode TestViscousUniformField(void)
195{
196 SimCtx *simCtx = NULL;
197 UserCtx *user = NULL;
198 Vec visc = NULL;
199
200 PetscFunctionBeginUser;
201 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
202 PetscCall(VecDuplicate(user->lUcont, &visc));
203 PetscCall(VecSet(user->Ucont, 1.0));
204 PetscCall(VecSet(user->Ucat, 1.0));
205 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
206 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
207 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
208 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
209
210 PetscCall(Viscous(user, user->lUcont, user->lUcat, visc));
211 PetscCall(PicurvAssertVecConstant(visc, 0.0, 1.0e-12, "Viscous should vanish for a uniform field"));
212
213 PetscCall(VecDestroy(&visc));
214 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
215 PetscFunctionReturn(0);
216}
217/**
218 * @brief Tests that the full RHS remains zero without forcing on a quiescent field.
219 */
220
221static PetscErrorCode TestComputeRHSZeroFieldNoForcing(void)
222{
223 SimCtx *simCtx = NULL;
224 UserCtx *user = NULL;
225
226 PetscFunctionBeginUser;
227 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
228 PetscCall(VecSet(user->Ucont, 0.0));
229 PetscCall(VecSet(user->Ucat, 0.0));
230 PetscCall(VecSet(user->Nvert, 0.0));
231 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
232 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
233 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
234 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
235 PetscCall(DMGlobalToLocalBegin(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
236 PetscCall(DMGlobalToLocalEnd(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
237
238 PetscCall(ComputeRHS(user, user->Rhs));
239 PetscCall(PicurvAssertVecConstant(user->Rhs, 0.0, 1.0e-12, "ComputeRHS should remain zero without forcing on a quiescent field"));
240
241 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
242 PetscFunctionReturn(0);
243}
244/**
245 * @brief Tests diffusivity-gradient computation on a constant field.
246 */
247
248static PetscErrorCode TestComputeEulerianDiffusivityGradientConstantField(void)
249{
250 SimCtx *simCtx = NULL;
251 UserCtx *user = NULL;
252
253 PetscFunctionBeginUser;
254 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
255 PetscCall(VecSet(user->Diffusivity, 0.75));
256 PetscCall(DMGlobalToLocalBegin(user->da, user->Diffusivity, INSERT_VALUES, user->lDiffusivity));
257 PetscCall(DMGlobalToLocalEnd(user->da, user->Diffusivity, INSERT_VALUES, user->lDiffusivity));
258
259 PetscCall(ComputeEulerianDiffusivityGradient(user));
260 PetscCall(PicurvAssertVecConstant(user->DiffusivityGradient, 0.0, 1.0e-12, "constant diffusivity should yield zero diffusivity gradient"));
261
262 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
263 PetscFunctionReturn(0);
264}
265/**
266 * @brief Tests verification-driven linear diffusivity override and its gradient.
267 */
268
269static PetscErrorCode TestComputeEulerianDiffusivityVerificationLinearX(void)
270{
271 SimCtx *simCtx = NULL;
272 UserCtx *user = NULL;
273 PetscReal ***diff_arr = NULL;
274 Cmpnts ***grad_arr = NULL;
275 Cmpnts ***cent = NULL;
276 const PetscReal gamma0 = 0.5;
277 const PetscReal slope_x = 0.25;
278
279 PetscFunctionBeginUser;
280 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
281 PetscCall(PetscStrncpy(simCtx->eulerianSource, "analytical", sizeof(simCtx->eulerianSource)));
282 simCtx->verificationDiffusivity.enabled = PETSC_TRUE;
283 PetscCall(PetscStrncpy(simCtx->verificationDiffusivity.mode, "analytical", sizeof(simCtx->verificationDiffusivity.mode)));
284 PetscCall(PetscStrncpy(simCtx->verificationDiffusivity.profile, "LINEAR_X", sizeof(simCtx->verificationDiffusivity.profile)));
285 simCtx->verificationDiffusivity.gamma0 = gamma0;
286 simCtx->verificationDiffusivity.slope_x = slope_x;
287
288 PetscCall(PicurvAssertBool(VerificationDiffusivityOverrideActive(simCtx),
289 "verification diffusivity override should report as active"));
290 PetscCall(ComputeEulerianDiffusivity(user));
291 PetscCall(ComputeEulerianDiffusivityGradient(user));
292
293 PetscCall(DMDAVecGetArrayRead(user->da, user->Diffusivity, &diff_arr));
294 PetscCall(DMDAVecGetArrayRead(user->fda, user->Cent, &cent));
295 PetscCall(PicurvAssertRealNear(gamma0 + slope_x * cent[2][2][2].x, diff_arr[2][2][2], 1.0e-12,
296 "verification override should populate the linear diffusivity field"));
297 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Cent, &cent));
298 PetscCall(DMDAVecRestoreArrayRead(user->da, user->Diffusivity, &diff_arr));
299
300 PetscCall(DMDAVecGetArrayRead(user->fda, user->DiffusivityGradient, &grad_arr));
301 PetscCall(PicurvAssertRealNear(slope_x / user->IM, grad_arr[2][2][2].x, 1.0e-12,
302 "linear verification diffusivity should yield the current finite-difference x gradient"));
303 PetscCall(PicurvAssertRealNear(0.0, grad_arr[2][2][2].y, 1.0e-12,
304 "linear verification diffusivity should keep y gradient zero"));
305 PetscCall(PicurvAssertRealNear(0.0, grad_arr[2][2][2].z, 1.0e-12,
306 "linear verification diffusivity should keep z gradient zero"));
307 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->DiffusivityGradient, &grad_arr));
308
309 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
310 PetscFunctionReturn(0);
311}
312/**
313 * @brief Tests that the Poisson null-space operator removes the mean from a constant field.
314 */
315
316static PetscErrorCode TestPoissonNullSpaceFunctionRemovesMean(void)
317{
318 SimCtx *simCtx = NULL;
319 UserCtx *user = NULL;
320 Vec x = NULL;
321 PetscReal sum = 0.0;
322 PetscReal ***x_arr = NULL;
323
324 PetscFunctionBeginUser;
325 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
326 PetscCall(VecDuplicate(user->P, &x));
327 PetscCall(VecSet(x, 3.0));
328 PetscCall(VecSet(user->Nvert, 0.0));
329 PetscCall(DMGlobalToLocalBegin(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
330 PetscCall(DMGlobalToLocalEnd(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
331
332 PetscCall(PoissonNullSpaceFunction(NULL, x, user));
333 PetscCall(VecSum(x, &sum));
334 PetscCall(PicurvAssertRealNear(0.0, sum, 1.0e-10, "PoissonNullSpaceFunction should remove the global mean"));
335 PetscCall(DMDAVecGetArrayRead(user->da, x, &x_arr));
336 PetscCall(PicurvAssertRealNear(0.0, x_arr[1][1][1], 1.0e-10, "PoissonNullSpaceFunction should zero a uniform interior field"));
337 PetscCall(DMDAVecRestoreArrayRead(user->da, x, &x_arr));
338
339 PetscCall(VecDestroy(&x));
340 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
341 PetscFunctionReturn(0);
342}
343/**
344 * @brief Tests that Poisson matrix assembly produces a populated operator on a tiny Cartesian grid.
345 */
346static PetscErrorCode TestPoissonLHSNewAssemblesOperator(void)
347{
348 SimCtx *simCtx = NULL;
349 UserCtx *user = NULL;
350 PetscInt rows = 0, cols = 0;
351 PetscInt ncols = 0;
352 const PetscInt *col_idx = NULL;
353 const PetscScalar *values = NULL;
354 PetscInt interior_row = 0;
355
356 PetscFunctionBeginUser;
357 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
358 PetscCall(VecSet(user->Nvert, 0.0));
359 PetscCall(DMGlobalToLocalBegin(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
360 PetscCall(DMGlobalToLocalEnd(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
361
362 PetscCall(PoissonLHSNew(user));
363 PetscCall(PicurvAssertBool((PetscBool)(user->A != NULL), "PoissonLHSNew should allocate the Poisson operator"));
364
365 PetscCall(MatGetSize(user->A, &rows, &cols));
366 PetscCall(PicurvAssertIntEqual(user->info.mx * user->info.my * user->info.mz, rows, "Poisson operator row count should match the DA node count"));
367 PetscCall(PicurvAssertIntEqual(rows, cols, "Poisson operator should be square"));
368
369 interior_row = (2 * user->info.my + 2) * user->info.mx + 2;
370 PetscCall(MatGetRow(user->A, interior_row, &ncols, &col_idx, &values));
371 PetscCall(PicurvAssertBool((PetscBool)(ncols > 1), "Interior Poisson row should couple to neighboring nodes"));
372 PetscCall(PicurvAssertBool((PetscBool)(values != NULL), "Interior Poisson row should expose non-null coefficients"));
373 PetscCall(MatRestoreRow(user->A, interior_row, &ncols, &col_idx, &values));
374
375 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
376 PetscFunctionReturn(0);
377}
378/**
379 * @brief Tests that projection leaves a zero pressure-correction field unchanged.
380 */
381
382static PetscErrorCode TestProjectionZeroPhiLeavesVelocityUnchanged(void)
383{
384 SimCtx *simCtx = NULL;
385 UserCtx *user = NULL;
386
387 PetscFunctionBeginUser;
388 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
389 PetscCall(VecSet(user->Ucont, 0.0));
390 PetscCall(VecSet(user->Phi, 0.0));
391 PetscCall(VecSet(user->Nvert, 0.0));
392 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
393 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
394 PetscCall(DMGlobalToLocalBegin(user->da, user->Phi, INSERT_VALUES, user->lPhi));
395 PetscCall(DMGlobalToLocalEnd(user->da, user->Phi, INSERT_VALUES, user->lPhi));
396 PetscCall(DMGlobalToLocalBegin(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
397 PetscCall(DMGlobalToLocalEnd(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
398
399 PetscCall(Projection(user));
400 PetscCall(PicurvAssertVecConstant(user->Ucont, 0.0, 1.0e-12, "Projection should leave a zero-velocity field unchanged when Phi is zero"));
401
402 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
403 PetscFunctionReturn(0);
404}
405/**
406 * @brief Tests that projection applies the expected x-direction correction for a linear pressure field.
407 */
408static PetscErrorCode TestProjectionLinearPhiCorrectsVelocity(void)
409{
410 SimCtx *simCtx = NULL;
411 UserCtx *user = NULL;
412 PetscReal ***phi = NULL;
413 Cmpnts ***ucont = NULL;
414
415 PetscFunctionBeginUser;
416 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
417 simCtx->dt = COEF_TIME_ACCURACY;
418 PetscCall(VecSet(user->Ucont, 0.0));
419 PetscCall(VecSet(user->Nvert, 0.0));
420 PetscCall(DMDAVecGetArray(user->da, user->Phi, &phi));
421 for (PetscInt k = user->info.zs; k < user->info.zs + user->info.zm; ++k) {
422 for (PetscInt j = user->info.ys; j < user->info.ys + user->info.ym; ++j) {
423 for (PetscInt i = user->info.xs; i < user->info.xs + user->info.xm; ++i) {
424 phi[k][j][i] = (PetscReal)i;
425 }
426 }
427 }
428 PetscCall(DMDAVecRestoreArray(user->da, user->Phi, &phi));
429 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
430 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
431 PetscCall(DMGlobalToLocalBegin(user->da, user->Phi, INSERT_VALUES, user->lPhi));
432 PetscCall(DMGlobalToLocalEnd(user->da, user->Phi, INSERT_VALUES, user->lPhi));
433 PetscCall(DMGlobalToLocalBegin(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
434 PetscCall(DMGlobalToLocalEnd(user->da, user->Nvert, INSERT_VALUES, user->lNvert));
435
436 PetscCall(Projection(user));
437
438 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
439 PetscCall(PicurvAssertRealNear(-1.0, ucont[2][2][2].x, 1.0e-10, "Projection should subtract the x pressure gradient under identity metrics"));
440 PetscCall(PicurvAssertRealNear(0.0, ucont[2][2][2].y, 1.0e-10, "Projection should leave the y component unchanged for an x-only gradient"));
441 PetscCall(PicurvAssertRealNear(0.0, ucont[2][2][2].z, 1.0e-10, "Projection should leave the z component unchanged for an x-only gradient"));
442 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
443
444 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
445 PetscFunctionReturn(0);
446}
447/**
448 * @brief Runs the unit-poisson-rhs PETSc test binary.
449 */
450
451int main(int argc, char **argv)
452{
453 PetscErrorCode ierr;
454 const PicurvTestCase cases[] = {
455 {"update-pressure-adds-phi", TestUpdatePressureAddsPhi},
456 {"poisson-rhs-zero-divergence", TestPoissonRHSZeroDivergence},
457 {"compute-body-forces-dispatcher", TestComputeBodyForcesDispatcher},
458 {"compute-eulerian-diffusivity-molecular-only", TestComputeEulerianDiffusivityMolecularOnly},
459 {"convection-zero-field", TestConvectionZeroField},
460 {"viscous-uniform-field", TestViscousUniformField},
461 {"compute-rhs-zero-field-no-forcing", TestComputeRHSZeroFieldNoForcing},
462 {"compute-eulerian-diffusivity-gradient-constant-field", TestComputeEulerianDiffusivityGradientConstantField},
463 {"compute-eulerian-diffusivity-verification-linear-x", TestComputeEulerianDiffusivityVerificationLinearX},
464 {"poisson-null-space-function-removes-mean", TestPoissonNullSpaceFunctionRemovesMean},
465 {"poisson-lhs-new-assembles-operator", TestPoissonLHSNewAssemblesOperator},
466 {"projection-zero-phi-leaves-velocity-unchanged", TestProjectionZeroPhiLeavesVelocityUnchanged},
467 {"projection-linear-phi-corrects-velocity", TestProjectionLinearPhiCorrectsVelocity},
468 };
469
470 ierr = PetscInitialize(&argc, &argv, NULL, "PICurv Poisson/RHS tests");
471 if (ierr) {
472 return (int)ierr;
473 }
474
475 ierr = PicurvRunTests("unit-poisson-rhs", cases, sizeof(cases) / sizeof(cases[0]));
476 if (ierr) {
477 PetscFinalize();
478 return (int)ierr;
479 }
480
481 ierr = PetscFinalize();
482 return (int)ierr;
483}
PoissonNullSpaceFunction
PetscErrorCode PoissonNullSpaceFunction(MatNullSpace nullsp, Vec X, void *ctx)
The callback function for PETSc's MatNullSpace object.
Definition poisson.c:1031
PoissonLHSNew
PetscErrorCode PoissonLHSNew(UserCtx *user)
Assembles the Left-Hand-Side (LHS) matrix (Laplacian operator) for the Poisson equation on a single g...
Definition poisson.c:1532
Projection
PetscErrorCode Projection(UserCtx *user)
Corrects the contravariant velocity field Ucont to be divergence-free using the gradient of the press...
Definition poisson.c:328
UpdatePressure
PetscErrorCode UpdatePressure(UserCtx *user)
Updates the pressure field P with the pressure correction Phi computed by the Poisson solver.
Definition poisson.c:854
PoissonRHS
PetscErrorCode PoissonRHS(UserCtx *user, Vec B)
Computes the Right-Hand-Side (RHS) of the Poisson equation, which is the divergence of the intermedia...
Definition poisson.c:2141
ComputeBodyForces
PetscErrorCode ComputeBodyForces(UserCtx *user, Vec Rct)
General dispatcher for applying all active body forces (momentum sources).
Definition rhs.c:1162
Viscous
PetscErrorCode Viscous(UserCtx *user, Vec Ucont, Vec Ucat, Vec Visc)
Computes the viscous contribution to the contravariant momentum RHS.
Definition rhs.c:519
ComputeEulerianDiffusivity
PetscErrorCode ComputeEulerianDiffusivity(UserCtx *user)
Computes the effective diffusivity scalar field (Gamma_eff) on the Eulerian grid.
Definition rhs.c:1959
ComputeEulerianDiffusivityGradient
PetscErrorCode ComputeEulerianDiffusivityGradient(UserCtx *user)
Computes the Eulerian gradient of the effective diffusivity field.
Definition rhs.c:2089
Convection
PetscErrorCode Convection(UserCtx *user, Vec Ucont, Vec Ucat, Vec Conv)
Computes the convective contribution to the contravariant momentum RHS.
Definition rhs.c:13
ComputeRHS
PetscErrorCode ComputeRHS(UserCtx *user, Vec Rhs)
Computes the Right-Hand Side (RHS) of the momentum equations.
Definition rhs.c:1185
TestComputeEulerianDiffusivityMolecularOnly
static PetscErrorCode TestComputeEulerianDiffusivityMolecularOnly(void)
Tests Eulerian diffusivity for the molecular-only configuration.
Definition test_poisson_rhs.c:139
TestComputeRHSZeroFieldNoForcing
static PetscErrorCode TestComputeRHSZeroFieldNoForcing(void)
Tests that the full RHS remains zero without forcing on a quiescent field.
Definition test_poisson_rhs.c:221
TestViscousUniformField
static PetscErrorCode TestViscousUniformField(void)
Tests that viscous terms vanish for a uniform field.
Definition test_poisson_rhs.c:194
main
int main(int argc, char **argv)
Runs the unit-poisson-rhs PETSc test binary.
Definition test_poisson_rhs.c:451
TestConvectionZeroField
static PetscErrorCode TestConvectionZeroField(void)
Tests that convection vanishes for a quiescent field.
Definition test_poisson_rhs.c:167
TestComputeBodyForcesDispatcher
static PetscErrorCode TestComputeBodyForcesDispatcher(void)
Tests the body-force dispatcher across supported source terms.
Definition test_poisson_rhs.c:85
TestPoissonLHSNewAssemblesOperator
static PetscErrorCode TestPoissonLHSNewAssemblesOperator(void)
Tests that Poisson matrix assembly produces a populated operator on a tiny Cartesian grid.
Definition test_poisson_rhs.c:346
TestUpdatePressureAddsPhi
static PetscErrorCode TestUpdatePressureAddsPhi(void)
Tests that pressure updates add the correction potential.
Definition test_poisson_rhs.c:30
TestProjectionZeroPhiLeavesVelocityUnchanged
static PetscErrorCode TestProjectionZeroPhiLeavesVelocityUnchanged(void)
Tests that projection leaves a zero pressure-correction field unchanged.
Definition test_poisson_rhs.c:382
TestPoissonRHSZeroDivergence
static PetscErrorCode TestPoissonRHSZeroDivergence(void)
Tests that the Poisson RHS is zero for zero-divergence velocity.
Definition test_poisson_rhs.c:51
TestPoissonNullSpaceFunctionRemovesMean
static PetscErrorCode TestPoissonNullSpaceFunctionRemovesMean(void)
Tests that the Poisson null-space operator removes the mean from a constant field.
Definition test_poisson_rhs.c:316
TestComputeEulerianDiffusivityVerificationLinearX
static PetscErrorCode TestComputeEulerianDiffusivityVerificationLinearX(void)
Tests verification-driven linear diffusivity override and its gradient.
Definition test_poisson_rhs.c:269
EnsurePoissonAndRhsVectors
static PetscErrorCode EnsurePoissonAndRhsVectors(UserCtx *user)
Allocates Poisson/RHS support vectors required by the tests.
Definition test_poisson_rhs.c:15
TestProjectionLinearPhiCorrectsVelocity
static PetscErrorCode TestProjectionLinearPhiCorrectsVelocity(void)
Tests that projection applies the expected x-direction correction for a linear pressure field.
Definition test_poisson_rhs.c:408
TestComputeEulerianDiffusivityGradientConstantField
static PetscErrorCode TestComputeEulerianDiffusivityGradientConstantField(void)
Tests diffusivity-gradient computation on a constant field.
Definition test_poisson_rhs.c:248
PicurvCreateMinimalContexts
PetscErrorCode PicurvCreateMinimalContexts(SimCtx **simCtx_out, UserCtx **user_out, PetscInt mx, PetscInt my, PetscInt mz)
Builds minimal SimCtx and UserCtx fixtures for C unit tests.
Definition test_support.c:530
PicurvAssertRealNear
PetscErrorCode PicurvAssertRealNear(PetscReal expected, PetscReal actual, PetscReal tol, const char *context)
Asserts that two real values agree within tolerance.
Definition test_support.c:902
PicurvDestroyMinimalContexts
PetscErrorCode PicurvDestroyMinimalContexts(SimCtx **simCtx_ptr, UserCtx **user_ptr)
Destroys minimal SimCtx/UserCtx fixtures and all owned PETSc objects.
Definition test_support.c:783
PicurvRunTests
PetscErrorCode PicurvRunTests(const char *suite_name, const PicurvTestCase *cases, size_t case_count)
Runs a named C test suite and prints pass/fail progress markers.
Definition test_support.c:138
PicurvAssertVecConstant
PetscErrorCode PicurvAssertVecConstant(Vec vec, PetscScalar expected, PetscReal tol, const char *context)
Asserts that a PETSc vector is spatially constant within tolerance.
Definition test_support.c:961
PicurvAssertIntEqual
PetscErrorCode PicurvAssertIntEqual(PetscInt expected, PetscInt actual, const char *context)
Asserts that two integer values are equal.
Definition test_support.c:917
PicurvAssertBool
PetscErrorCode PicurvAssertBool(PetscBool value, const char *context)
Asserts that one boolean condition is true.
Definition test_support.c:932
test_support.h
Shared declarations for the PICurv C test fixture and assertion layer.
PicurvTestCase
Named test case descriptor consumed by PicurvRunTests.
Definition test_support.h:17
PERIODIC
@ PERIODIC
Definition variables.h:260
UserCtx::Rhs
Vec Rhs
Definition variables.h:845
SimCtx::schmidt_number
PetscReal schmidt_number
Definition variables.h:709
UserCtx::boundary_faces
BoundaryFaceConfig boundary_faces[6]
Definition variables.h:829
UserCtx::lNvert
Vec lNvert
Definition variables.h:837
UserCtx::Phi
Vec Phi
Definition variables.h:837
SimCtx::forceScalingFactor
PetscReal forceScalingFactor
Definition variables.h:723
SimCtx::rans
PetscInt rans
Definition variables.h:732
UserCtx::lZet
Vec lZet
Definition variables.h:858
BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX
@ BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX
Definition variables.h:286
SimCtx::ren
PetscReal ren
Definition variables.h:691
BoundaryFaceConfig::handler_type
BCHandlerType handler_type
Definition variables.h:337
SimCtx::dt
PetscReal dt
Definition variables.h:658
SimCtx::bulkVelocityCorrection
PetscReal bulkVelocityCorrection
Definition variables.h:725
UserCtx::DiffusivityGradient
Vec DiffusivityGradient
Definition variables.h:841
UserCtx::Ucont
Vec Ucont
Definition variables.h:837
Cmpnts::x
PetscScalar x
Definition variables.h:101
UserCtx::lPhi
Vec lPhi
Definition variables.h:837
UserCtx::lCsi
Vec lCsi
Definition variables.h:858
SimCtx::eulerianSource
char eulerianSource[PETSC_MAX_PATH_LEN]
Definition variables.h:663
Cmpnts::z
PetscScalar z
Definition variables.h:101
UserCtx::Ucat
Vec Ucat
Definition variables.h:837
UserCtx::lUcont
Vec lUcont
Definition variables.h:837
UserCtx::Diffusivity
Vec Diffusivity
Definition variables.h:840
UserCtx::lAj
Vec lAj
Definition variables.h:858
UserCtx::info
DMDALocalInfo info
Definition variables.h:818
UserCtx::lUcat
Vec lUcat
Definition variables.h:837
Cmpnts::y
PetscScalar y
Definition variables.h:101
UserCtx::IM
PetscInt IM
Definition variables.h:820
UserCtx::lEta
Vec lEta
Definition variables.h:858
UserCtx::Cent
Vec Cent
Definition variables.h:858
SimCtx::les
PetscInt les
Definition variables.h:732
UserCtx::Nvert
Vec Nvert
Definition variables.h:837
UserCtx::lDiffusivity
Vec lDiffusivity
Definition variables.h:840
BoundaryFaceConfig::mathematical_type
BCType mathematical_type
Definition variables.h:336
SimCtx::summationRHS
PetscReal summationRHS
Definition variables.h:770
COEF_TIME_ACCURACY
#define COEF_TIME_ACCURACY
Coefficient controlling the temporal accuracy scheme (e.g., 1.5 for 2nd Order Backward Difference).
Definition variables.h:57
SimCtx::verificationDiffusivity
VerificationDiffusivityConfig verificationDiffusivity
Definition variables.h:699
SimCtx::drivingForceMagnitude
PetscReal drivingForceMagnitude
Definition variables.h:723
BC_FACE_NEG_X
@ BC_FACE_NEG_X
Definition variables.h:245
BC_FACE_POS_X
@ BC_FACE_POS_X
Definition variables.h:245
Cmpnts
A 3D point or vector with PetscScalar components.
Definition variables.h:100
SimCtx
The master context for the entire simulation.
Definition variables.h:643
UserCtx
User-defined context containing data specific to a single computational grid level.
Definition variables.h:811
verification_sources.h
VerificationDiffusivityOverrideActive
PetscBool VerificationDiffusivityOverrideActive(const SimCtx *simCtx)
Reports whether a verification-only diffusivity override is active.
Definition verification_sources.c:11
Generated by doxygen 1.9.8