PICurv 0.1.0
A Parallel Particle-In-Cell Solver for Curvilinear LES
 
Loading...
Searching...
No Matches
test_runtime_kernels.c
Go to the documentation of this file.
1/**
2 * @file test_runtime_kernels.c
3 * @brief C unit tests for runtime, particle, wall, and walltime-guard helpers.
4 */
5
6#include "test_support.h"
7
8#include "BC_Handlers.h"
9#include "ParticleMotion.h"
10#include "ParticlePhysics.h"
11#include "ParticleSwarm.h"
12#include "interpolation.h"
13#include "initialcondition.h"
14#include "les.h"
15#include "runloop.h"
16#include "setup.h"
17#include "wallfunction.h"
18#include "walkingsearch.h"
19
20/**
21 * @brief Synchronizes the minimal runtime fixture's global fields into their persistent local ghosts.
22 */
23static PetscErrorCode SyncRuntimeFieldGhosts(UserCtx *user)
24{
25 PetscFunctionBeginUser;
26 PetscCall(DMGlobalToLocalBegin(user->da, user->P, INSERT_VALUES, user->lP));
27 PetscCall(DMGlobalToLocalEnd(user->da, user->P, INSERT_VALUES, user->lP));
28 PetscCall(DMGlobalToLocalBegin(user->da, user->Psi, INSERT_VALUES, user->lPsi));
29 PetscCall(DMGlobalToLocalEnd(user->da, user->Psi, INSERT_VALUES, user->lPsi));
30 PetscCall(DMGlobalToLocalBegin(user->da, user->Diffusivity, INSERT_VALUES, user->lDiffusivity));
31 PetscCall(DMGlobalToLocalEnd(user->da, user->Diffusivity, INSERT_VALUES, user->lDiffusivity));
32 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
33 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
34 PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
35 PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucont, INSERT_VALUES, user->lUcont));
36 PetscCall(DMGlobalToLocalBegin(user->fda, user->DiffusivityGradient, INSERT_VALUES, user->lDiffusivityGradient));
37 PetscCall(DMGlobalToLocalEnd(user->fda, user->DiffusivityGradient, INSERT_VALUES, user->lDiffusivityGradient));
38 PetscFunctionReturn(0);
39}
40
41/**
42 * @brief Seeds one localized swarm particle with the cell, position, weight, and status data used by runtime tests.
43 */
44static PetscErrorCode SeedSingleParticle(UserCtx *user,
45 PetscInt ci,
46 PetscInt cj,
47 PetscInt ck,
48 PetscReal x,
49 PetscReal y,
50 PetscReal z,
51 PetscReal wx,
52 PetscReal wy,
53 PetscReal wz,
54 PetscInt status_value)
55{
56 PetscReal *positions = NULL;
57 PetscReal *weights = NULL;
58 PetscInt *cell_ids = NULL;
59 PetscInt *status = NULL;
60
61 PetscFunctionBeginUser;
62 PetscCall(DMSwarmGetField(user->swarm, "position", NULL, NULL, (void **)&positions));
63 PetscCall(DMSwarmGetField(user->swarm, "weight", NULL, NULL, (void **)&weights));
64 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
65 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
66
67 positions[0] = x;
68 positions[1] = y;
69 positions[2] = z;
70 weights[0] = wx;
71 weights[1] = wy;
72 weights[2] = wz;
73 cell_ids[0] = ci;
74 cell_ids[1] = cj;
75 cell_ids[2] = ck;
76 status[0] = status_value;
77
78 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
79 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
80 PetscCall(DMSwarmRestoreField(user->swarm, "weight", NULL, NULL, (void **)&weights));
81 PetscCall(DMSwarmRestoreField(user->swarm, "position", NULL, NULL, (void **)&positions));
82 PetscFunctionReturn(0);
83}
84/**
85 * @brief Tests particle distribution remainder handling across ranks.
86 */
87
88static PetscErrorCode TestDistributeParticlesRemainderHandling(void)
89{
90 PetscInt particles_per_rank = -1;
91 PetscInt remainder = -1;
92
93 PetscFunctionBeginUser;
94 PetscCall(DistributeParticles(10, 0, 3, &particles_per_rank, &remainder));
95 PetscCall(PicurvAssertIntEqual(4, particles_per_rank, "rank 0 should receive one remainder particle"));
96 PetscCall(PicurvAssertIntEqual(1, remainder, "remainder should be reported correctly"));
97
98 PetscCall(DistributeParticles(10, 2, 3, &particles_per_rank, &remainder));
99 PetscCall(PicurvAssertIntEqual(3, particles_per_rank, "last rank should receive base particle count"));
100 PetscCall(PicurvAssertIntEqual(1, remainder, "remainder should remain unchanged across ranks"));
101 PetscFunctionReturn(0);
102}
103/**
104 * @brief Tests basic particle-inside-bounding-box classification cases.
105 */
106
107static PetscErrorCode TestIsParticleInsideBoundingBoxBasicCases(void)
108{
109 BoundingBox bbox;
110 Particle particle;
111
112 PetscFunctionBeginUser;
113 PetscCall(PetscMemzero(&bbox, sizeof(bbox)));
114 PetscCall(PetscMemzero(&particle, sizeof(particle)));
115
116 bbox.min_coords.x = 0.0;
117 bbox.min_coords.y = 0.0;
118 bbox.min_coords.z = 0.0;
119 bbox.max_coords.x = 1.0;
120 bbox.max_coords.y = 2.0;
121 bbox.max_coords.z = 3.0;
122
123 particle.loc.x = 0.25;
124 particle.loc.y = 1.0;
125 particle.loc.z = 2.5;
126 PetscCall(PicurvAssertBool(IsParticleInsideBoundingBox(&bbox, &particle), "particle should be inside bounding box"));
127
128 particle.loc.x = 1.5;
129 PetscCall(PicurvAssertBool((PetscBool)!IsParticleInsideBoundingBox(&bbox, &particle), "particle should be outside bounding box"));
130 PetscFunctionReturn(0);
131}
132/**
133 * @brief Tests particle weight updates against expected ratios.
134 */
135
136static PetscErrorCode TestUpdateParticleWeightsComputesExpectedRatios(void)
137{
138 Particle particle;
139 PetscReal distances[NUM_FACES] = {1.0, 3.0, 2.0, 2.0, 4.0, 1.0};
140 PetscReal clamped[NUM_FACES] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
141
142 PetscFunctionBeginUser;
143 PetscCall(PetscMemzero(&particle, sizeof(particle)));
144 PetscCall(UpdateParticleWeights(distances, &particle));
145 PetscCall(PicurvAssertRealNear(0.25, particle.weights.x, 1.0e-12, "x interpolation weight"));
146 PetscCall(PicurvAssertRealNear(0.5, particle.weights.y, 1.0e-12, "y interpolation weight"));
147 PetscCall(PicurvAssertRealNear(0.2, particle.weights.z, 1.0e-12, "z interpolation weight"));
148
149 PetscCall(UpdateParticleWeights(clamped, &particle));
150 PetscCall(PicurvAssertRealNear(0.5, particle.weights.x, 1.0e-12, "clamped x weight remains centered"));
151 PetscCall(PicurvAssertRealNear(0.5, particle.weights.y, 1.0e-12, "clamped y weight remains centered"));
152 PetscCall(PicurvAssertRealNear(0.5, particle.weights.z, 1.0e-12, "clamped z weight remains centered"));
153 PetscFunctionReturn(0);
154}
155/**
156 * @brief Tests particle position updates without Brownian forcing.
157 */
158
159static PetscErrorCode TestUpdateParticlePositionWithoutBrownianContribution(void)
160{
161 SimCtx *simCtx = NULL;
162 UserCtx *user = NULL;
163 Particle particle;
164
165 PetscFunctionBeginUser;
166 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
167 simCtx->dt = 0.25;
168
169 PetscCall(PetscMemzero(&particle, sizeof(particle)));
170 particle.loc.x = 1.0;
171 particle.loc.y = -2.0;
172 particle.loc.z = 3.0;
173 particle.vel.x = 0.5;
174 particle.vel.y = -1.0;
175 particle.vel.z = 2.0;
176 particle.diffusivitygradient.x = 0.1;
177 particle.diffusivitygradient.y = 0.2;
178 particle.diffusivitygradient.z = -0.3;
179 particle.diffusivity = 0.0;
180
181 PetscCall(UpdateParticlePosition(user, &particle));
182 PetscCall(PicurvAssertRealNear(1.15, particle.loc.x, 1.0e-12, "x position update"));
183 PetscCall(PicurvAssertRealNear(-2.2, particle.loc.y, 1.0e-12, "y position update"));
184 PetscCall(PicurvAssertRealNear(3.425, particle.loc.z, 1.0e-12, "z position update"));
185
186 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
187 PetscFunctionReturn(0);
188}
189/**
190 * @brief Tests particle position updates driven only by diffusivity-gradient drift.
191 */
192
193static PetscErrorCode TestUpdateParticlePositionDiffusivityGradientOnly(void)
194{
195 SimCtx *simCtx = NULL;
196 UserCtx *user = NULL;
197 Particle particle;
198
199 PetscFunctionBeginUser;
200 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
201 simCtx->dt = 0.5;
202
203 PetscCall(PetscMemzero(&particle, sizeof(particle)));
204 particle.loc.x = 0.25;
205 particle.loc.y = 0.5;
206 particle.loc.z = 0.75;
207 particle.vel.x = 0.0;
208 particle.vel.y = 0.0;
209 particle.vel.z = 0.0;
210 particle.diffusivitygradient.x = 0.2;
211 particle.diffusivitygradient.y = -0.1;
212 particle.diffusivitygradient.z = 0.3;
213 particle.diffusivity = 0.0;
214
215 PetscCall(UpdateParticlePosition(user, &particle));
216 PetscCall(PicurvAssertRealNear(0.35, particle.loc.x, 1.0e-12, "drift-only x position update"));
217 PetscCall(PicurvAssertRealNear(0.45, particle.loc.y, 1.0e-12, "drift-only y position update"));
218 PetscCall(PicurvAssertRealNear(0.90, particle.loc.z, 1.0e-12, "drift-only z position update"));
219
220 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
221 PetscFunctionReturn(0);
222}
223/**
224 * @brief Tests IEM relaxation updates for particle-carried fields.
225 */
226
227static PetscErrorCode TestUpdateParticleFieldIEMRelaxation(void)
228{
229 PetscReal psi = 1.0;
230 const PetscReal dt = 0.5;
231 const PetscReal diffusivity = 0.2;
232 const PetscReal mean_val = 3.0;
233 const PetscReal cell_vol = 8.0;
234 const PetscReal c_model = 2.0;
235 PetscReal unchanged = 7.0;
236
237 PetscFunctionBeginUser;
238 PetscCall(UpdateParticleField(PARTICLE_FIELD_ID_PSI, dt, &psi, diffusivity, mean_val, cell_vol, c_model));
239 PetscCall(PicurvAssertRealNear(
240 mean_val + (1.0 - mean_val) * PetscExpReal(-(c_model * diffusivity / PetscPowReal(cell_vol, 0.6666667)) * dt),
241 psi,
242 1.0e-12,
243 "IEM update should match analytical relaxation"));
244
245 PetscCall(UpdateParticleField(PARTICLE_FIELD_ID_VELOCITY, dt, &unchanged, diffusivity, mean_val, cell_vol, c_model));
246 PetscCall(PicurvAssertRealNear(7.0, unchanged, 1.0e-12, "unknown field should remain unchanged"));
247 PetscFunctionReturn(0);
248}
249/**
250 * @brief Tests that non-Ucont requests do not modify interior field initialization.
251 */
252
253static PetscErrorCode TestSetInitialInteriorFieldIgnoresNonUcontRequest(void)
254{
255 SimCtx *simCtx = NULL;
256 UserCtx *user = NULL;
257
258 PetscFunctionBeginUser;
259 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
260 PetscCall(VecSet(user->Ucont, 7.0));
261
262 PetscCall(SetInitialInteriorField(user, FIELD_ID_P));
263 PetscCall(PicurvAssertVecConstant(user->Ucont, 7.0, 1.0e-12, "non-Ucont request should not modify Ucont"));
264
265 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
266 PetscFunctionReturn(0);
267}
268/**
269 * @brief Tests cartesian Constant IC: Cart2Contra sets contravariant flux via metric dot product.
270 *
271 * The minimal runtime fixture uses identity face metrics. UniformCart2Contra
272 * therefore maps (u,v,w)=(0,0,2) directly to a Zeta flux of 2.
273 */
274static PetscErrorCode TestSetInitialInteriorFieldCartesianConstantSetsContravariantFlux(void)
275{
276 SimCtx *simCtx = NULL;
277 UserCtx *user = NULL;
278 Cmpnts ***ucont = NULL;
279
280 PetscFunctionBeginUser;
281 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
282 simCtx->initialConditionMode = IC_MODE_CONSTANT_CARTESIAN;
283 simCtx->InitialConstantContra.x = 0.0;
284 simCtx->InitialConstantContra.y = 0.0;
285 simCtx->InitialConstantContra.z = 2.0;
286 PetscCall(VecSet(user->Ucont, 0.0));
287
288 PetscCall(SetInitialInteriorField(user, FIELD_ID_UCONT));
289 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
290 PetscCall(PicurvAssertRealNear(0.0, ucont[1][1][1].x, 1.0e-10, "Xi flux is zero for pure-z velocity on Cartesian grid"));
291 PetscCall(PicurvAssertRealNear(0.0, ucont[1][1][1].y, 1.0e-10, "Eta flux is zero for pure-z velocity on Cartesian grid"));
292 PetscCall(PicurvAssertRealNear(2.0, ucont[1][1][1].z, 1.0e-10, "Zeta flux follows the identity fixture metric"));
293 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
294
295 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
296 PetscFunctionReturn(0);
297}
298
299/**
300 * @brief Tests curvilinear Constant IC: flow_direction selects the streamwise axis.
301 *
302 * icVelocityPhysical=2.0, FLOW_DIR_POS_ZETA. The minimal runtime fixture uses
303 * identity face metrics, so ucont.z = 2.0.
304 * Xi and Eta components are left at zero regardless of InitialConstantContra.
305 */
306static PetscErrorCode TestSetInitialInteriorFieldCurvilinearConstantViaFlowDirection(void)
307{
308 SimCtx *simCtx = NULL;
309 UserCtx *user = NULL;
310 Cmpnts ***ucont = NULL;
311
312 PetscFunctionBeginUser;
313 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
314 simCtx->initialConditionMode = IC_MODE_CONSTANT_STREAMWISE;
315 simCtx->icVelocityPhysical = 2.0;
316 simCtx->flowDirection = FLOW_DIR_POS_ZETA;
317 user->GridOrientation = 1;
318 /* InitialConstantContra intentionally non-zero to confirm only z-component is set */
319 simCtx->InitialConstantContra.x = 99.0;
320 simCtx->InitialConstantContra.y = 99.0;
321 simCtx->InitialConstantContra.z = 99.0;
322 PetscCall(VecSet(user->Ucont, 0.0));
323
324 PetscCall(SetInitialInteriorField(user, FIELD_ID_UCONT));
325 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
326 PetscCall(PicurvAssertRealNear(0.0, ucont[1][1][1].x, 1.0e-10, "Xi flux is zero in curvilinear Zeta mode"));
327 PetscCall(PicurvAssertRealNear(0.0, ucont[1][1][1].y, 1.0e-10, "Eta flux is zero in curvilinear Zeta mode"));
328 PetscCall(PicurvAssertRealNear(2.0, ucont[1][1][1].z, 1.0e-10, "Zeta flux follows the identity fixture metric"));
329 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
330
331 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
332 PetscFunctionReturn(0);
333}
334
335/**
336 * @brief Tests zero IC clears physical-cell contravariant velocity.
337 */
338static PetscErrorCode TestSetInitialInteriorFieldZeroClearsInterior(void)
339{
340 SimCtx *simCtx = NULL;
341 UserCtx *user = NULL;
342 Cmpnts ***ucont = NULL;
343
344 PetscFunctionBeginUser;
345 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
346 simCtx->initialConditionMode = IC_MODE_ZERO;
347 PetscCall(VecSet(user->Ucont, 7.0));
348
349 PetscCall(SetInitialInteriorField(user, FIELD_ID_UCONT));
350 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
351 PetscCall(PicurvAssertRealNear(0.0, ucont[2][2][2].x, 1.0e-12, "zero IC clears interior Xi flux"));
352 PetscCall(PicurvAssertRealNear(0.0, ucont[2][2][2].y, 1.0e-12, "zero IC clears interior Eta flux"));
353 PetscCall(PicurvAssertRealNear(0.0, ucont[2][2][2].z, 1.0e-12, "zero IC clears interior Zeta flux"));
354 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
355
356 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
357 PetscFunctionReturn(0);
358}
359
360/**
361 * @brief Tests Poiseuille IC is the product parabola that vanishes on the walls.
362 * @details Four cells per cross-stream axis: cell centres sit at logical 0.5..3.5 with the
363 * walls at 0 and 4, so the normalized coordinate is -0.75, -0.25, 0.25, 0.75.
364 */
365static PetscErrorCode TestSetInitialInteriorFieldPoiseuilleProfile(void)
366{
367 SimCtx *simCtx = NULL;
368 UserCtx *user = NULL;
369 Cmpnts ***ucont = NULL;
370
371 PetscFunctionBeginUser;
372 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
373 simCtx->initialConditionMode = IC_MODE_POISEUILLE;
374 simCtx->icVelocityPhysical = 3.0;
375 simCtx->flowDirection = FLOW_DIR_POS_ZETA;
376 user->GridOrientation = 1;
377 PetscCall(VecZeroEntries(user->Ucont));
378
379 PetscCall(SetInitialInteriorField(user, FIELD_ID_UCONT));
380 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
381 PetscCall(PicurvAssertRealNear(3.0 * 225.0 / 256.0, ucont[2][2][2].z, 1.0e-12,
382 "Poiseuille IC follows the parabola at a centre-adjacent cell"));
383 PetscCall(PicurvAssertRealNear(3.0 * 49.0 / 256.0, ucont[2][1][1].z, 1.0e-12,
384 "Poiseuille IC vanishes on the wall, not at the first cell centre"));
385 PetscCall(PicurvAssertRealNear(0.0, ucont[2][2][2].x, 1.0e-12, "Poiseuille IC leaves Xi flux zero"));
386 PetscCall(PicurvAssertRealNear(0.0, ucont[2][2][2].y, 1.0e-12, "Poiseuille IC leaves Eta flux zero"));
387 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
388
389 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
390 PetscFunctionReturn(0);
391}
392
393/**
394 * @brief Tests Poiseuille IC is uniform along a periodic cross-stream axis.
395 * @details A channel periodic in i has walls only in j, so its start is the exact
396 * one-dimensional parabola in j.
397 */
398static PetscErrorCode TestSetInitialInteriorFieldPoiseuilleUniformAlongPeriodicAxis(void)
399{
400 SimCtx *simCtx = NULL;
401 UserCtx *user = NULL;
402 Cmpnts ***ucont = NULL;
403
404 PetscFunctionBeginUser;
405 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
406 simCtx->initialConditionMode = IC_MODE_POISEUILLE;
407 simCtx->icVelocityPhysical = 3.0;
408 simCtx->flowDirection = FLOW_DIR_POS_ZETA;
409 simCtx->i_periodic = 1;
410 user->GridOrientation = 1;
411 PetscCall(VecZeroEntries(user->Ucont));
412
413 PetscCall(SetInitialInteriorField(user, FIELD_ID_UCONT));
414 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
415 PetscCall(PicurvAssertRealNear(3.0 * 15.0 / 16.0, ucont[2][2][1].z, 1.0e-12,
416 "periodic axis carries no parabola factor"));
417 PetscCall(PicurvAssertRealNear(ucont[2][2][1].z, ucont[2][2][3].z, 1.0e-12,
418 "profile is uniform along the periodic axis"));
419 PetscCall(PicurvAssertRealNear(3.0 * 7.0 / 16.0, ucont[2][1][2].z, 1.0e-12,
420 "wall-bounded axis keeps its parabola"));
421 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
422
423 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
424 PetscFunctionReturn(0);
425}
426
427/**
428 * @brief Tests spatially varying Cartesian-field conversion to contravariant fluxes.
429 */
430static PetscErrorCode TestCart2ContraConvertsCartesianField(void)
431{
432 SimCtx *simCtx = NULL;
433 UserCtx *user = NULL;
434 Cmpnts ***ucat = NULL;
435 Cmpnts ***ucont = NULL;
436
437 PetscFunctionBeginUser;
438 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
439 PetscCall(DMDAVecGetArray(user->fda, user->Ucat, &ucat));
440 for (PetscInt k = 0; k < 6; k++)
441 for (PetscInt j = 0; j < 6; j++)
442 for (PetscInt i = 0; i < 6; i++)
443 ucat[k][j][i] = (Cmpnts){(PetscReal)i, 2.0, 3.0};
444 PetscCall(DMDAVecRestoreArray(user->fda, user->Ucat, &ucat));
445 PetscCall(UpdateLocalGhosts(user, FIELD_ID_UCAT));
446 PetscCall(Cart2Contra(user));
447
448 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
449 PetscCall(PicurvAssertRealNear(1.5, ucont[1][1][1].x, 1.0e-10,
450 "Xi flux uses face-interpolated Ucat"));
451 PetscCall(PicurvAssertRealNear(2.0, ucont[1][1][1].y, 1.0e-10,
452 "Eta flux uses Cartesian y velocity"));
453 PetscCall(PicurvAssertRealNear(3.0, ucont[1][1][1].z, 1.0e-10,
454 "Zeta flux uses Cartesian z velocity"));
455 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
456
457 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
458 PetscFunctionReturn(0);
459}
460
461/** @brief Verifies periodic Ucat endpoints are repaired before Cart2Contra terminal faces. */
462static PetscErrorCode TestCart2ContraUsesFinalizedPeriodicUcat(void)
463{
464 SimCtx *simCtx=NULL; UserCtx *user=NULL; Cmpnts ***ucat=NULL,***ucont=NULL; PetscInt mx;
465 const FieldId cell_fields[] = {FIELD_ID_UCAT}, staggered_fields[] = {FIELD_ID_UCONT};
466 PetscFunctionBeginUser;
467 PetscCall(PicurvCreateMinimalContextsWithPeriodicity(&simCtx,&user,4,4,4,PETSC_TRUE,PETSC_FALSE,PETSC_FALSE));
468 user->boundary_faces[BC_FACE_NEG_X].mathematical_type=PERIODIC;
469 user->boundary_faces[BC_FACE_POS_X].mathematical_type=PERIODIC;
470 mx=user->info.mx;
471 PetscCall(VecSet(user->Ucat,0.0));
472 PetscCall(DMDAVecGetArray(user->fda,user->Ucat,&ucat));
473 for(PetscInt k=1;k<user->info.mz-1;k++) for(PetscInt j=1;j<user->info.my-1;j++)
474 for(PetscInt i=1;i<mx-1;i++) ucat[k][j][i]=(Cmpnts){10.0+(PetscReal)i,2.0,3.0};
475 ucat[2][2][0]=(Cmpnts){-1000.0,-1000.0,-1000.0};
476 ucat[2][2][mx-1]=(Cmpnts){-2000.0,-2000.0,-2000.0};
477 PetscCall(DMDAVecRestoreArray(user->fda,user->Ucat,&ucat));
478 PetscCall(SynchronizePeriodicCellFields(user,1,cell_fields));
479 PetscCall(UpdateLocalGhosts(user, FIELD_ID_UCAT));
480 PetscCall(Cart2Contra(user));
481 PetscCall(SynchronizePeriodicStaggeredFields(user,1,staggered_fields));
482 PetscCall(DMDAVecGetArrayRead(user->fda,user->Ucont,&ucont));
483 PetscCall(PicurvAssertRealNear(0.5*((10.0+mx-2)+(10.0+1)),ucont[2][2][mx-2].x,1e-12,"terminal Xi face uses corrected periodic Ucat endpoint"));
484 PetscCall(DMDAVecRestoreArrayRead(user->fda,user->Ucont,&ucont));
485 PetscCall(PicurvDestroyMinimalContexts(&simCtx,&user));
486 PetscFunctionReturn(0);
487}
488
489/**
490 * @brief Tests loading a staged Ucat file IC and converting it to Ucont.
491 */
492static PetscErrorCode TestPopulateInitialUcontLoadsStagedUcat(void)
493{
494 SimCtx *simCtx = NULL;
495 UserCtx *user = NULL;
496 Cmpnts ***ucont = NULL;
497 char tmpdir[PETSC_MAX_PATH_LEN];
498
499 PetscFunctionBeginUser;
500 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
501 PetscCall(PicurvMakeTempDir(tmpdir, sizeof(tmpdir)));
502 PetscCall(VecSet(user->Ucat, 2.5));
503 PetscCall(PetscStrncpy(simCtx->_io_context_buffer, tmpdir, sizeof(simCtx->_io_context_buffer)));
504 simCtx->current_io_directory = simCtx->_io_context_buffer;
505 PetscCall(WriteFieldData(user, "ufield00000_0", user->Ucat, "dat"));
506 simCtx->current_io_directory = NULL;
507
508 PetscCall(VecZeroEntries(user->Ucat));
509 PetscCall(VecZeroEntries(user->Ucont));
510 simCtx->initialConditionMode = IC_MODE_FILE;
511 simCtx->initialConditionField = IC_FIELD_UCAT;
512 PetscCall(PetscStrncpy(simCtx->initialConditionDirectory, tmpdir, sizeof(simCtx->initialConditionDirectory)));
513 PetscCall(PopulateInitialUcont(user));
514
515 PetscCall(PicurvAssertVecConstant(user->Ucat, 2.5, 1.0e-12, "file IC should restore Ucat"));
516 PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
517 PetscCall(PicurvAssertRealNear(2.5, ucont[1][1][1].x, 1.0e-12, "file Ucat IC should populate interior Xi flux"));
518 PetscCall(PicurvAssertRealNear(2.5, ucont[1][1][1].y, 1.0e-12, "file Ucat IC should populate interior Eta flux"));
519 PetscCall(PicurvAssertRealNear(2.5, ucont[1][1][1].z, 1.0e-12, "file Ucat IC should populate interior Zeta flux"));
520 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
521 PetscCall(PicurvRemoveTempDir(tmpdir));
522 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
523 PetscFunctionReturn(0);
524}
525
526/** @brief Tests LES output staging from local vectors on a multiply periodic DMDA. */
527static PetscErrorCode TestWriteLESFieldsCopiesOwnedPeriodicValues(void)
528{
529 SimCtx *simCtx = NULL;
530 UserCtx *user = NULL;
531 char tmpdir[PETSC_MAX_PATH_LEN];
532 char euler_dir[PETSC_MAX_PATH_LEN];
533 char block_dir[PETSC_MAX_PATH_LEN];
534
535 PetscFunctionBeginUser;
536 PetscCall(PicurvCreateMinimalContextsWithPeriodicity(
537 &simCtx, &user, 4, 4, 4, PETSC_TRUE, PETSC_TRUE, PETSC_TRUE));
538 PetscCall(DMCreateGlobalVector(user->da, &user->CS));
539 PetscCall(DMCreateLocalVector(user->da, &user->lCs));
540 PetscCall(DMCreateGlobalVector(user->da, &user->Nu_t));
541 PetscCall(DMCreateLocalVector(user->da, &user->lNu_t));
542 PetscCall(VecSet(user->CS, 0.0));
543 PetscCall(VecSet(user->lCs, 2.5));
544 PetscCall(VecSet(user->Nu_t, 0.0));
545 PetscCall(VecSet(user->lNu_t, 7.5));
546 simCtx->les = DYNAMIC_SMAGORINSKY;
547
548 PetscCall(PicurvMakeTempDir(tmpdir, sizeof(tmpdir)));
549 PetscCall(PetscSNPrintf(euler_dir, sizeof(euler_dir), "%s/eulerian", tmpdir));
550 PetscCall(PetscSNPrintf(block_dir, sizeof(block_dir), "%s/block_0000", euler_dir));
551 PetscCall(PicurvEnsureDir(euler_dir));
552 PetscCall(PicurvEnsureDir(block_dir));
553 PetscCall(WriteSimulationFields(user, tmpdir));
554
555 PetscCall(PicurvAssertVecConstant(user->CS, 2.5, 1.0e-12,
556 "checkpoint output should copy owned periodic lCs values"));
557 PetscCall(PicurvAssertVecConstant(user->Nu_t, 7.5, 1.0e-12,
558 "checkpoint output should copy owned periodic lNu_t values"));
559
560 PetscCall(PicurvRemoveTempDir(tmpdir));
561 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
562 PetscFunctionReturn(0);
563}
564
565/**
566 * @brief Stages a real Cartesian Fourier mode through the production file-IC path.
567 *
568 * The returned maximum is the native face-interpolation/backward-difference
569 * divergence. On the uniform periodic fixture its symbol is
570 * (sin(theta_x), sin(theta_y), sin(theta_z)).
571 */
572static PetscErrorCode StagedFourierModeMaxDivergence(UserCtx *user,
573 PetscInt mode_x,
574 PetscInt mode_y,
575 PetscInt mode_z,
576 Cmpnts amplitude,
577 const char *tmpdir,
578 PetscReal *max_divergence)
579{
580 SimCtx *simCtx = user->simCtx;
581 Cmpnts ***ucat = NULL, ***ucont = NULL;
582 PetscInt n_x = user->info.mx - 2, n_y = user->info.my - 2, n_z = user->info.mz - 2;
583 const FieldId staggered_fields[] = {FIELD_ID_UCONT};
584
585 PetscFunctionBeginUser;
586 PetscCall(VecZeroEntries(user->Ucat));
587 PetscCall(DMDAVecGetArray(user->fda, user->Ucat, &ucat));
588 for (PetscInt k = 0; k < user->info.mz; ++k) {
589 for (PetscInt j = 0; j < user->info.my; ++j) {
590 for (PetscInt i = 0; i < user->info.mx; ++i) {
591 PetscReal phase = 2.0 * PETSC_PI * (
592 (PetscReal)(mode_x * (i - 1)) / (PetscReal)n_x +
593 (PetscReal)(mode_y * (j - 1)) / (PetscReal)n_y +
594 (PetscReal)(mode_z * (k - 1)) / (PetscReal)n_z);
595 PetscReal wave = PetscCosReal(phase);
596 ucat[k][j][i] = (Cmpnts){amplitude.x * wave, amplitude.y * wave, amplitude.z * wave};
597 }
598 }
599 }
600 PetscCall(DMDAVecRestoreArray(user->fda, user->Ucat, &ucat));
601
602 PetscCall(PetscStrncpy(simCtx->_io_context_buffer, tmpdir, sizeof(simCtx->_io_context_buffer)));
603 simCtx->current_io_directory = simCtx->_io_context_buffer;
604 PetscCall(WriteFieldData(user, "ufield00000_0", user->Ucat, "dat"));
605 simCtx->current_io_directory = NULL;
606 PetscCall(VecZeroEntries(user->Ucat));
607 PetscCall(VecZeroEntries(user->Ucont));
608 PetscCall(PopulateInitialUcont(user));
609 PetscCall(SynchronizePeriodicStaggeredFields(user, 1, staggered_fields));
610 PetscCall(UpdateLocalGhosts(user, FIELD_ID_UCONT));
611
612 *max_divergence = 0.0;
613 PetscCall(DMDAVecGetArrayRead(user->fda, user->lUcont, &ucont));
614 for (PetscInt k = 1; k < user->info.mz - 1; ++k) {
615 for (PetscInt j = 1; j < user->info.my - 1; ++j) {
616 for (PetscInt i = 1; i < user->info.mx - 1; ++i) {
617 PetscReal divergence =
618 ucont[k][j][i].x - ucont[k][j][i - 1].x +
619 ucont[k][j][i].y - ucont[k][j - 1][i].y +
620 ucont[k][j][i].z - ucont[k - 1][j][i].z;
621 *max_divergence = PetscMax(*max_divergence, PetscAbsReal(divergence));
622 }
623 }
624 }
625 PetscCall(DMDAVecRestoreArrayRead(user->fda, user->lUcont, &ucont));
626 PetscFunctionReturn(0);
627}
628
629/** @brief Cross-checks selected exact Fourier modes against PICurv's native discrete symbol. */
630static PetscErrorCode TestStagedUcatFourierModesUsePicurvDiscreteSymbol(void)
631{
632 SimCtx *simCtx = NULL;
633 UserCtx *user = NULL;
634 char tmpdir[PETSC_MAX_PATH_LEN];
635 PetscReal max_divergence = 0.0, sx, sy, sz;
636 PetscInt n;
637
638 PetscFunctionBeginUser;
639 PetscCall(PicurvCreateMinimalContextsWithPeriodicity(
640 &simCtx, &user, 9, 9, 9, PETSC_TRUE, PETSC_TRUE, PETSC_TRUE));
641 for (PetscInt face = 0; face < NUM_FACES; ++face)
642 user->boundary_faces[face].mathematical_type = PERIODIC;
643 PetscCall(PicurvMakeTempDir(tmpdir, sizeof(tmpdir)));
644 simCtx->initialConditionMode = IC_MODE_FILE;
645 simCtx->initialConditionField = IC_FIELD_UCAT;
646 PetscCall(PetscStrncpy(simCtx->initialConditionDirectory, tmpdir, sizeof(simCtx->initialConditionDirectory)));
647
648 PetscCall(StagedFourierModeMaxDivergence(
649 user, 1, 0, 0, (Cmpnts){0.0, 1.0, 0.0}, tmpdir, &max_divergence));
650 PetscCall(PicurvAssertRealNear(0.0, max_divergence, 1.0e-12,
651 "axis-aligned transverse staged mode"));
652
653 n = user->info.mx - 2;
654 sx = PetscSinReal(2.0 * PETSC_PI / (PetscReal)n);
655 sy = PetscSinReal(4.0 * PETSC_PI / (PetscReal)n);
656 sz = PetscSinReal(6.0 * PETSC_PI / (PetscReal)n);
657 PetscCall(StagedFourierModeMaxDivergence(
658 user, 1, 2, 0, (Cmpnts){sy, -sx, 0.0}, tmpdir, &max_divergence));
659 PetscCall(PicurvAssertRealNear(0.0, max_divergence, 1.0e-12,
660 "oblique first discrete transverse polarization"));
661 PetscCall(StagedFourierModeMaxDivergence(
662 user, 1, 2, 3,
663 (Cmpnts){sz * sx, sz * sy, -(sx * sx + sy * sy)},
664 tmpdir, &max_divergence));
665 PetscCall(PicurvAssertRealNear(0.0, max_divergence, 1.0e-12,
666 "oblique second discrete transverse polarization"));
667
668 PetscCall(PicurvRemoveTempDir(tmpdir));
669 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
670 PetscFunctionReturn(0);
671}
672
673/**
674 * @brief Tests loading a staged Ucont file IC without Cartesian conversion.
675 */
676static PetscErrorCode TestPopulateInitialUcontLoadsStagedUcont(void)
677{
678 SimCtx *simCtx = NULL;
679 UserCtx *user = NULL;
680 char tmpdir[PETSC_MAX_PATH_LEN];
681
682 PetscFunctionBeginUser;
683 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
684 PetscCall(PicurvMakeTempDir(tmpdir, sizeof(tmpdir)));
685 PetscCall(VecSet(user->Ucont, 4.5));
686 PetscCall(PetscStrncpy(simCtx->_io_context_buffer, tmpdir, sizeof(simCtx->_io_context_buffer)));
687 simCtx->current_io_directory = simCtx->_io_context_buffer;
688 PetscCall(WriteFieldData(user, "vfield00000_0", user->Ucont, "dat"));
689 simCtx->current_io_directory = NULL;
690
691 PetscCall(VecZeroEntries(user->Ucont));
692 simCtx->initialConditionMode = IC_MODE_FILE;
693 simCtx->initialConditionField = IC_FIELD_UCONT;
694 PetscCall(PetscStrncpy(simCtx->initialConditionDirectory, tmpdir, sizeof(simCtx->initialConditionDirectory)));
695 PetscCall(PopulateInitialUcont(user));
696
697 PetscCall(PicurvAssertVecConstant(user->Ucont, 4.5, 1.0e-12, "file Ucont IC should restore Ucont directly"));
698 PetscCall(PicurvRemoveTempDir(tmpdir));
699 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
700 PetscFunctionReturn(0);
701}
702/**
703 * @brief Tests direct interpolation from Eulerian fields to one localized swarm particle.
704 */
705
706static PetscErrorCode TestInterpolateAllFieldsToSwarmConstantFields(void)
707{
708 SimCtx *simCtx = NULL;
709 UserCtx *user = NULL;
710 Cmpnts ***grad = NULL;
711 PetscReal *velocity = NULL;
712 PetscReal *diffusivity = NULL;
713 PetscReal *diffusivity_gradient = NULL;
714
715 PetscFunctionBeginUser;
716 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
717 PetscCall(PicurvCreateSwarmPair(user, 1, "ske"));
718 PetscCall(VecSet(user->Ucat, 2.0));
719 PetscCall(VecSet(user->Diffusivity, 0.25));
720
721 PetscCall(DMDAVecGetArray(user->fda, user->DiffusivityGradient, &grad));
722 for (PetscInt k = user->info.zs; k < user->info.zs + user->info.zm; ++k) {
723 for (PetscInt j = user->info.ys; j < user->info.ys + user->info.ym; ++j) {
724 for (PetscInt i = user->info.xs; i < user->info.xs + user->info.xm; ++i) {
725 grad[k][j][i].x = 0.1;
726 grad[k][j][i].y = 0.2;
727 grad[k][j][i].z = 0.3;
728 }
729 }
730 }
731 PetscCall(DMDAVecRestoreArray(user->fda, user->DiffusivityGradient, &grad));
732 PetscCall(SyncRuntimeFieldGhosts(user));
733 PetscCall(SeedSingleParticle(user, 0, 0, 0, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, ACTIVE_AND_LOCATED));
734
735 PetscCall(InterpolateAllFieldsToSwarm(user));
736
737 PetscCall(DMSwarmGetField(user->swarm, "velocity", NULL, NULL, (void **)&velocity));
738 PetscCall(DMSwarmGetField(user->swarm, "Diffusivity", NULL, NULL, (void **)&diffusivity));
739 PetscCall(DMSwarmGetField(user->swarm, "DiffusivityGradient", NULL, NULL, (void **)&diffusivity_gradient));
740 PetscCall(PicurvAssertRealNear(2.0, velocity[0], 1.0e-12, "Interpolated velocity x should match constant Eulerian field"));
741 PetscCall(PicurvAssertRealNear(2.0, velocity[1], 1.0e-12, "Interpolated velocity y should match constant Eulerian field"));
742 PetscCall(PicurvAssertRealNear(2.0, velocity[2], 1.0e-12, "Interpolated velocity z should match constant Eulerian field"));
743 PetscCall(PicurvAssertRealNear(0.25, diffusivity[0], 1.0e-12, "Interpolated scalar diffusivity should match constant Eulerian field"));
744 PetscCall(PicurvAssertRealNear(0.1, diffusivity_gradient[0], 1.0e-12, "Interpolated diffusivity-gradient x component"));
745 PetscCall(PicurvAssertRealNear(0.2, diffusivity_gradient[1], 1.0e-12, "Interpolated diffusivity-gradient y component"));
746 PetscCall(PicurvAssertRealNear(0.3, diffusivity_gradient[2], 1.0e-12, "Interpolated diffusivity-gradient z component"));
747 PetscCall(DMSwarmRestoreField(user->swarm, "DiffusivityGradient", NULL, NULL, (void **)&diffusivity_gradient));
748 PetscCall(DMSwarmRestoreField(user->swarm, "Diffusivity", NULL, NULL, (void **)&diffusivity));
749 PetscCall(DMSwarmRestoreField(user->swarm, "velocity", NULL, NULL, (void **)&velocity));
750
751 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
752 PetscFunctionReturn(0);
753}
754/**
755 * @brief Tests the corner-averaged (legacy) interpolation path on constant fields.
756 */
757
758static PetscErrorCode TestInterpolateAllFieldsToSwarmCornerAveragedConstantFields(void)
759{
760 SimCtx *simCtx = NULL;
761 UserCtx *user = NULL;
762 Cmpnts ***grad = NULL;
763 PetscReal *velocity = NULL;
764 PetscReal *diffusivity = NULL;
765 PetscReal *diffusivity_gradient = NULL;
766
767 PetscFunctionBeginUser;
768 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
769 simCtx->interpolationMethod = INTERP_CORNER_AVERAGED;
770 PetscCall(PicurvCreateSwarmPair(user, 1, "ske"));
771 PetscCall(VecSet(user->Ucat, 2.0));
772 PetscCall(VecSet(user->Diffusivity, 0.25));
773
774 PetscCall(DMDAVecGetArray(user->fda, user->DiffusivityGradient, &grad));
775 for (PetscInt k = user->info.zs; k < user->info.zs + user->info.zm; ++k) {
776 for (PetscInt j = user->info.ys; j < user->info.ys + user->info.ym; ++j) {
777 for (PetscInt i = user->info.xs; i < user->info.xs + user->info.xm; ++i) {
778 grad[k][j][i].x = 0.1;
779 grad[k][j][i].y = 0.2;
780 grad[k][j][i].z = 0.3;
781 }
782 }
783 }
784 PetscCall(DMDAVecRestoreArray(user->fda, user->DiffusivityGradient, &grad));
785 PetscCall(SyncRuntimeFieldGhosts(user));
786 PetscCall(SeedSingleParticle(user, 0, 0, 0, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, ACTIVE_AND_LOCATED));
787
788 PetscCall(InterpolateAllFieldsToSwarm(user));
789
790 PetscCall(DMSwarmGetField(user->swarm, "velocity", NULL, NULL, (void **)&velocity));
791 PetscCall(DMSwarmGetField(user->swarm, "Diffusivity", NULL, NULL, (void **)&diffusivity));
792 PetscCall(DMSwarmGetField(user->swarm, "DiffusivityGradient", NULL, NULL, (void **)&diffusivity_gradient));
793 PetscCall(PicurvAssertRealNear(2.0, velocity[0], 1.0e-12, "CornerAveraged: interpolated velocity x should match constant Eulerian field"));
794 PetscCall(PicurvAssertRealNear(2.0, velocity[1], 1.0e-12, "CornerAveraged: interpolated velocity y should match constant Eulerian field"));
795 PetscCall(PicurvAssertRealNear(2.0, velocity[2], 1.0e-12, "CornerAveraged: interpolated velocity z should match constant Eulerian field"));
796 PetscCall(PicurvAssertRealNear(0.25, diffusivity[0], 1.0e-12, "CornerAveraged: interpolated scalar diffusivity should match constant Eulerian field"));
797 PetscCall(PicurvAssertRealNear(0.1, diffusivity_gradient[0], 1.0e-12, "CornerAveraged: interpolated diffusivity-gradient x component"));
798 PetscCall(PicurvAssertRealNear(0.2, diffusivity_gradient[1], 1.0e-12, "CornerAveraged: interpolated diffusivity-gradient y component"));
799 PetscCall(PicurvAssertRealNear(0.3, diffusivity_gradient[2], 1.0e-12, "CornerAveraged: interpolated diffusivity-gradient z component"));
800 PetscCall(DMSwarmRestoreField(user->swarm, "DiffusivityGradient", NULL, NULL, (void **)&diffusivity_gradient));
801 PetscCall(DMSwarmRestoreField(user->swarm, "Diffusivity", NULL, NULL, (void **)&diffusivity));
802 PetscCall(DMSwarmRestoreField(user->swarm, "velocity", NULL, NULL, (void **)&velocity));
803
804 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
805 PetscFunctionReturn(0);
806}
807/**
808 * @brief Tests particle-to-grid scattering using known cell occupancy and scalar values.
809 */
810
811static PetscErrorCode TestScatterAllParticleFieldsToEulerFieldsAveragesPsi(void)
812{
813 SimCtx *simCtx = NULL;
814 UserCtx *user = NULL;
815 PetscInt *cell_ids = NULL;
816 PetscReal *psi = NULL;
817 PetscReal ***psi_grid = NULL;
818
819 PetscFunctionBeginUser;
820 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
821 PetscCall(PicurvCreateSwarmPair(user, 2, "ske"));
822
823 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
824 cell_ids[0] = 0; cell_ids[1] = 0; cell_ids[2] = 0;
825 cell_ids[3] = 0; cell_ids[4] = 0; cell_ids[5] = 0;
826 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
827
828 PetscCall(DMSwarmGetField(user->swarm, "Psi", NULL, NULL, (void **)&psi));
829 psi[0] = 1.0;
830 psi[1] = 3.0;
831 PetscCall(DMSwarmRestoreField(user->swarm, "Psi", NULL, NULL, (void **)&psi));
832
833 PetscCall(ScatterAllParticleFieldsToEulerFields(user));
834
835 PetscCall(DMDAVecGetArrayRead(user->da, user->Psi, &psi_grid));
836 PetscCall(PicurvAssertRealNear(2.0, psi_grid[1][1][1], 1.0e-12, "Scatter should average particle Psi values into the owning cell"));
837 PetscCall(DMDAVecRestoreArrayRead(user->da, user->Psi, &psi_grid));
838
839 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
840 PetscFunctionReturn(0);
841}
842/**
843 * @brief Tests particle counting by geometric cell IDs using the production +1 storage shift.
844 */
845
846static PetscErrorCode TestCalculateParticleCountPerCellCountsGlobalCellIDs(void)
847{
848 SimCtx *simCtx = NULL;
849 UserCtx *user = NULL;
850 PetscInt *cell_ids = NULL;
851 PetscReal ***counts = NULL;
852
853 PetscFunctionBeginUser;
854 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
855 PetscCall(PicurvCreateSwarmPair(user, 3, "ske"));
856
857 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
858 cell_ids[0] = 0; cell_ids[1] = 0; cell_ids[2] = 0;
859 cell_ids[3] = 0; cell_ids[4] = 0; cell_ids[5] = 0;
860 cell_ids[6] = 1; cell_ids[7] = 0; cell_ids[8] = 0;
861 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
862
863 PetscCall(CalculateParticleCountPerCell(user));
864
865 PetscCall(DMDAVecGetArrayRead(user->da, user->ParticleCount, &counts));
866 PetscCall(PicurvAssertRealNear(2.0, counts[1][1][1], 1.0e-12, "Two particles in cell (0,0,0) should accumulate at shifted index (1,1,1)"));
867 PetscCall(PicurvAssertRealNear(1.0, counts[1][1][2], 1.0e-12, "One particle in cell (1,0,0) should accumulate at shifted index (2,1,1)"));
868 PetscCall(DMDAVecRestoreArrayRead(user->da, user->ParticleCount, &counts));
869
870 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
871 PetscFunctionReturn(0);
872}
873/**
874 * @brief Tests localized particle-status reset behavior for restart of the location workflow.
875 */
876
877static PetscErrorCode TestResetAllParticleStatusesLeavesLostParticlesUntouched(void)
878{
879 SimCtx *simCtx = NULL;
880 UserCtx *user = NULL;
881 PetscInt *status = NULL;
882
883 PetscFunctionBeginUser;
884 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
885 PetscCall(PicurvCreateSwarmPair(user, 3, "ske"));
886
887 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
888 status[0] = ACTIVE_AND_LOCATED;
889 status[1] = LOST;
890 status[2] = NEEDS_LOCATION;
891 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
892
893 PetscCall(ResetAllParticleStatuses(user));
894
895 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
896 PetscCall(PicurvAssertIntEqual(NEEDS_LOCATION, status[0], "ACTIVE_AND_LOCATED particles should be reset to NEEDS_LOCATION"));
897 PetscCall(PicurvAssertIntEqual(LOST, status[1], "LOST particles should remain LOST"));
898 PetscCall(PicurvAssertIntEqual(NEEDS_LOCATION, status[2], "NEEDS_LOCATION particles should remain unchanged"));
899 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
900
901 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
902 PetscFunctionReturn(0);
903}
904/**
905 * @brief Tests direct removal of particles that leave every rank bounding box.
906 */
907
908static PetscErrorCode TestCheckAndRemoveOutOfBoundsParticlesRemovesEscapedParticle(void)
909{
910 SimCtx *simCtx = NULL;
911 UserCtx *user = NULL;
912 PetscReal *positions = NULL;
913 PetscInt removed_local = 0;
914 PetscInt removed_global = 0;
915 PetscInt nlocal = 0;
916
917 PetscFunctionBeginUser;
918 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
919 PetscCall(PicurvCreateSwarmPair(user, 2, "ske"));
920
921 PetscCall(DMSwarmGetField(user->swarm, "position", NULL, NULL, (void **)&positions));
922 positions[0] = 0.5; positions[1] = 0.5; positions[2] = 0.5;
923 positions[3] = 9.0; positions[4] = 9.0; positions[5] = 9.0;
924 PetscCall(DMSwarmRestoreField(user->swarm, "position", NULL, NULL, (void **)&positions));
925
926 PetscCall(CheckAndRemoveOutOfBoundsParticles(user, &removed_local, &removed_global, simCtx->bboxlist));
927 PetscCall(DMSwarmGetLocalSize(user->swarm, &nlocal));
928 PetscCall(PicurvAssertIntEqual(1, removed_local, "Exactly one particle should be removed as out-of-bounds on a single rank"));
929 PetscCall(PicurvAssertIntEqual(1, removed_global, "Global out-of-bounds removal count should match the local single-rank result"));
930 PetscCall(PicurvAssertIntEqual(1, nlocal, "One in-bounds particle should remain after out-of-bounds removal"));
931
932 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
933 PetscFunctionReturn(0);
934}
935/**
936 * @brief Tests direct removal of particles already marked LOST by the location workflow.
937 */
938
939static PetscErrorCode TestCheckAndRemoveLostParticlesRemovesLostEntries(void)
940{
941 SimCtx *simCtx = NULL;
942 UserCtx *user = NULL;
943 PetscInt *status = NULL;
944 PetscInt removed_local = 0;
945 PetscInt removed_global = 0;
946 PetscInt nlocal = 0;
947 PetscReal *psi = NULL;
948 PetscReal removed_psi = 0.0;
949
950 PetscFunctionBeginUser;
951 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
952 PetscCall(PicurvCreateSwarmPair(user, 3, "ske"));
953 PetscCall(DMSwarmGetField(user->swarm, "Psi", NULL, NULL, (void **)&psi));
954 psi[0] = 0.5;
955 psi[1] = 2.0;
956 psi[2] = 3.0;
957 PetscCall(DMSwarmRestoreField(user->swarm, "Psi", NULL, NULL, (void **)&psi));
958
959 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
960 status[0] = ACTIVE_AND_LOCATED;
961 status[1] = LOST;
962 status[2] = LOST;
963 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
964
965 PetscCall(CheckAndRemoveLostParticles(user, &removed_local, &removed_global, &removed_psi));
966 PetscCall(DMSwarmGetLocalSize(user->swarm, &nlocal));
967 PetscCall(PicurvAssertRealNear(5.0, removed_psi, 0.0, "the removed particles carried Psi 2 and 3"));
968 PetscCall(PicurvAssertIntEqual(2, removed_local, "Two LOST particles should be removed locally"));
969 PetscCall(PicurvAssertIntEqual(2, removed_global, "Global LOST-particle removal count should match the local single-rank result"));
970 PetscCall(PicurvAssertIntEqual(1, nlocal, "One non-LOST particle should remain"));
971
972 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
973 PetscFunctionReturn(0);
974}
975/**
976 * @brief Tests Brownian displacement generation against a duplicated seeded RNG stream.
977 */
978
979static PetscErrorCode TestCalculateBrownianDisplacementDeterministicSeed(void)
980{
981 SimCtx *simCtx = NULL;
982 UserCtx *user = NULL;
983 char tmpdir[PETSC_MAX_PATH_LEN];
984 Cmpnts first;
985 Cmpnts second;
986
987 PetscFunctionBeginUser;
988 PetscCall(PicurvBuildTinyRuntimeContext(NULL, PETSC_FALSE, &simCtx, &user, tmpdir, sizeof(tmpdir)));
989 simCtx->dt = 0.25;
990 PetscCall(PicurvAssertBool((PetscBool)(simCtx->BrownianMotionRNG != NULL),
991 "runtime setup path should initialize the Brownian RNG"));
992
993 PetscCall(PetscRandomSetSeed(simCtx->BrownianMotionRNG, 12345));
994 PetscCall(PetscRandomSeed(simCtx->BrownianMotionRNG));
995
996 PetscCall(CalculateBrownianDisplacement(user, 0.5, &first));
997 PetscCall(PetscRandomSetSeed(simCtx->BrownianMotionRNG, 12345));
998 PetscCall(PetscRandomSeed(simCtx->BrownianMotionRNG));
999 PetscCall(CalculateBrownianDisplacement(user, 0.5, &second));
1000
1001 PetscCall(PicurvAssertRealNear(first.x, second.x, 1.0e-12, "Resetting the Brownian RNG seed should reproduce the x displacement"));
1002 PetscCall(PicurvAssertRealNear(first.y, second.y, 1.0e-12, "Resetting the Brownian RNG seed should reproduce the y displacement"));
1003 PetscCall(PicurvAssertRealNear(first.z, second.z, 1.0e-12, "Resetting the Brownian RNG seed should reproduce the z displacement"));
1004
1005 PetscCall(PicurvDestroyRuntimeContext(&simCtx));
1006 PetscCall(PicurvRemoveTempDir(tmpdir));
1007 PetscFunctionReturn(0);
1008}
1009/**
1010 * @brief Tests swarm-wide particle position updates using the same transport path as the runtime loop.
1011 */
1012static PetscErrorCode TestUpdateAllParticlePositionsMovesSwarmEntries(void)
1013{
1014 SimCtx *simCtx = NULL;
1015 UserCtx *user = NULL;
1016 PetscReal *positions = NULL;
1017 PetscReal *velocities = NULL;
1018 PetscReal *diffusivity = NULL;
1019 Cmpnts *diffusivity_gradient = NULL;
1020 PetscReal *psi = NULL;
1021 PetscReal *weights = NULL;
1022 PetscInt *cell_ids = NULL;
1023 PetscInt *status = NULL;
1024 PetscInt64 *pid = NULL;
1025
1026 PetscFunctionBeginUser;
1027 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
1028 PetscCall(PicurvCreateSwarmPair(user, 1, "ske"));
1029 simCtx->dt = 0.25;
1030
1031 PetscCall(DMSwarmGetField(user->swarm, "position", NULL, NULL, (void **)&positions));
1032 PetscCall(DMSwarmGetField(user->swarm, "velocity", NULL, NULL, (void **)&velocities));
1033 PetscCall(DMSwarmGetField(user->swarm, "Diffusivity", NULL, NULL, (void **)&diffusivity));
1034 PetscCall(DMSwarmGetField(user->swarm, "DiffusivityGradient", NULL, NULL, (void **)&diffusivity_gradient));
1035 PetscCall(DMSwarmGetField(user->swarm, "Psi", NULL, NULL, (void **)&psi));
1036 PetscCall(DMSwarmGetField(user->swarm, "weight", NULL, NULL, (void **)&weights));
1037 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1038 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1039 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_pid", NULL, NULL, (void **)&pid));
1040
1041 positions[0] = 0.20; positions[1] = 0.30; positions[2] = 0.40;
1042 velocities[0] = 0.40; velocities[1] = -0.20; velocities[2] = 0.10;
1043 diffusivity[0] = 0.0;
1044 diffusivity_gradient[0].x = 0.10;
1045 diffusivity_gradient[0].y = 0.20;
1046 diffusivity_gradient[0].z = -0.10;
1047 psi[0] = 0.5;
1048 weights[0] = 0.5; weights[1] = 0.5; weights[2] = 0.5;
1049 cell_ids[0] = 0; cell_ids[1] = 0; cell_ids[2] = 0;
1050 status[0] = ACTIVE_AND_LOCATED;
1051 pid[0] = 7;
1052
1053 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_pid", NULL, NULL, (void **)&pid));
1054 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1055 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1056 PetscCall(DMSwarmRestoreField(user->swarm, "weight", NULL, NULL, (void **)&weights));
1057 PetscCall(DMSwarmRestoreField(user->swarm, "Psi", NULL, NULL, (void **)&psi));
1058 PetscCall(DMSwarmRestoreField(user->swarm, "DiffusivityGradient", NULL, NULL, (void **)&diffusivity_gradient));
1059 PetscCall(DMSwarmRestoreField(user->swarm, "Diffusivity", NULL, NULL, (void **)&diffusivity));
1060 PetscCall(DMSwarmRestoreField(user->swarm, "velocity", NULL, NULL, (void **)&velocities));
1061 PetscCall(DMSwarmRestoreField(user->swarm, "position", NULL, NULL, (void **)&positions));
1062
1063 PetscCall(UpdateAllParticlePositions(user));
1064
1065 PetscCall(DMSwarmGetField(user->swarm, "position", NULL, NULL, (void **)&positions));
1066 PetscCall(PicurvAssertRealNear(0.325, positions[0], 1.0e-12, "UpdateAllParticlePositions should advect x"));
1067 PetscCall(PicurvAssertRealNear(0.300, positions[1], 1.0e-12, "UpdateAllParticlePositions should advect y"));
1068 PetscCall(PicurvAssertRealNear(0.400, positions[2], 1.0e-12, "UpdateAllParticlePositions should advect z"));
1069 PetscCall(DMSwarmRestoreField(user->swarm, "position", NULL, NULL, (void **)&positions));
1070
1071 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1072 PetscFunctionReturn(0);
1073}
1074/**
1075 * @brief Tests the location orchestrator fast path when a particle already carries a valid prior cell.
1076 */
1077static PetscErrorCode TestLocateAllParticlesInGridPriorCellFastPath(void)
1078{
1079 SimCtx *simCtx = NULL;
1080 UserCtx *user = NULL;
1081 PetscReal *positions = NULL;
1082 PetscReal *weights = NULL;
1083 PetscInt *cell_ids = NULL;
1084 PetscInt *status = NULL;
1085 PetscInt64 *pid = NULL;
1086
1087 PetscFunctionBeginUser;
1088 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
1089 PetscCall(PetscMalloc1(simCtx->size, &user->RankCellInfoMap));
1090 PetscCall(GetOwnedCellRange(&user->info, 0, &user->RankCellInfoMap[0].xs_cell, &user->RankCellInfoMap[0].xm_cell));
1091 PetscCall(GetOwnedCellRange(&user->info, 1, &user->RankCellInfoMap[0].ys_cell, &user->RankCellInfoMap[0].ym_cell));
1092 PetscCall(GetOwnedCellRange(&user->info, 2, &user->RankCellInfoMap[0].zs_cell, &user->RankCellInfoMap[0].zm_cell));
1093 PetscCall(PicurvCreateSwarmPair(user, 1, "ske"));
1094
1095 PetscCall(DMSwarmGetField(user->swarm, "position", NULL, NULL, (void **)&positions));
1096 PetscCall(DMSwarmGetField(user->swarm, "weight", NULL, NULL, (void **)&weights));
1097 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1098 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1099 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_pid", NULL, NULL, (void **)&pid));
1100 positions[0] = 0.375; positions[1] = 0.375; positions[2] = 0.375;
1101 weights[0] = 0.5; weights[1] = 0.5; weights[2] = 0.5;
1102 cell_ids[0] = 1; cell_ids[1] = 1; cell_ids[2] = 1;
1103 status[0] = NEEDS_LOCATION;
1104 pid[0] = 11;
1105 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_pid", NULL, NULL, (void **)&pid));
1106 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1107 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1108 PetscCall(DMSwarmRestoreField(user->swarm, "weight", NULL, NULL, (void **)&weights));
1109 PetscCall(DMSwarmRestoreField(user->swarm, "position", NULL, NULL, (void **)&positions));
1110
1111 PetscCall(LocateAllParticlesInGrid(user, simCtx->bboxlist));
1112
1113 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1114 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1115 PetscCall(PicurvAssertIntEqual(1, cell_ids[0], "prior-cell fast path should preserve the i cell id"));
1116 PetscCall(PicurvAssertIntEqual(1, cell_ids[1], "prior-cell fast path should preserve the j cell id"));
1117 PetscCall(PicurvAssertIntEqual(1, cell_ids[2], "prior-cell fast path should preserve the k cell id"));
1118 PetscCall(PicurvAssertIntEqual(ACTIVE_AND_LOCATED, status[0], "prior-cell fast path should mark the particle ACTIVE_AND_LOCATED"));
1119 PetscCall(PicurvAssertIntEqual(1, simCtx->searchMetrics.searchAttempts, "prior-cell fast path should record one search attempt"));
1120 PetscCall(PicurvAssertIntEqual(1, (PetscInt)simCtx->searchMetrics.searchPopulation, "prior-cell fast path should record one input particle"));
1121 PetscCall(PicurvAssertIntEqual(1, (PetscInt)simCtx->searchMetrics.searchLocatedCount, "prior-cell fast path should count one located particle"));
1122 PetscCall(PicurvAssertIntEqual(0, (PetscInt)simCtx->searchMetrics.searchLostCount, "prior-cell fast path should not lose the particle"));
1123 PetscCall(PicurvAssertIntEqual(0, (PetscInt)simCtx->searchMetrics.reSearchCount, "prior-cell fast path should not re-search on later passes"));
1124 PetscCall(PicurvAssertBool((PetscBool)(simCtx->searchMetrics.traversalStepsSum > 0), "prior-cell fast path should accumulate traversal steps"));
1125 PetscCall(PicurvAssertIntEqual(1, simCtx->searchMetrics.maxParticlePassDepth, "prior-cell fast path should report one settlement pass"));
1126 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1127 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1128
1129 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1130 PetscFunctionReturn(0);
1131}
1132/**
1133 * @brief Tests the guess-then-verify orchestrator path for a local particle with an unknown prior cell.
1134 */
1135static PetscErrorCode TestLocateAllParticlesInGridGuessPathResolvesLocalParticle(void)
1136{
1137 SimCtx *simCtx = NULL;
1138 UserCtx *user = NULL;
1139 PetscReal *positions = NULL;
1140 PetscReal *weights = NULL;
1141 PetscInt *cell_ids = NULL;
1142 PetscInt *status = NULL;
1143 PetscInt64 *pid = NULL;
1144
1145 PetscFunctionBeginUser;
1146 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
1147 PetscCall(PetscMalloc1(simCtx->size, &user->RankCellInfoMap));
1148 PetscCall(GetOwnedCellRange(&user->info, 0, &user->RankCellInfoMap[0].xs_cell, &user->RankCellInfoMap[0].xm_cell));
1149 PetscCall(GetOwnedCellRange(&user->info, 1, &user->RankCellInfoMap[0].ys_cell, &user->RankCellInfoMap[0].ym_cell));
1150 PetscCall(GetOwnedCellRange(&user->info, 2, &user->RankCellInfoMap[0].zs_cell, &user->RankCellInfoMap[0].zm_cell));
1151 PetscCall(PicurvCreateSwarmPair(user, 1, "ske"));
1152
1153 PetscCall(DMSwarmGetField(user->swarm, "position", NULL, NULL, (void **)&positions));
1154 PetscCall(DMSwarmGetField(user->swarm, "weight", NULL, NULL, (void **)&weights));
1155 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1156 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1157 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_pid", NULL, NULL, (void **)&pid));
1158 positions[0] = 0.625; positions[1] = 0.625; positions[2] = 0.625;
1159 weights[0] = 0.5; weights[1] = 0.5; weights[2] = 0.5;
1160 cell_ids[0] = -1; cell_ids[1] = -1; cell_ids[2] = -1;
1161 status[0] = NEEDS_LOCATION;
1162 pid[0] = 22;
1163 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_pid", NULL, NULL, (void **)&pid));
1164 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1165 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1166 PetscCall(DMSwarmRestoreField(user->swarm, "weight", NULL, NULL, (void **)&weights));
1167 PetscCall(DMSwarmRestoreField(user->swarm, "position", NULL, NULL, (void **)&positions));
1168
1169 PetscCall(LocateAllParticlesInGrid(user, simCtx->bboxlist));
1170
1171 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1172 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1173 PetscCall(PicurvAssertIntEqual(2, cell_ids[0], "guess-path location should resolve the i cell id"));
1174 PetscCall(PicurvAssertIntEqual(2, cell_ids[1], "guess-path location should resolve the j cell id"));
1175 PetscCall(PicurvAssertIntEqual(2, cell_ids[2], "guess-path location should resolve the k cell id"));
1176 PetscCall(PicurvAssertIntEqual(ACTIVE_AND_LOCATED, status[0], "guess-path location should mark the particle ACTIVE_AND_LOCATED"));
1177 PetscCall(PicurvAssertIntEqual(1, simCtx->searchMetrics.searchAttempts, "guess-path location should perform one robust search"));
1178 PetscCall(PicurvAssertIntEqual(1, (PetscInt)simCtx->searchMetrics.searchPopulation, "guess-path location should count one input particle"));
1179 PetscCall(PicurvAssertIntEqual(1, (PetscInt)simCtx->searchMetrics.searchLocatedCount, "guess-path location should count one located particle"));
1180 PetscCall(PicurvAssertIntEqual(0, (PetscInt)simCtx->searchMetrics.searchLostCount, "guess-path location should not lose the particle"));
1181 PetscCall(PicurvAssertIntEqual(0, (PetscInt)simCtx->searchMetrics.reSearchCount, "guess-path location should not count later-pass re-searches"));
1182 PetscCall(PicurvAssertIntEqual(1, simCtx->searchMetrics.bboxGuessFallbackCount, "guess-path location should record one bbox fallback"));
1183 PetscCall(PicurvAssertIntEqual(0, simCtx->searchMetrics.bboxGuessSuccessCount, "guess-path local resolution should not count as remote bbox success"));
1184 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", NULL, NULL, (void **)&status));
1185 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", NULL, NULL, (void **)&cell_ids));
1186
1187 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1188 PetscFunctionReturn(0);
1189}
1190/**
1191 * @brief Verifies that later settlement passes increment re-search metrics.
1192 */
1193static PetscErrorCode TestLocateParticleOrFindMigrationTargetCountsReSearch(void)
1194{
1195 SimCtx *simCtx = NULL;
1196 UserCtx *user = NULL;
1197 Particle particle;
1198 ParticleLocationStatus status = NEEDS_LOCATION;
1199
1200 PetscFunctionBeginUser;
1201 PetscCall(PetscMemzero(&particle, sizeof(particle)));
1202 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
1203 PetscCall(PetscMalloc1(simCtx->size, &user->RankCellInfoMap));
1204 PetscCall(GetOwnedCellRange(&user->info, 0, &user->RankCellInfoMap[0].xs_cell, &user->RankCellInfoMap[0].xm_cell));
1205 PetscCall(GetOwnedCellRange(&user->info, 1, &user->RankCellInfoMap[0].ys_cell, &user->RankCellInfoMap[0].ym_cell));
1206 PetscCall(GetOwnedCellRange(&user->info, 2, &user->RankCellInfoMap[0].zs_cell, &user->RankCellInfoMap[0].zm_cell));
1207
1208 simCtx->searchMetrics.currentSettlementPass = 2;
1209 particle.PID = 33;
1210 particle.cell[0] = 1;
1211 particle.cell[1] = 1;
1212 particle.cell[2] = 1;
1213 particle.loc.x = 0.375;
1214 particle.loc.y = 0.375;
1215 particle.loc.z = 0.375;
1216 particle.weights.x = 0.5;
1217 particle.weights.y = 0.5;
1218 particle.weights.z = 0.5;
1219
1220 PetscCall(LocateParticleOrFindMigrationTarget(user, &particle, &status));
1221
1222 PetscCall(PicurvAssertIntEqual(ACTIVE_AND_LOCATED, status, "direct re-search test should locate the particle"));
1223 PetscCall(PicurvAssertIntEqual(1, (PetscInt)simCtx->searchMetrics.searchAttempts, "direct re-search test should record one robust walk"));
1224 PetscCall(PicurvAssertIntEqual(1, (PetscInt)simCtx->searchMetrics.reSearchCount, "direct re-search test should increment re_search_count on later passes"));
1225 PetscCall(PicurvAssertIntEqual(0, (PetscInt)simCtx->searchMetrics.maxTraversalFailCount, "direct re-search test should not hit MAX_TRAVERSAL"));
1226
1227 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1228 PetscFunctionReturn(0);
1229}
1230/**
1231 * @brief Tests no-slip and free-slip wall helper kernels.
1232 */
1233
1234static PetscErrorCode TestWallNoSlipAndFreeSlipHelpers(void)
1235{
1236 Cmpnts wall_velocity = {0.0, 0.0, 0.0};
1237 Cmpnts reference_velocity = {2.0, 4.0, 6.0};
1238 Cmpnts boundary_velocity = {0.0, 0.0, 0.0};
1239 Cmpnts free_slip_reference = {2.0, 3.0, 4.0};
1240
1241 PetscFunctionBeginUser;
1242 noslip(NULL, 2.0, 1.0, wall_velocity, reference_velocity, &boundary_velocity, 1.0, 0.0, 0.0);
1243 PetscCall(PicurvAssertRealNear(1.0, boundary_velocity.x, 1.0e-12, "no-slip interpolated x"));
1244 PetscCall(PicurvAssertRealNear(2.0, boundary_velocity.y, 1.0e-12, "no-slip interpolated y"));
1245 PetscCall(PicurvAssertRealNear(3.0, boundary_velocity.z, 1.0e-12, "no-slip interpolated z"));
1246
1247 freeslip(NULL, 2.0, 1.0, wall_velocity, free_slip_reference, &boundary_velocity, 1.0, 0.0, 0.0);
1248 PetscCall(PicurvAssertRealNear(1.0, boundary_velocity.x, 1.0e-12, "free-slip interpolated normal component"));
1249 PetscCall(PicurvAssertRealNear(3.0, boundary_velocity.y, 1.0e-12, "free-slip tangential y preserved"));
1250 PetscCall(PicurvAssertRealNear(4.0, boundary_velocity.z, 1.0e-12, "free-slip tangential z preserved"));
1251 PetscFunctionReturn(0);
1252}
1253/**
1254 * @brief Tests wall-model scalar helper kernels.
1255 */
1256
1257static PetscErrorCode TestWallModelScalarHelpers(void)
1258{
1259 const PetscReal expected_smooth_e = PetscExpReal(0.41 * 5.5);
1260 PetscReal e_coeff = 0.0;
1261 PetscReal utau = 0.0;
1262 PetscReal residual = 0.0;
1263
1264 PetscFunctionBeginUser;
1265 e_coeff = E_coeff(0.1, 0.0, 1.0e-3);
1266 PetscCall(PicurvAssertRealNear(expected_smooth_e, e_coeff, 1.0e-10, "smooth-wall E coefficient"));
1267
1268 utau = find_utau_hydset(1.0e-3, 1.0, 1.0e-2, 0.1, 0.0);
1269 PetscCall(PicurvAssertBool((PetscBool)(utau > 0.0), "friction velocity should remain positive"));
1270 residual = f_hydset(1.0e-3, 1.0, 1.0e-2, utau, 0.0);
1271 PetscCall(PicurvAssertBool((PetscBool)(PetscAbsReal(residual) < 1.0e-5), "Newton solve residual should be small"));
1272
1273 PetscCall(PicurvAssertRealNear(0.0, nu_t(0.0), 1.0e-12, "eddy viscosity ratio at wall"));
1274 PetscCall(PicurvAssertBool((PetscBool)(integrate_1(1.0e-3, 1.0e-2, 0.1, 0) > 0.0), "integral helper should be positive"));
1275 PetscFunctionReturn(0);
1276}
1277/**
1278 * @brief Tests closed-form and iterative wall-model velocity helpers against inverse reconstructions.
1279 */
1280
1281static PetscErrorCode TestWallModelVelocityHelpers(void)
1282{
1283 const PetscReal kinematic_viscosity = 1.0e-3;
1284 const PetscReal wall_distance = 2.0e-2;
1285 const PetscReal target_velocity = 1.0;
1286 const PetscReal roughness_length = 1.0e-4;
1287 PetscReal utau_loglaw = 0.0;
1288 PetscReal utau_werner = 0.0;
1289 PetscReal utau_cabot = 0.0;
1290 PetscReal wall_shear_velocity = 0.0;
1291 PetscReal wall_shear_normal = 0.0;
1292
1293 PetscFunctionBeginUser;
1294 utau_loglaw = find_utau_loglaw(target_velocity, wall_distance, roughness_length);
1295 PetscCall(PicurvAssertRealNear(target_velocity, u_loglaw(wall_distance, utau_loglaw, roughness_length), 1.0e-12,
1296 "simple log-law inversion should reconstruct the target velocity"));
1297
1298 utau_werner = find_utau_Werner(kinematic_viscosity, target_velocity, wall_distance, 0.1);
1299 PetscCall(PicurvAssertBool((PetscBool)(utau_werner > 0.0), "Werner-Wengle friction velocity should remain positive"));
1300 PetscCall(PicurvAssertRealNear(target_velocity, u_Werner(kinematic_viscosity, wall_distance, utau_werner), 1.0e-12,
1301 "Werner-Wengle inversion should reconstruct the target velocity"));
1302
1303 /* The pointwise pair the solver uses must invert exactly on BOTH sides of
1304 y+ = 11.81. The previous inverse round-tripped only in the sublayer, which is the
1305 one range the sample above happens to sit in, and so passed over the mismatch. */
1306 {
1307 const PetscReal yplus[5] = {2.0, 11.0, 11.81, 30.0, 500.0};
1308 const PetscReal u_tau = 0.05;
1309 const PetscReal nu = kinematic_viscosity;
1310
1311 for (PetscInt i = 0; i < 5; ++i) {
1312 const PetscReal y = yplus[i] * nu / u_tau;
1313 const PetscReal u = u_Werner(nu, y, u_tau);
1314 PetscCall(PicurvAssertRealNear(u_tau, find_utau_Werner(nu, u, y, 0.0), 1.0e-12 * u_tau,
1315 "find_utau_Werner must invert u_Werner exactly on both branches"));
1316 }
1317 }
1318
1319 /* The cell-integral pair: correct for a cell-averaged input, not wired into the
1320 solver. Its round trip must still hold on both sides of its own switch. */
1321 {
1322 const PetscReal A = 8.3, Bexp = 1.0 / 7.0;
1323 const PetscReal u_m = 0.5 * kinematic_viscosity / wall_distance *
1324 PetscPowReal(A, 2.0 / (1.0 - Bexp));
1325 const PetscReal samples[2] = {0.25 * u_m, 4.0 * u_m};
1326
1327 for (PetscInt i = 0; i < 2; ++i) {
1328 const PetscReal tau = taw_Werner(kinematic_viscosity, samples[i], wall_distance);
1329
1330 PetscCall(PicurvAssertBool((PetscBool)(tau > 0.0),
1331 "the explicit Werner-Wengle stress must be positive"));
1332 PetscCall(PicurvAssertRealNear(
1333 samples[i], u_Werner_explicit(kinematic_viscosity, wall_distance, tau),
1334 1.0e-9 + 1.0e-9 * samples[i],
1335 "the explicit Werner-Wengle pair must round-trip on both branches"));
1336 }
1337
1338 /* The branches meet at the switch, so approaching it from either side must give
1339 the same stress. A threshold written through the friction velocity, as the
1340 iterative residual's is, cannot satisfy this. */
1341 PetscCall(PicurvAssertRealNear(
1342 taw_Werner(kinematic_viscosity, u_m * (1.0 - 1.0e-9), wall_distance),
1343 taw_Werner(kinematic_viscosity, u_m * (1.0 + 1.0e-9), wall_distance), 1.0e-9,
1344 "the Werner-Wengle branches must be continuous at the transition speed"));
1345
1346 /* Nearer the wall is slower, whatever the branch. This is the property the
1347 previous mismatched pairing violated. */
1348 {
1349 const PetscReal tau = taw_Werner(kinematic_viscosity, 4.0 * u_m, wall_distance);
1350
1351 PetscCall(PicurvAssertBool(
1352 (PetscBool)(u_Werner_explicit(kinematic_viscosity, 0.5 * wall_distance, tau) <
1353 u_Werner_explicit(kinematic_viscosity, wall_distance, tau)),
1354 "a cell nearer the wall must be assigned the slower speed"));
1355 }
1356 }
1357
1358 find_utau_Cabot(kinematic_viscosity, target_velocity, wall_distance, 0.1, 0.0, 0.0,
1359 &utau_cabot, &wall_shear_velocity, &wall_shear_normal);
1360 PetscCall(PicurvAssertBool((PetscBool)(utau_cabot > 0.0), "Cabot friction velocity should remain positive"));
1361 PetscCall(PicurvAssertRealNear(target_velocity, u_Cabot(kinematic_viscosity, wall_distance, utau_cabot, 0.0, wall_shear_velocity), 1.0e-6,
1362 "Cabot inversion should reconstruct the target velocity when pressure gradient is zero"));
1363 PetscCall(PicurvAssertRealNear(0.0, wall_shear_normal, 1.0e-10,
1364 "zero normal pressure gradient should keep Cabot normal wall shear at zero"));
1365 PetscFunctionReturn(0);
1366}
1367/**
1368 * @brief Tests the vector wall-function wrappers on a tangential reference flow.
1369 */
1370static PetscErrorCode TestWallFunctionVectorWrappers(void)
1371{
1372 SimCtx *simCtx = NULL;
1373 UserCtx *user = NULL;
1374 Cmpnts wall_velocity = {0.0, 0.0, 0.0};
1375 Cmpnts reference_velocity = {0.0, 1.0, 0.0};
1376 Cmpnts boundary_velocity = {0.0, 0.0, 0.0};
1377 PetscReal friction_velocity = 0.0;
1378
1379 PetscFunctionBeginUser;
1380 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
1381 simCtx->ren = 1000.0;
1382
1383 wall_function(user, 2.0e-2, 1.0e-2, wall_velocity, reference_velocity, &boundary_velocity, &friction_velocity, 1.0, 0.0, 0.0);
1384 PetscCall(PicurvAssertRealNear(0.0, boundary_velocity.x, 1.0e-12, "Werner wall function should preserve zero normal velocity"));
1385 PetscCall(PicurvAssertBool((PetscBool)(boundary_velocity.y > 0.0 && boundary_velocity.y < 1.0), "Werner wall function should damp tangential velocity"));
1386 PetscCall(PicurvAssertBool((PetscBool)(friction_velocity > 0.0), "Werner wall function should compute positive friction velocity"));
1387
1388 /* A reference speed taken from the pointwise profile itself must return that
1389 profile's friction velocity, and the boundary cell must land on the same profile.
1390 The reference sits in the power-law branch (y+ = 20), the boundary cell in the
1391 sublayer (y+ = 10), so both branches are exercised. The cell-integral relation
1392 read this point value as a cell average and returned u_tau 14% high. */
1393 {
1394 const PetscReal nu = 1.0 / simCtx->ren, u_tau = 1.0;
1395 Cmpnts profile_reference = {0.0, (PetscReal)u_Werner(nu, 2.0e-2, u_tau), 0.0};
1396
1397 wall_function(user, 2.0e-2, 1.0e-2, wall_velocity, profile_reference, &boundary_velocity, &friction_velocity, 1.0, 0.0, 0.0);
1398 PetscCall(PicurvAssertRealNear(u_tau, friction_velocity, 1.0e-12,
1399 "Werner wall function must return the friction velocity of the profile it is given"));
1400 PetscCall(PicurvAssertRealNear(u_Werner(nu, 1.0e-2, u_tau), boundary_velocity.y, 1.0e-12,
1401 "Werner wall function must put the boundary cell on the same pointwise profile"));
1402 }
1403
1404 wall_function_loglaw(user, 1.0e-4, 2.0e-2, 1.0e-2, wall_velocity, reference_velocity, &boundary_velocity, &friction_velocity, 1.0, 0.0, 0.0);
1405 PetscCall(PicurvAssertRealNear(0.0, boundary_velocity.x, 1.0e-12, "log-law wall function should preserve zero normal velocity"));
1406 PetscCall(PicurvAssertBool((PetscBool)(boundary_velocity.y > 0.0 && boundary_velocity.y <= 1.0), "log-law wall function should keep tangential velocity bounded"));
1407 PetscCall(PicurvAssertBool((PetscBool)(friction_velocity > 0.0), "log-law wall function should compute positive friction velocity"));
1408
1409 wall_function_Cabot(user, 1.0e-4, 2.0e-2, 1.0e-2, wall_velocity, reference_velocity, &boundary_velocity, &friction_velocity,
1410 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 10);
1411 PetscCall(PicurvAssertRealNear(0.0, boundary_velocity.x, 1.0e-12, "Cabot wall function should preserve zero normal velocity"));
1412 PetscCall(PicurvAssertBool((PetscBool)(boundary_velocity.y > 0.0 && boundary_velocity.y <= 1.0), "Cabot wall function should keep tangential velocity bounded"));
1413 PetscCall(PicurvAssertBool((PetscBool)(friction_velocity > 0.0), "Cabot wall function should compute positive friction velocity"));
1414
1415 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1416 PetscFunctionReturn(0);
1417}
1418/**
1419 * @brief Tests driven-flow validation when no driven handlers are present.
1420 */
1421
1422static PetscErrorCode TestValidateDrivenFlowConfigurationNoDrivenHandlers(void)
1423{
1424 SimCtx *simCtx = NULL;
1425 UserCtx *user = NULL;
1426
1427 PetscFunctionBeginUser;
1428 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
1429 PetscCall(Validate_DrivenFlowConfiguration(user));
1430 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1431 PetscFunctionReturn(0);
1432}
1433/**
1434 * @brief Tests that the shared minimal fixture mirrors the production DA contract.
1435 */
1436
1437static PetscErrorCode TestMinimalFixtureMirrorsProductionDMLayout(void)
1438{
1439 SimCtx *simCtx = NULL;
1440 UserCtx *user = NULL;
1441 DM coord_dm = NULL;
1442 PetscInt mx = 0, my = 0, mz = 0;
1443
1444 PetscFunctionBeginUser;
1445 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 8, 6, 4));
1446
1447 PetscCall(DMDAGetInfo(user->da, NULL, &mx, &my, &mz, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL));
1448 PetscCall(PicurvAssertIntEqual(user->IM + 1, mx, "minimal fixture should size da with IM+1 nodes"));
1449 PetscCall(PicurvAssertIntEqual(user->JM + 1, my, "minimal fixture should size da with JM+1 nodes"));
1450 PetscCall(PicurvAssertIntEqual(user->KM + 1, mz, "minimal fixture should size da with KM+1 nodes"));
1451 PetscCall(PicurvAssertIntEqual(mx, user->info.mx, "user->info should be sourced from the production da"));
1452 PetscCall(PicurvAssertIntEqual(my, user->info.my, "user->info my should match the da dimensions"));
1453 PetscCall(PicurvAssertIntEqual(mz, user->info.mz, "user->info mz should match the da dimensions"));
1454
1455 PetscCall(DMGetCoordinateDM(user->da, &coord_dm));
1456 PetscCall(PicurvAssertBool((PetscBool)(coord_dm == user->fda),
1457 "minimal fixture should derive fda from the coordinate-DM path"));
1458
1459 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1460 PetscFunctionReturn(0);
1461}
1462/**
1463 * @brief Tests that the shared swarm fixture registers the production field set.
1464 */
1465
1466static PetscErrorCode TestMinimalFixtureRegistersProductionSwarmFields(void)
1467{
1468 SimCtx *simCtx = NULL;
1469 UserCtx *user = NULL;
1470 PetscInt bs = 0;
1471 void *field_ptr = NULL;
1472
1473 PetscFunctionBeginUser;
1474 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
1475 PetscCall(PicurvCreateSwarmPair(user, 2, "ske"));
1476
1477 PetscCall(DMSwarmGetField(user->swarm, "position", &bs, NULL, &field_ptr));
1478 PetscCall(PicurvAssertIntEqual(3, bs, "solver swarm should register particle position"));
1479 PetscCall(PicurvAssertBool((PetscBool)(field_ptr != NULL), "position field should be retrievable"));
1480 PetscCall(DMSwarmRestoreField(user->swarm, "position", &bs, NULL, &field_ptr));
1481
1482 PetscCall(DMSwarmGetField(user->swarm, "velocity", &bs, NULL, &field_ptr));
1483 PetscCall(PicurvAssertIntEqual(3, bs, "solver swarm should register particle velocity"));
1484 PetscCall(DMSwarmRestoreField(user->swarm, "velocity", &bs, NULL, &field_ptr));
1485
1486 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_CellID", &bs, NULL, &field_ptr));
1487 PetscCall(PicurvAssertIntEqual(3, bs, "solver swarm should register DMSwarm_CellID"));
1488 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_CellID", &bs, NULL, &field_ptr));
1489
1490 PetscCall(DMSwarmGetField(user->swarm, "weight", &bs, NULL, &field_ptr));
1491 PetscCall(PicurvAssertIntEqual(3, bs, "solver swarm should register particle weight"));
1492 PetscCall(DMSwarmRestoreField(user->swarm, "weight", &bs, NULL, &field_ptr));
1493
1494 PetscCall(DMSwarmGetField(user->swarm, "Diffusivity", &bs, NULL, &field_ptr));
1495 PetscCall(PicurvAssertIntEqual(1, bs, "solver swarm should register particle diffusivity"));
1496 PetscCall(DMSwarmRestoreField(user->swarm, "Diffusivity", &bs, NULL, &field_ptr));
1497
1498 PetscCall(DMSwarmGetField(user->swarm, "DiffusivityGradient", &bs, NULL, &field_ptr));
1499 PetscCall(PicurvAssertIntEqual(3, bs, "solver swarm should register particle diffusivity gradients"));
1500 PetscCall(DMSwarmRestoreField(user->swarm, "DiffusivityGradient", &bs, NULL, &field_ptr));
1501
1502 PetscCall(DMSwarmGetField(user->swarm, "Psi", &bs, NULL, &field_ptr));
1503 PetscCall(PicurvAssertIntEqual(1, bs, "solver swarm should register particle scalar Psi"));
1504 PetscCall(DMSwarmRestoreField(user->swarm, "Psi", &bs, NULL, &field_ptr));
1505
1506 PetscCall(DMSwarmGetField(user->swarm, "DMSwarm_location_status", &bs, NULL, &field_ptr));
1507 PetscCall(PicurvAssertIntEqual(1, bs, "solver swarm should register particle location status"));
1508 PetscCall(DMSwarmRestoreField(user->swarm, "DMSwarm_location_status", &bs, NULL, &field_ptr));
1509
1510 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1511 PetscFunctionReturn(0);
1512}
1513/**
1514 * @brief Tests solver history-vector shifting between time levels.
1515 */
1516
1517static PetscErrorCode TestUpdateSolverHistoryVectorsShiftsStates(void)
1518{
1519 SimCtx *simCtx = NULL;
1520 UserCtx *user = NULL;
1521
1522 PetscFunctionBeginUser;
1523 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
1524
1525 PetscCall(VecSet(user->Ucont, 11.0));
1526 PetscCall(VecSet(user->Ucont_o, 7.0));
1527 PetscCall(VecSet(user->Ucont_rm1, 3.0));
1528 PetscCall(VecSet(user->Ucat, 5.0));
1529 PetscCall(VecSet(user->Ucat_o, -1.0));
1530 PetscCall(VecSet(user->P, 9.0));
1531 PetscCall(VecSet(user->P_o, -2.0));
1532
1533 PetscCall(UpdateSolverHistoryVectors(user, PETSC_FALSE));
1534
1535 PetscCall(PicurvAssertVecConstant(user->Ucont_o, 11.0, 1.0e-12, "Ucont_o should receive current Ucont"));
1536 PetscCall(PicurvAssertVecConstant(user->Ucont_rm1, 7.0, 1.0e-12, "Ucont_rm1 should receive prior Ucont_o"));
1537 PetscCall(PicurvAssertVecConstant(user->Ucat_o, 5.0, 1.0e-12, "Ucat_o should receive current Ucat"));
1538 PetscCall(PicurvAssertVecConstant(user->P_o, 9.0, 1.0e-12, "P_o should receive current P"));
1539 PetscCall(PicurvAssertVecConstant(user->lUcont_o, 11.0, 1.0e-12, "lUcont_o ghost sync should match Ucont_o"));
1540 PetscCall(PicurvAssertVecConstant(user->lUcont_rm1, 7.0, 1.0e-12, "lUcont_rm1 ghost sync should match Ucont_rm1"));
1541
1542 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1543 PetscFunctionReturn(0);
1544}
1545/**
1546 * @brief Tests owned-cell range accounting on a single MPI rank.
1547 */
1548
1549static PetscErrorCode TestGetOwnedCellRangeSingleRankAccounting(void)
1550{
1551 SimCtx *simCtx = NULL;
1552 UserCtx *user = NULL;
1553 PetscInt xs = -1, ys = -1, zs = -1;
1554 PetscInt xm = -1, ym = -1, zm = -1;
1555
1556 PetscFunctionBeginUser;
1557 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 8, 6, 4));
1558
1559 PetscCall(GetOwnedCellRange(&user->info, 0, &xs, &xm));
1560 PetscCall(GetOwnedCellRange(&user->info, 1, &ys, &ym));
1561 PetscCall(GetOwnedCellRange(&user->info, 2, &zs, &zm));
1562
1563 PetscCall(PicurvAssertIntEqual(0, xs, "single-rank x cell-start index"));
1564 PetscCall(PicurvAssertIntEqual(0, ys, "single-rank y cell-start index"));
1565 PetscCall(PicurvAssertIntEqual(0, zs, "single-rank z cell-start index"));
1566 PetscCall(PicurvAssertIntEqual(user->info.mx - 2, xm, "single-rank x owned cell count"));
1567 PetscCall(PicurvAssertIntEqual(user->info.my - 2, ym, "single-rank y owned cell count"));
1568 PetscCall(PicurvAssertIntEqual(user->info.mz - 2, zm, "single-rank z owned cell count"));
1569
1570 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1571 PetscFunctionReturn(0);
1572}
1573/**
1574 * @brief Tests neighbor-rank discovery on a single MPI rank.
1575 */
1576
1577static PetscErrorCode TestComputeAndStoreNeighborRanksSingleRank(void)
1578{
1579 SimCtx *simCtx = NULL;
1580 UserCtx *user = NULL;
1581
1582 PetscFunctionBeginUser;
1583 PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
1584 PetscCall(ComputeAndStoreNeighborRanks(user));
1585
1586 PetscCall(PicurvAssertIntEqual(MPI_PROC_NULL, user->neighbors.rank_xm, "single-rank xm neighbor should be null"));
1587 PetscCall(PicurvAssertIntEqual(MPI_PROC_NULL, user->neighbors.rank_xp, "single-rank xp neighbor should be null"));
1588 PetscCall(PicurvAssertIntEqual(MPI_PROC_NULL, user->neighbors.rank_ym, "single-rank ym neighbor should be null"));
1589 PetscCall(PicurvAssertIntEqual(MPI_PROC_NULL, user->neighbors.rank_yp, "single-rank yp neighbor should be null"));
1590 PetscCall(PicurvAssertIntEqual(MPI_PROC_NULL, user->neighbors.rank_zm, "single-rank zm neighbor should be null"));
1591 PetscCall(PicurvAssertIntEqual(MPI_PROC_NULL, user->neighbors.rank_zp, "single-rank zp neighbor should be null"));
1592
1593 PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
1594 PetscFunctionReturn(0);
1595}
1596/**
1597 * @brief Tests parsing of positive runtime walltime metadata values.
1598 */
1599
1600static PetscErrorCode TestRuntimeWalltimeGuardParsesPositiveSeconds(void)
1601{
1602 PetscReal seconds = 0.0;
1603
1604 PetscFunctionBeginUser;
1605 PetscCall(PicurvAssertBool(RuntimeWalltimeGuardParsePositiveSeconds("300", &seconds), "plain positive seconds should parse"));
1606 PetscCall(PicurvAssertRealNear(300.0, seconds, 1.0e-12, "parsed walltime seconds"));
1607 PetscCall(PicurvAssertBool(RuntimeWalltimeGuardParsePositiveSeconds(" 42.5 ", &seconds), "whitespace-wrapped decimal seconds should parse"));
1608 PetscCall(PicurvAssertRealNear(42.5, seconds, 1.0e-12, "parsed decimal walltime seconds"));
1609 PetscCall(PicurvAssertBool((PetscBool)!RuntimeWalltimeGuardParsePositiveSeconds("nope", &seconds), "non-numeric metadata should fail parsing"));
1610 PetscCall(PicurvAssertBool((PetscBool)!RuntimeWalltimeGuardParsePositiveSeconds("-10", &seconds), "negative metadata should fail parsing"));
1611 PetscFunctionReturn(0);
1612}
1613/**
1614 * @brief Tests walltime-guard estimator helper calculations.
1615 */
1616
1617static PetscErrorCode TestRuntimeWalltimeGuardEstimatorHelpers(void)
1618{
1619 PetscReal ewma_fast = 0.0;
1620 PetscReal ewma_slow = 0.0;
1621 PetscReal conservative_fast = 0.0;
1622 PetscReal conservative_slow = 0.0;
1623 PetscReal required_headroom = 0.0;
1624
1625 PetscFunctionBeginUser;
1626 ewma_fast = RuntimeWalltimeGuardUpdateEWMA(PETSC_TRUE, 4.0, 6.0, 0.5);
1627 ewma_slow = RuntimeWalltimeGuardUpdateEWMA(PETSC_TRUE, ewma_fast, 12.0, 0.5);
1628 conservative_fast = RuntimeWalltimeGuardConservativeEstimate(5.0, ewma_fast, 6.0);
1629 conservative_slow = RuntimeWalltimeGuardConservativeEstimate(5.0, ewma_slow, 12.0);
1630 required_headroom = RuntimeWalltimeGuardRequiredHeadroom(8.0, 2.0, conservative_slow);
1631
1632 PetscCall(PicurvAssertRealNear(5.0, ewma_fast, 1.0e-12, "EWMA after moderate step"));
1633 PetscCall(PicurvAssertRealNear(8.5, ewma_slow, 1.0e-12, "EWMA after newer slow step"));
1634 PetscCall(PicurvAssertRealNear(6.0, conservative_fast, 1.0e-12, "conservative estimate tracks latest moderate step"));
1635 PetscCall(PicurvAssertRealNear(12.0, conservative_slow, 1.0e-12, "conservative estimate tracks newest slow step"));
1636 PetscCall(PicurvAssertRealNear(24.0, required_headroom, 1.0e-12, "required headroom scales with conservative estimate"));
1637 PetscFunctionReturn(0);
1638}
1639/**
1640 * @brief Tests runtime walltime-guard shutdown trigger decisions.
1641 */
1642
1643static PetscErrorCode TestRuntimeWalltimeGuardTriggerDecision(void)
1644{
1645 PetscBool should_trigger = PETSC_FALSE;
1646 PetscReal required_headroom = 0.0;
1647
1648 PetscFunctionBeginUser;
1649 should_trigger = RuntimeWalltimeGuardShouldTrigger(9, 10, 15.0, 5.0, 2.0, 6.0, 6.0, 6.0, &required_headroom);
1650 PetscCall(PicurvAssertBool((PetscBool)!should_trigger, "guard should not trigger before warmup completes"));
1651
1652 should_trigger = RuntimeWalltimeGuardShouldTrigger(10, 10, 40.0, 5.0, 2.0, 10.0, 12.0, 14.0, &required_headroom);
1653 PetscCall(PicurvAssertBool((PetscBool)!should_trigger, "guard should not trigger when remaining walltime exceeds required headroom"));
1654 PetscCall(PicurvAssertRealNear(28.0, required_headroom, 1.0e-12, "required headroom after warmup"));
1655
1656 should_trigger = RuntimeWalltimeGuardShouldTrigger(10, 10, 28.0, 5.0, 2.0, 10.0, 12.0, 14.0, &required_headroom);
1657 PetscCall(PicurvAssertBool(should_trigger, "guard should trigger when remaining walltime reaches required headroom"));
1658 PetscCall(PicurvAssertRealNear(28.0, required_headroom, 1.0e-12, "required headroom remains unchanged at trigger threshold"));
1659 PetscFunctionReturn(0);
1660}
1661/**
1662 * @brief Runs the unit-runtime PETSc test binary.
1663 */
1664
1665int main(int argc, char **argv)
1666{
1667 PetscErrorCode ierr;
1668 const PicurvTestCase cases[] = {
1669 {"distribute-particles-remainder-handling", TestDistributeParticlesRemainderHandling},
1670 {"is-particle-inside-bbox-basic-cases", TestIsParticleInsideBoundingBoxBasicCases},
1671 {"update-particle-weights-computes-expected-ratios", TestUpdateParticleWeightsComputesExpectedRatios},
1672 {"update-particle-position-without-brownian-contribution", TestUpdateParticlePositionWithoutBrownianContribution},
1673 {"update-particle-position-diffusivity-gradient-only", TestUpdateParticlePositionDiffusivityGradientOnly},
1674 {"update-particle-field-iem-relaxation", TestUpdateParticleFieldIEMRelaxation},
1675 {"set-initial-interior-field-ignores-non-ucont-request", TestSetInitialInteriorFieldIgnoresNonUcontRequest},
1676 {"set-initial-interior-field-cartesian-constant-sets-contravariant-flux", TestSetInitialInteriorFieldCartesianConstantSetsContravariantFlux},
1677 {"set-initial-interior-field-curvilinear-constant-via-flow-direction", TestSetInitialInteriorFieldCurvilinearConstantViaFlowDirection},
1678 {"set-initial-interior-field-zero-clears-interior", TestSetInitialInteriorFieldZeroClearsInterior},
1679 {"set-initial-interior-field-poiseuille-profile", TestSetInitialInteriorFieldPoiseuilleProfile},
1680 {"set-initial-interior-field-poiseuille-uniform-along-periodic-axis", TestSetInitialInteriorFieldPoiseuilleUniformAlongPeriodicAxis},
1681 {"cart2contra-converts-cartesian-field", TestCart2ContraConvertsCartesianField},
1682 {"cart2contra-uses-finalized-periodic-ucat", TestCart2ContraUsesFinalizedPeriodicUcat},
1683 {"populate-initial-ucont-loads-staged-ucat", TestPopulateInitialUcontLoadsStagedUcat},
1684 {"write-les-fields-copies-owned-periodic-values", TestWriteLESFieldsCopiesOwnedPeriodicValues},
1685 {"staged-ucat-fourier-modes-use-picurv-discrete-symbol", TestStagedUcatFourierModesUsePicurvDiscreteSymbol},
1686 {"populate-initial-ucont-loads-staged-ucont", TestPopulateInitialUcontLoadsStagedUcont},
1687 {"interpolate-all-fields-to-swarm-constant-fields", TestInterpolateAllFieldsToSwarmConstantFields},
1688 {"interpolate-all-fields-to-swarm-corner-averaged-constant-fields", TestInterpolateAllFieldsToSwarmCornerAveragedConstantFields},
1689 {"scatter-all-particle-fields-to-euler-fields-averages-psi", TestScatterAllParticleFieldsToEulerFieldsAveragesPsi},
1690 {"calculate-particle-count-per-cell-counts-global-cell-ids", TestCalculateParticleCountPerCellCountsGlobalCellIDs},
1691 {"reset-all-particle-statuses-leaves-lost-particles-untouched", TestResetAllParticleStatusesLeavesLostParticlesUntouched},
1692 {"check-and-remove-out-of-bounds-particles-removes-escaped-particle", TestCheckAndRemoveOutOfBoundsParticlesRemovesEscapedParticle},
1693 {"check-and-remove-lost-particles-removes-lost-entries", TestCheckAndRemoveLostParticlesRemovesLostEntries},
1694 {"calculate-brownian-displacement-deterministic-seed", TestCalculateBrownianDisplacementDeterministicSeed},
1695 {"update-all-particle-positions-moves-swarm-entries", TestUpdateAllParticlePositionsMovesSwarmEntries},
1696 {"locate-all-particles-in-grid-prior-cell-fast-path", TestLocateAllParticlesInGridPriorCellFastPath},
1697 {"locate-all-particles-in-grid-guess-path-resolves-local-particle", TestLocateAllParticlesInGridGuessPathResolvesLocalParticle},
1698 {"locate-particle-or-find-migration-target-counts-research", TestLocateParticleOrFindMigrationTargetCountsReSearch},
1699 {"wall-noslip-and-freeslip-helpers", TestWallNoSlipAndFreeSlipHelpers},
1700 {"wall-model-scalar-helpers", TestWallModelScalarHelpers},
1701 {"wall-model-velocity-helpers", TestWallModelVelocityHelpers},
1702 {"wall-function-vector-wrappers", TestWallFunctionVectorWrappers},
1703 {"validate-driven-flow-configuration-no-driven-handlers", TestValidateDrivenFlowConfigurationNoDrivenHandlers},
1704 {"minimal-fixture-mirrors-production-dm-layout", TestMinimalFixtureMirrorsProductionDMLayout},
1705 {"minimal-fixture-registers-production-swarm-fields", TestMinimalFixtureRegistersProductionSwarmFields},
1706 {"update-solver-history-vectors-shifts-states", TestUpdateSolverHistoryVectorsShiftsStates},
1707 {"get-owned-cell-range-single-rank-accounting", TestGetOwnedCellRangeSingleRankAccounting},
1708 {"compute-and-store-neighbor-ranks-single-rank", TestComputeAndStoreNeighborRanksSingleRank},
1709 {"runtime-walltime-guard-parses-positive-seconds", TestRuntimeWalltimeGuardParsesPositiveSeconds},
1710 {"runtime-walltime-guard-estimator-helpers", TestRuntimeWalltimeGuardEstimatorHelpers},
1711 {"runtime-walltime-guard-trigger-decision", TestRuntimeWalltimeGuardTriggerDecision},
1712 };
1713
1714 ierr = PetscInitialize(&argc, &argv, NULL, "PICurv runtime-kernel tests");
1715 if (ierr) {
1716 return (int)ierr;
1717 }
1718
1719 ierr = PicurvRunTests("unit-runtime", cases, sizeof(cases) / sizeof(cases[0]));
1720 if (ierr) {
1721 PetscFinalize();
1722 return (int)ierr;
1723 }
1724
1725 ierr = PetscFinalize();
1726 return (int)ierr;
1727}
BC_Handlers.h
Validate_DrivenFlowConfiguration
PetscErrorCode Validate_DrivenFlowConfiguration(UserCtx *user)
(Private) Validates all consistency rules for a driven flow (channel/pipe) setup.
Definition BC_Handlers.c:15
SynchronizePeriodicStaggeredFields
PetscErrorCode SynchronizePeriodicStaggeredFields(UserCtx *user, PetscInt num_fields, const FieldId field_ids[])
Synchronizes persistent component-staggered vector fields.
Definition Boundaries.c:2015
SynchronizePeriodicCellFields
PetscErrorCode SynchronizePeriodicCellFields(UserCtx *user, PetscInt num_fields, const FieldId field_ids[])
Synchronizes periodic endpoint cells for a list of cell-centered fields.
Definition Boundaries.c:1647
ParticleMotion.h
Header file for Particle Motion and migration related functions.
CheckAndRemoveLostParticles
PetscErrorCode CheckAndRemoveLostParticles(UserCtx *user, PetscInt *removedCountLocal, PetscInt *removedCountGlobal, PetscReal *removedScalarSumGlobal)
Removes particles that have been definitively flagged as LOST by the location algorithm.
Definition ParticleMotion.c:336
CheckAndRemoveOutOfBoundsParticles
PetscErrorCode CheckAndRemoveOutOfBoundsParticles(UserCtx *user, PetscInt *removedCountLocal, PetscInt *removedCountGlobal, const BoundingBox *bboxlist)
Checks for particles outside the physical domain boundaries and removes them using DMSwarmRemovePoint...
Definition ParticleMotion.c:228
UpdateAllParticlePositions
PetscErrorCode UpdateAllParticlePositions(UserCtx *user)
Loops over all local particles in the DMSwarm, updating their positions based on velocity and the glo...
Definition ParticleMotion.c:134
CalculateParticleCountPerCell
PetscErrorCode CalculateParticleCountPerCell(UserCtx *user)
Counts particles in each cell of the DMDA 'da' and stores the result in user->ParticleCount.
Definition ParticleMotion.c:521
CalculateBrownianDisplacement
PetscErrorCode CalculateBrownianDisplacement(UserCtx *user, PetscReal diff_eff, Cmpnts *displacement)
Calculates the stochastic displacement vector (Brownian motion) for a single particle.
Definition ParticleMotion.c:57
LocateAllParticlesInGrid
PetscErrorCode LocateAllParticlesInGrid(UserCtx *user, BoundingBox *bboxlist)
Orchestrates the complete particle location and migration process for one timestep.
Definition ParticleMotion.c:1399
UpdateParticlePosition
PetscErrorCode UpdateParticlePosition(UserCtx *user, Particle *particle)
Updates a particle's position based on its velocity and the timestep dt (stored in user->dt).
Definition ParticleMotion.c:105
ResetAllParticleStatuses
PetscErrorCode ResetAllParticleStatuses(UserCtx *user)
Marks all local particles as NEEDS_LOCATION for the next settlement pass.
Definition ParticleMotion.c:1633
ParticlePhysics.h
Header file for Particle related physics modules.
UpdateParticleField
PetscErrorCode UpdateParticleField(ParticleFieldId field_id, PetscReal dt, PetscReal *psi_io, PetscReal diffusivity, PetscReal mean_val, PetscReal cell_vol, PetscReal C_model)
Updates a single particle's field based on its state and physics model.
Definition ParticlePhysics.c:14
ParticleSwarm.h
Header file for Particle Swarm management functions.
UpdateParticleWeights
PetscErrorCode UpdateParticleWeights(PetscReal *d, Particle *particle)
Updates a particle's interpolation weights based on distances to cell faces.
Definition ParticleSwarm.c:1012
DistributeParticles
PetscErrorCode DistributeParticles(PetscInt numParticles, PetscMPIInt rank, PetscMPIInt size, PetscInt *particlesPerProcess, PetscInt *remainder)
Distributes particles evenly across MPI processes, handling any remainders.
Definition ParticleSwarm.c:577
IsParticleInsideBoundingBox
PetscBool IsParticleInsideBoundingBox(const BoundingBox *bbox, const Particle *particle)
Checks if a particle's location is within a specified bounding box.
Definition ParticleSwarm.c:961
FieldId
FieldId
Compile-time identity for a catalogued Eulerian field.
Definition field_catalog.h:23
FIELD_ID_UCAT
@ FIELD_ID_UCAT
Definition field_catalog.h:26
FIELD_ID_UCONT
@ FIELD_ID_UCONT
Definition field_catalog.h:27
FIELD_ID_P
@ FIELD_ID_P
Definition field_catalog.h:30
ScatterAllParticleFieldsToEulerFields
PetscErrorCode ScatterAllParticleFieldsToEulerFields(UserCtx *user)
Scatters a predefined set of particle fields to their corresponding Eulerian fields.
Definition interpolation.c:1888
SetInitialInteriorField
PetscErrorCode SetInitialInteriorField(UserCtx *user, FieldId field_id)
Sets the initial values for the INTERIOR of a specified Eulerian field.
Definition initialcondition.c:14
PopulateInitialUcont
PetscErrorCode PopulateInitialUcont(UserCtx *user)
Populate Ucont for one fresh-start block from the configured IC mode.
Definition initialcondition.c:231
InterpolateAllFieldsToSwarm
PetscErrorCode InterpolateAllFieldsToSwarm(UserCtx *user)
Interpolates all relevant fields from the DMDA to the DMSwarm.
Definition interpolation.c:1400
WriteSimulationFields
PetscErrorCode WriteSimulationFields(UserCtx *user, const char *checkpoint_directory)
Writes simulation fields to files.
Definition io.c:2021
WriteFieldData
PetscErrorCode WriteFieldData(UserCtx *user, const char *field_name, Vec field_vec, const char *ext)
Writes data from a specific PETSc vector to a file.
Definition io.c:1966
PARTICLE_FIELD_ID_PSI
@ PARTICLE_FIELD_ID_PSI
Definition particle_field_catalog.h:31
PARTICLE_FIELD_ID_VELOCITY
@ PARTICLE_FIELD_ID_VELOCITY
Definition particle_field_catalog.h:26
RuntimeWalltimeGuardUpdateEWMA
PetscReal RuntimeWalltimeGuardUpdateEWMA(PetscBool has_previous, PetscReal previous_ewma_seconds, PetscReal latest_step_seconds, PetscReal alpha)
Update an EWMA estimate for timestep wall-clock duration.
Definition runloop.c:153
RuntimeWalltimeGuardConservativeEstimate
PetscReal RuntimeWalltimeGuardConservativeEstimate(PetscReal warmup_average_seconds, PetscReal ewma_seconds, PetscReal latest_step_seconds)
Return the conservative timestep estimate used by the walltime guard.
Definition runloop.c:165
RuntimeWalltimeGuardRequiredHeadroom
PetscReal RuntimeWalltimeGuardRequiredHeadroom(PetscReal min_seconds, PetscReal multiplier, PetscReal conservative_estimate_seconds)
Compute the required shutdown headroom from timestep estimate and floor.
Definition runloop.c:176
RuntimeWalltimeGuardShouldTrigger
PetscBool RuntimeWalltimeGuardShouldTrigger(PetscInt completed_steps, PetscInt warmup_steps, PetscReal remaining_seconds, PetscReal min_seconds, PetscReal multiplier, PetscReal warmup_average_seconds, PetscReal ewma_seconds, PetscReal latest_step_seconds, PetscReal *required_headroom_seconds_out)
Decide whether the runtime walltime guard should stop before another step.
Definition runloop.c:187
UpdateSolverHistoryVectors
PetscErrorCode UpdateSolverHistoryVectors(UserCtx *user, PetscBool preserve_previous_state)
Copies the current time step's solution fields into history vectors (e.g., U(t_n) -> U_o,...
Definition runloop.c:314
GetOwnedCellRange
PetscErrorCode GetOwnedCellRange(const DMDALocalInfo *info_nodes, PetscInt dim, PetscInt *xs_cell_global_out, PetscInt *xm_cell_local_out)
Determines the global starting index and number of CELLS owned by the current processor in a specifie...
Definition setup.c:2936
ComputeAndStoreNeighborRanks
PetscErrorCode ComputeAndStoreNeighborRanks(UserCtx *user)
Computes and stores the Cartesian neighbor ranks for the DMDA decomposition.
Definition setup.c:3033
UpdateLocalGhosts
PetscErrorCode UpdateLocalGhosts(UserCtx *user, FieldId field_id)
Updates the local vector (including ghost points) from its corresponding global vector.
Definition setup.c:2489
Cart2Contra
PetscErrorCode Cart2Contra(UserCtx *user)
Convert the ghosted Cartesian velocity field to contravariant face fluxes.
Definition setup.c:3435
RuntimeWalltimeGuardParsePositiveSeconds
PetscBool RuntimeWalltimeGuardParsePositiveSeconds(const char *text, PetscReal *seconds_out)
Parse a positive floating-point seconds value from runtime metadata.
Definition setup.c:67
TestUpdateParticleFieldIEMRelaxation
static PetscErrorCode TestUpdateParticleFieldIEMRelaxation(void)
Tests IEM relaxation updates for particle-carried fields.
Definition test_runtime_kernels.c:227
TestMinimalFixtureMirrorsProductionDMLayout
static PetscErrorCode TestMinimalFixtureMirrorsProductionDMLayout(void)
Tests that the shared minimal fixture mirrors the production DA contract.
Definition test_runtime_kernels.c:1437
TestUpdateParticlePositionDiffusivityGradientOnly
static PetscErrorCode TestUpdateParticlePositionDiffusivityGradientOnly(void)
Tests particle position updates driven only by diffusivity-gradient drift.
Definition test_runtime_kernels.c:193
TestLocateAllParticlesInGridGuessPathResolvesLocalParticle
static PetscErrorCode TestLocateAllParticlesInGridGuessPathResolvesLocalParticle(void)
Tests the guess-then-verify orchestrator path for a local particle with an unknown prior cell.
Definition test_runtime_kernels.c:1135
TestWallModelVelocityHelpers
static PetscErrorCode TestWallModelVelocityHelpers(void)
Tests closed-form and iterative wall-model velocity helpers against inverse reconstructions.
Definition test_runtime_kernels.c:1281
TestSetInitialInteriorFieldIgnoresNonUcontRequest
static PetscErrorCode TestSetInitialInteriorFieldIgnoresNonUcontRequest(void)
Tests that non-Ucont requests do not modify interior field initialization.
Definition test_runtime_kernels.c:253
TestLocateParticleOrFindMigrationTargetCountsReSearch
static PetscErrorCode TestLocateParticleOrFindMigrationTargetCountsReSearch(void)
Verifies that later settlement passes increment re-search metrics.
Definition test_runtime_kernels.c:1193
TestRuntimeWalltimeGuardParsesPositiveSeconds
static PetscErrorCode TestRuntimeWalltimeGuardParsesPositiveSeconds(void)
Tests parsing of positive runtime walltime metadata values.
Definition test_runtime_kernels.c:1600
TestCart2ContraUsesFinalizedPeriodicUcat
static PetscErrorCode TestCart2ContraUsesFinalizedPeriodicUcat(void)
Verifies periodic Ucat endpoints are repaired before Cart2Contra terminal faces.
Definition test_runtime_kernels.c:462
TestSetInitialInteriorFieldCurvilinearConstantViaFlowDirection
static PetscErrorCode TestSetInitialInteriorFieldCurvilinearConstantViaFlowDirection(void)
Tests curvilinear Constant IC: flow_direction selects the streamwise axis.
Definition test_runtime_kernels.c:306
main
int main(int argc, char **argv)
Runs the unit-runtime PETSc test binary.
Definition test_runtime_kernels.c:1665
TestInterpolateAllFieldsToSwarmConstantFields
static PetscErrorCode TestInterpolateAllFieldsToSwarmConstantFields(void)
Tests direct interpolation from Eulerian fields to one localized swarm particle.
Definition test_runtime_kernels.c:706
SyncRuntimeFieldGhosts
static PetscErrorCode SyncRuntimeFieldGhosts(UserCtx *user)
Synchronizes the minimal runtime fixture's global fields into their persistent local ghosts.
Definition test_runtime_kernels.c:23
TestStagedUcatFourierModesUsePicurvDiscreteSymbol
static PetscErrorCode TestStagedUcatFourierModesUsePicurvDiscreteSymbol(void)
Cross-checks selected exact Fourier modes against PICurv's native discrete symbol.
Definition test_runtime_kernels.c:630
TestDistributeParticlesRemainderHandling
static PetscErrorCode TestDistributeParticlesRemainderHandling(void)
Tests particle distribution remainder handling across ranks.
Definition test_runtime_kernels.c:88
TestResetAllParticleStatusesLeavesLostParticlesUntouched
static PetscErrorCode TestResetAllParticleStatusesLeavesLostParticlesUntouched(void)
Tests localized particle-status reset behavior for restart of the location workflow.
Definition test_runtime_kernels.c:877
SeedSingleParticle
static PetscErrorCode SeedSingleParticle(UserCtx *user, PetscInt ci, PetscInt cj, PetscInt ck, PetscReal x, PetscReal y, PetscReal z, PetscReal wx, PetscReal wy, PetscReal wz, PetscInt status_value)
Seeds one localized swarm particle with the cell, position, weight, and status data used by runtime t...
Definition test_runtime_kernels.c:44
TestSetInitialInteriorFieldCartesianConstantSetsContravariantFlux
static PetscErrorCode TestSetInitialInteriorFieldCartesianConstantSetsContravariantFlux(void)
Tests cartesian Constant IC: Cart2Contra sets contravariant flux via metric dot product.
Definition test_runtime_kernels.c:274
TestPopulateInitialUcontLoadsStagedUcat
static PetscErrorCode TestPopulateInitialUcontLoadsStagedUcat(void)
Tests loading a staged Ucat file IC and converting it to Ucont.
Definition test_runtime_kernels.c:492
TestPopulateInitialUcontLoadsStagedUcont
static PetscErrorCode TestPopulateInitialUcontLoadsStagedUcont(void)
Tests loading a staged Ucont file IC without Cartesian conversion.
Definition test_runtime_kernels.c:676
TestInterpolateAllFieldsToSwarmCornerAveragedConstantFields
static PetscErrorCode TestInterpolateAllFieldsToSwarmCornerAveragedConstantFields(void)
Tests the corner-averaged (legacy) interpolation path on constant fields.
Definition test_runtime_kernels.c:758
TestGetOwnedCellRangeSingleRankAccounting
static PetscErrorCode TestGetOwnedCellRangeSingleRankAccounting(void)
Tests owned-cell range accounting on a single MPI rank.
Definition test_runtime_kernels.c:1549
TestUpdateParticlePositionWithoutBrownianContribution
static PetscErrorCode TestUpdateParticlePositionWithoutBrownianContribution(void)
Tests particle position updates without Brownian forcing.
Definition test_runtime_kernels.c:159
TestComputeAndStoreNeighborRanksSingleRank
static PetscErrorCode TestComputeAndStoreNeighborRanksSingleRank(void)
Tests neighbor-rank discovery on a single MPI rank.
Definition test_runtime_kernels.c:1577
TestCalculateParticleCountPerCellCountsGlobalCellIDs
static PetscErrorCode TestCalculateParticleCountPerCellCountsGlobalCellIDs(void)
Tests particle counting by geometric cell IDs using the production +1 storage shift.
Definition test_runtime_kernels.c:846
TestWriteLESFieldsCopiesOwnedPeriodicValues
static PetscErrorCode TestWriteLESFieldsCopiesOwnedPeriodicValues(void)
Tests LES output staging from local vectors on a multiply periodic DMDA.
Definition test_runtime_kernels.c:527
TestValidateDrivenFlowConfigurationNoDrivenHandlers
static PetscErrorCode TestValidateDrivenFlowConfigurationNoDrivenHandlers(void)
Tests driven-flow validation when no driven handlers are present.
Definition test_runtime_kernels.c:1422
TestWallNoSlipAndFreeSlipHelpers
static PetscErrorCode TestWallNoSlipAndFreeSlipHelpers(void)
Tests no-slip and free-slip wall helper kernels.
Definition test_runtime_kernels.c:1234
TestMinimalFixtureRegistersProductionSwarmFields
static PetscErrorCode TestMinimalFixtureRegistersProductionSwarmFields(void)
Tests that the shared swarm fixture registers the production field set.
Definition test_runtime_kernels.c:1466
TestWallFunctionVectorWrappers
static PetscErrorCode TestWallFunctionVectorWrappers(void)
Tests the vector wall-function wrappers on a tangential reference flow.
Definition test_runtime_kernels.c:1370
TestSetInitialInteriorFieldZeroClearsInterior
static PetscErrorCode TestSetInitialInteriorFieldZeroClearsInterior(void)
Tests zero IC clears physical-cell contravariant velocity.
Definition test_runtime_kernels.c:338
TestLocateAllParticlesInGridPriorCellFastPath
static PetscErrorCode TestLocateAllParticlesInGridPriorCellFastPath(void)
Tests the location orchestrator fast path when a particle already carries a valid prior cell.
Definition test_runtime_kernels.c:1077
TestCalculateBrownianDisplacementDeterministicSeed
static PetscErrorCode TestCalculateBrownianDisplacementDeterministicSeed(void)
Tests Brownian displacement generation against a duplicated seeded RNG stream.
Definition test_runtime_kernels.c:979
TestScatterAllParticleFieldsToEulerFieldsAveragesPsi
static PetscErrorCode TestScatterAllParticleFieldsToEulerFieldsAveragesPsi(void)
Tests particle-to-grid scattering using known cell occupancy and scalar values.
Definition test_runtime_kernels.c:811
TestUpdateSolverHistoryVectorsShiftsStates
static PetscErrorCode TestUpdateSolverHistoryVectorsShiftsStates(void)
Tests solver history-vector shifting between time levels.
Definition test_runtime_kernels.c:1517
TestUpdateParticleWeightsComputesExpectedRatios
static PetscErrorCode TestUpdateParticleWeightsComputesExpectedRatios(void)
Tests particle weight updates against expected ratios.
Definition test_runtime_kernels.c:136
TestCheckAndRemoveOutOfBoundsParticlesRemovesEscapedParticle
static PetscErrorCode TestCheckAndRemoveOutOfBoundsParticlesRemovesEscapedParticle(void)
Tests direct removal of particles that leave every rank bounding box.
Definition test_runtime_kernels.c:908
StagedFourierModeMaxDivergence
static PetscErrorCode StagedFourierModeMaxDivergence(UserCtx *user, PetscInt mode_x, PetscInt mode_y, PetscInt mode_z, Cmpnts amplitude, const char *tmpdir, PetscReal *max_divergence)
Stages a real Cartesian Fourier mode through the production file-IC path.
Definition test_runtime_kernels.c:572
TestSetInitialInteriorFieldPoiseuilleProfile
static PetscErrorCode TestSetInitialInteriorFieldPoiseuilleProfile(void)
Tests Poiseuille IC is the product parabola that vanishes on the walls.
Definition test_runtime_kernels.c:365
TestSetInitialInteriorFieldPoiseuilleUniformAlongPeriodicAxis
static PetscErrorCode TestSetInitialInteriorFieldPoiseuilleUniformAlongPeriodicAxis(void)
Tests Poiseuille IC is uniform along a periodic cross-stream axis.
Definition test_runtime_kernels.c:398
TestCart2ContraConvertsCartesianField
static PetscErrorCode TestCart2ContraConvertsCartesianField(void)
Tests spatially varying Cartesian-field conversion to contravariant fluxes.
Definition test_runtime_kernels.c:430
TestCheckAndRemoveLostParticlesRemovesLostEntries
static PetscErrorCode TestCheckAndRemoveLostParticlesRemovesLostEntries(void)
Tests direct removal of particles already marked LOST by the location workflow.
Definition test_runtime_kernels.c:939
TestRuntimeWalltimeGuardEstimatorHelpers
static PetscErrorCode TestRuntimeWalltimeGuardEstimatorHelpers(void)
Tests walltime-guard estimator helper calculations.
Definition test_runtime_kernels.c:1617
TestIsParticleInsideBoundingBoxBasicCases
static PetscErrorCode TestIsParticleInsideBoundingBoxBasicCases(void)
Tests basic particle-inside-bounding-box classification cases.
Definition test_runtime_kernels.c:107
TestRuntimeWalltimeGuardTriggerDecision
static PetscErrorCode TestRuntimeWalltimeGuardTriggerDecision(void)
Tests runtime walltime-guard shutdown trigger decisions.
Definition test_runtime_kernels.c:1643
TestUpdateAllParticlePositionsMovesSwarmEntries
static PetscErrorCode TestUpdateAllParticlePositionsMovesSwarmEntries(void)
Tests swarm-wide particle position updates using the same transport path as the runtime loop.
Definition test_runtime_kernels.c:1012
TestWallModelScalarHelpers
static PetscErrorCode TestWallModelScalarHelpers(void)
Tests wall-model scalar helper kernels.
Definition test_runtime_kernels.c:1257
PicurvMakeTempDir
PetscErrorCode PicurvMakeTempDir(char *path, size_t path_len)
Creates a unique temporary directory for one test case.
Definition test_support.c:168
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:580
PicurvEnsureDir
PetscErrorCode PicurvEnsureDir(const char *path)
Ensures a directory exists for test output.
Definition test_support.c:156
PicurvAssertRealNear
PetscErrorCode PicurvAssertRealNear(PetscReal expected, PetscReal actual, PetscReal tol, const char *context)
Asserts that two real values agree within tolerance.
Definition test_support.c:1063
PicurvDestroyMinimalContexts
PetscErrorCode PicurvDestroyMinimalContexts(SimCtx **simCtx_ptr, UserCtx **user_ptr)
Destroys minimal SimCtx/UserCtx fixtures and all owned PETSc objects.
Definition test_support.c:927
PicurvCreateMinimalContextsWithPeriodicity
PetscErrorCode PicurvCreateMinimalContextsWithPeriodicity(SimCtx **simCtx_out, UserCtx **user_out, PetscInt mx, PetscInt my, PetscInt mz, PetscBool x_periodic, PetscBool y_periodic, PetscBool z_periodic)
Builds minimal SimCtx and UserCtx fixtures for C unit tests with configurable periodicity.
Definition test_support.c:339
PicurvCreateSwarmPair
PetscErrorCode PicurvCreateSwarmPair(UserCtx *user, PetscInt nlocal, const char *post_field_name)
Creates matched solver and post-processing swarms for tests.
Definition test_support.c:748
PicurvDestroyRuntimeContext
PetscErrorCode PicurvDestroyRuntimeContext(SimCtx **simCtx_ptr)
Finalizes and frees a runtime context built by PicurvBuildTinyRuntimeContext.
Definition test_support.c:913
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
PicurvBuildTinyRuntimeContext
PetscErrorCode PicurvBuildTinyRuntimeContext(const char *bcs_contents, PetscBool enable_particles, SimCtx **simCtx_out, UserCtx **user_out, char *tmpdir, size_t tmpdir_len)
Builds a tiny runtime context through the real setup path for behavior-level tests.
Definition test_support.c:785
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:1124
PicurvAssertIntEqual
PetscErrorCode PicurvAssertIntEqual(PetscInt expected, PetscInt actual, const char *context)
Asserts that two integer values are equal.
Definition test_support.c:1080
PicurvAssertBool
PetscErrorCode PicurvAssertBool(PetscBool value, const char *context)
Asserts that one boolean condition is true.
Definition test_support.c:1095
PicurvRemoveTempDir
PetscErrorCode PicurvRemoveTempDir(const char *path)
Recursively removes a temporary directory created by PicurvMakeTempDir.
Definition test_support.c:185
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
RankNeighbors::rank_zm
PetscMPIInt rank_zm
Definition variables.h:225
DYNAMIC_SMAGORINSKY
@ DYNAMIC_SMAGORINSKY
Definition variables.h:551
SimCtx::icVelocityPhysical
PetscReal icVelocityPhysical
Definition variables.h:921
UserCtx::lDiffusivityGradient
Vec lDiffusivityGradient
Definition variables.h:1117
PERIODIC
@ PERIODIC
Definition variables.h:318
Particle::vel
Cmpnts vel
Definition variables.h:212
RankCellInfo::ys_cell
PetscInt ys_cell
Definition variables.h:230
RankCellInfo::xs_cell
PetscInt xs_cell
Definition variables.h:230
UserCtx::boundary_faces
BoundaryFaceConfig boundary_faces[6]
Definition variables.h:1099
RankNeighbors::rank_yp
PetscMPIInt rank_yp
Definition variables.h:224
SearchMetricsState::searchLocatedCount
PetscInt64 searchLocatedCount
Definition variables.h:267
SearchMetricsState::searchLostCount
PetscInt64 searchLostCount
Definition variables.h:268
Particle::cell
PetscInt cell[3]
Definition variables.h:210
SimCtx::initialConditionMode
InitialConditionMode initialConditionMode
Definition variables.h:916
UserCtx::simCtx
SimCtx * simCtx
Back-pointer to the master simulation context.
Definition variables.h:1077
ParticleLocationStatus
ParticleLocationStatus
Defines the state of a particle with respect to its location and migration status during the iterativ...
Definition variables.h:163
LOST
@ LOST
Definition variables.h:167
NEEDS_LOCATION
@ NEEDS_LOCATION
Definition variables.h:164
ACTIVE_AND_LOCATED
@ ACTIVE_AND_LOCATED
Definition variables.h:165
RankNeighbors::rank_ym
PetscMPIInt rank_ym
Definition variables.h:224
SimCtx::flowDirection
FlowDirection flowDirection
Definition variables.h:920
RankNeighbors::rank_xp
PetscMPIInt rank_xp
Definition variables.h:223
UserCtx::KM
PetscInt KM
Definition variables.h:1088
SearchMetricsState::traversalStepsSum
PetscInt64 traversalStepsSum
Definition variables.h:269
SimCtx::ren
PetscReal ren
Definition variables.h:906
UserCtx::lUcont_rm1
Vec lUcont_rm1
Definition variables.h:1121
RankCellInfo::zm_cell
PetscInt zm_cell
Definition variables.h:231
BoundingBox::max_coords
Cmpnts max_coords
Maximum x, y, z coordinates of the bounding box.
Definition variables.h:199
RankCellInfo::zs_cell
PetscInt zs_cell
Definition variables.h:230
Particle::diffusivitygradient
Cmpnts diffusivitygradient
Definition variables.h:217
SearchMetricsState::searchPopulation
PetscInt64 searchPopulation
Definition variables.h:266
SimCtx::dt
PetscReal dt
Definition variables.h:874
UserCtx::neighbors
RankNeighbors neighbors
Definition variables.h:1091
SearchMetricsState::currentSettlementPass
PetscInt currentSettlementPass
Definition variables.h:278
UserCtx::DiffusivityGradient
Vec DiffusivityGradient
Definition variables.h:1117
UserCtx::Ucont
Vec Ucont
Definition variables.h:1113
BoundingBox::min_coords
Cmpnts min_coords
Minimum x, y, z coordinates of the bounding box.
Definition variables.h:198
Cmpnts::x
PetscScalar x
Definition variables.h:122
Particle::loc
Cmpnts loc
Definition variables.h:211
SearchMetricsState::reSearchCount
PetscInt64 reSearchCount
Definition variables.h:270
SimCtx::current_io_directory
char * current_io_directory
Definition variables.h:885
RankCellInfo::xm_cell
PetscInt xm_cell
Definition variables.h:231
UserCtx::lUcont_o
Vec lUcont_o
Definition variables.h:1120
SearchMetricsState::bboxGuessFallbackCount
PetscInt64 bboxGuessFallbackCount
Definition variables.h:276
SimCtx::interpolationMethod
InterpolationMethod interpolationMethod
Definition variables.h:995
UserCtx::RankCellInfoMap
RankCellInfo * RankCellInfoMap
Definition variables.h:1170
RankCellInfo::ym_cell
PetscInt ym_cell
Definition variables.h:231
UserCtx::Ucat_o
Vec Ucat_o
Definition variables.h:1120
SearchMetricsState::bboxGuessSuccessCount
PetscInt64 bboxGuessSuccessCount
Definition variables.h:275
SimCtx::bboxlist
BoundingBox * bboxlist
Definition variables.h:993
UserCtx::lNu_t
Vec lNu_t
Definition variables.h:1156
RankNeighbors::rank_xm
PetscMPIInt rank_xm
Definition variables.h:223
SearchMetricsState::maxParticlePassDepth
PetscInt64 maxParticlePassDepth
Definition variables.h:277
Cmpnts::z
PetscScalar z
Definition variables.h:122
INTERP_CORNER_AVERAGED
@ INTERP_CORNER_AVERAGED
Definition variables.h:724
UserCtx::ParticleCount
Vec ParticleCount
Definition variables.h:1171
UserCtx::Ucont_o
Vec Ucont_o
Definition variables.h:1120
UserCtx::JM
PetscInt JM
Definition variables.h:1088
FLOW_DIR_POS_ZETA
@ FLOW_DIR_POS_ZETA
Definition variables.h:303
SimCtx::initialConditionDirectory
char initialConditionDirectory[PETSC_MAX_PATH_LEN]
Definition variables.h:918
IC_MODE_CONSTANT_CARTESIAN
@ IC_MODE_CONSTANT_CARTESIAN
Definition variables.h:179
IC_MODE_POISEUILLE
@ IC_MODE_POISEUILLE
Definition variables.h:180
IC_MODE_CONSTANT_STREAMWISE
@ IC_MODE_CONSTANT_STREAMWISE
Definition variables.h:181
IC_MODE_FILE
@ IC_MODE_FILE
Definition variables.h:182
IC_MODE_ZERO
@ IC_MODE_ZERO
Definition variables.h:178
SimCtx::InitialConstantContra
Cmpnts InitialConstantContra
Definition variables.h:919
RankNeighbors::rank_zp
PetscMPIInt rank_zp
Definition variables.h:225
UserCtx::Ucont_rm1
Vec Ucont_rm1
Definition variables.h:1121
SimCtx::searchMetrics
SearchMetricsState searchMetrics
Definition variables.h:1005
SimCtx::i_periodic
PetscInt i_periodic
Definition variables.h:953
Particle::diffusivity
PetscReal diffusivity
Definition variables.h:216
UserCtx::lUcont
Vec lUcont
Definition variables.h:1113
UserCtx::Diffusivity
Vec Diffusivity
Definition variables.h:1116
SimCtx::BrownianMotionRNG
PetscRandom BrownianMotionRNG
Definition variables.h:1006
UserCtx::GridOrientation
PetscInt GridOrientation
Definition variables.h:1092
UserCtx::info
DMDALocalInfo info
Definition variables.h:1086
UserCtx::lUcat
Vec lUcat
Definition variables.h:1113
Cmpnts::y
PetscScalar y
Definition variables.h:122
IC_FIELD_UCONT
@ IC_FIELD_UCONT
Definition variables.h:188
IC_FIELD_UCAT
@ IC_FIELD_UCAT
Definition variables.h:187
SimCtx::size
PetscMPIInt size
Definition variables.h:863
UserCtx::IM
PetscInt IM
Definition variables.h:1088
SimCtx::_io_context_buffer
char _io_context_buffer[PETSC_MAX_PATH_LEN]
Definition variables.h:884
Particle::weights
Cmpnts weights
Definition variables.h:213
NUM_FACES
@ NUM_FACES
Definition variables.h:173
SimCtx::les
PetscInt les
Active LES closure; an LESModelType value.
Definition variables.h:984
UserCtx::lDiffusivity
Vec lDiffusivity
Definition variables.h:1116
BoundaryFaceConfig::mathematical_type
BCType mathematical_type
Definition variables.h:392
SearchMetricsState::searchAttempts
PetscInt64 searchAttempts
Definition variables.h:265
SimCtx::initialConditionField
InitialConditionField initialConditionField
Definition variables.h:917
Particle::PID
PetscInt64 PID
Definition variables.h:209
SearchMetricsState::maxTraversalFailCount
PetscInt64 maxTraversalFailCount
Definition variables.h:272
BC_FACE_NEG_X
@ BC_FACE_NEG_X
Definition variables.h:288
BC_FACE_POS_X
@ BC_FACE_POS_X
Definition variables.h:288
BoundingBox
Defines a 3D axis-aligned bounding box.
Definition variables.h:197
Cmpnts
A 3D point or vector with PetscScalar components.
Definition variables.h:121
Particle
Defines a particle's core properties for Lagrangian tracking.
Definition variables.h:208
SimCtx
The master context for the entire simulation.
Definition variables.h:859
UserCtx
User-defined context containing data specific to a single computational grid level.
Definition variables.h:1074
walkingsearch.h
Header file for particle location functions using the walking search algorithm.
LocateParticleOrFindMigrationTarget
PetscErrorCode LocateParticleOrFindMigrationTarget(UserCtx *user, Particle *particle, ParticleLocationStatus *status_out)
Locates a particle's host cell or identifies its migration target using a robust walk search.
Definition walkingsearch.c:959
wall_function_loglaw
void wall_function_loglaw(UserCtx *user, double roughness_height, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts *velocity_boundary, PetscReal *friction_velocity, double normal_x, double normal_y, double normal_z)
Applies log-law wall function with roughness correction.
Definition wallfunction.c:711
find_utau_loglaw
double find_utau_loglaw(double velocity, double wall_distance, double roughness_length)
Solves for friction velocity using simple log-law (explicit formula)
Definition wallfunction.c:629
u_Werner
double u_Werner(double kinematic_viscosity, double wall_distance, double friction_velocity)
Computes velocity using Werner-Wengle wall function.
Definition wallfunction.c:490
wall_function
void wall_function(UserCtx *user, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts *velocity_boundary, PetscReal *friction_velocity, double normal_x, double normal_y, double normal_z)
Applies standard wall function with Werner-Wengle model.
Definition wallfunction.c:648
u_Werner_explicit
double u_Werner_explicit(double kinematic_viscosity, double wall_distance, double wall_shear_stress)
Wall-parallel velocity at a distance, from the Werner-Wengle wall stress.
Definition wallfunction.c:557
find_utau_hydset
double find_utau_hydset(double kinematic_viscosity, double known_velocity, double wall_distance, double initial_guess, double roughness_height)
Solves for friction velocity using Newton-Raphson iteration.
Definition wallfunction.c:256
u_Cabot
double u_Cabot(double kinematic_viscosity, double wall_distance, double friction_velocity, double pressure_gradient_tangent, double wall_shear_stress)
Computes velocity using Cabot wall function.
Definition wallfunction.c:386
u_loglaw
double u_loglaw(double wall_distance, double friction_velocity, double roughness_length)
Computes velocity using simple log-law (smooth wall with roughness offset)
Definition wallfunction.c:620
wall_function_Cabot
void wall_function_Cabot(UserCtx *user, double roughness_height, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts *velocity_boundary, PetscReal *friction_velocity, double normal_x, double normal_y, double normal_z, double pressure_gradient_x, double pressure_gradient_y, double pressure_gradient_z, int iteration_count)
Applies Cabot non-equilibrium wall function with pressure gradients.
Definition wallfunction.c:773
noslip
void noslip(UserCtx *user, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts *velocity_boundary, double normal_x, double normal_y, double normal_z)
Applies no-slip wall boundary condition with linear interpolation.
Definition wallfunction.c:68
find_utau_Cabot
void find_utau_Cabot(double kinematic_viscosity, double velocity, double wall_distance, double initial_guess, double pressure_gradient_tangent, double pressure_gradient_normal, double *friction_velocity, double *wall_shear_velocity, double *wall_shear_normal)
Solves for friction velocity using Cabot wall function.
Definition wallfunction.c:440
integrate_1
double integrate_1(double kinematic_viscosity, double wall_distance, double friction_velocity, int integration_mode)
Integrates eddy viscosity profile from wall to distance y.
Definition wallfunction.c:318
f_hydset
double f_hydset(double kinematic_viscosity, double known_velocity, double wall_distance, double friction_velocity_guess, double roughness_height)
Residual function for friction velocity equation (log-law with roughness)
Definition wallfunction.c:209
E_coeff
double E_coeff(double friction_velocity, double roughness_height, double kinematic_viscosity)
Computes roughness-modified log-law coefficient E.
Definition wallfunction.c:140
taw_Werner
double taw_Werner(double kinematic_viscosity, double velocity, double wall_distance)
Wall shear stress from the cell-integrated Werner-Wengle relation.
Definition wallfunction.c:523
nu_t
double nu_t(double yplus)
Computes turbulent eddy viscosity ratio (ν_t / ν)
Definition wallfunction.c:303
find_utau_Werner
double find_utau_Werner(double kinematic_viscosity, double velocity, double wall_distance, double initial_guess)
Friction velocity from a point velocity, inverting u_Werner() exactly.
Definition wallfunction.c:591
freeslip
void freeslip(UserCtx *user, double distance_reference, double distance_boundary, Cmpnts velocity_wall, Cmpnts velocity_reference, Cmpnts *velocity_boundary, double normal_x, double normal_y, double normal_z)
Applies free-slip wall boundary condition.
Definition wallfunction.c:101
Generated by doxygen 1.9.8