PICurv 0.1.0
A Parallel Particle-In-Cell Solver for Curvilinear LES
Loading...
Searching...
No Matches
logging.c
Go to the documentation of this file.
1// logging.c
2#include "logging.h"
3#include "statistics_window.h"
4#include "statistics_accumulator.h"
5#include "AnalyticalSolutions.h"
6#include "particle_field_catalog.h"
7#include "verification_sources.h"
8
9/* Maximum temporary buffer size for converting numbers to strings */
10#define TMP_BUF_SIZE 128
11
12// --------------------- Static Variable for Log Level ---------------------
13
14/**
15 * @brief Static variable to cache the current logging level.
16 *
17 * Initialized to -1 to indicate that the log level has not been set yet.
18 */
19static LogLevel current_log_level = -1;
20
21// --------------------- Static Variables for Allow-List -------------------
22
23/**
24 * @brief Global/static array of function names allowed to log.
25 */
26static char** gAllowedFunctions = NULL;
27
28/**
29 * @brief Number of entries in the gAllowedFunctions array.
30 */
31static int gNumAllowed = 0;
32
33enum {
34 SEARCH_METRIC_SUM_SEARCH_ATTEMPTS = 0,
35 SEARCH_METRIC_SUM_SEARCH_POPULATION,
36 SEARCH_METRIC_SUM_SEARCH_LOCATED,
37 SEARCH_METRIC_SUM_SEARCH_LOST,
38 SEARCH_METRIC_SUM_TRAVERSAL_STEPS,
39 SEARCH_METRIC_SUM_RESEARCH,
40 SEARCH_METRIC_SUM_TIE_BREAKS,
41 SEARCH_METRIC_SUM_BOUNDARY_CLAMPS,
42 SEARCH_METRIC_SUM_BBOX_GUESS_SUCCESS,
43 SEARCH_METRIC_SUM_BBOX_GUESS_FALLBACK,
44 SEARCH_METRIC_SUM_MAX_TRAVERSAL_FAILS,
45 SEARCH_METRIC_MAX_TRAVERSAL_STEPS,
46 SEARCH_METRIC_MAX_PASS_DEPTH,
47 SEARCH_METRIC_REDUCTION_LEN
48};
49
50/**
51 * @brief Internal reduction callback for packed search metrics.
52 * @details Local to this translation unit.
53 */
54static void SearchMetricsReduceOp(void *invec, void *inoutvec, int *len, MPI_Datatype *datatype)
55{
56 PetscReal *in = (PetscReal *)invec;
57 PetscReal *inout = (PetscReal *)inoutvec;
58 (void)datatype;
59
60 for (int idx = 0; idx < *len; idx += SEARCH_METRIC_REDUCTION_LEN) {
61 inout[idx + SEARCH_METRIC_SUM_SEARCH_ATTEMPTS] += in[idx + SEARCH_METRIC_SUM_SEARCH_ATTEMPTS];
62 inout[idx + SEARCH_METRIC_SUM_SEARCH_POPULATION] += in[idx + SEARCH_METRIC_SUM_SEARCH_POPULATION];
63 inout[idx + SEARCH_METRIC_SUM_SEARCH_LOCATED] += in[idx + SEARCH_METRIC_SUM_SEARCH_LOCATED];
64 inout[idx + SEARCH_METRIC_SUM_SEARCH_LOST] += in[idx + SEARCH_METRIC_SUM_SEARCH_LOST];
65 inout[idx + SEARCH_METRIC_SUM_TRAVERSAL_STEPS] += in[idx + SEARCH_METRIC_SUM_TRAVERSAL_STEPS];
66 inout[idx + SEARCH_METRIC_SUM_RESEARCH] += in[idx + SEARCH_METRIC_SUM_RESEARCH];
67 inout[idx + SEARCH_METRIC_SUM_TIE_BREAKS] += in[idx + SEARCH_METRIC_SUM_TIE_BREAKS];
68 inout[idx + SEARCH_METRIC_SUM_BOUNDARY_CLAMPS] += in[idx + SEARCH_METRIC_SUM_BOUNDARY_CLAMPS];
69 inout[idx + SEARCH_METRIC_SUM_BBOX_GUESS_SUCCESS] += in[idx + SEARCH_METRIC_SUM_BBOX_GUESS_SUCCESS];
70 inout[idx + SEARCH_METRIC_SUM_BBOX_GUESS_FALLBACK] += in[idx + SEARCH_METRIC_SUM_BBOX_GUESS_FALLBACK];
71 inout[idx + SEARCH_METRIC_SUM_MAX_TRAVERSAL_FAILS] += in[idx + SEARCH_METRIC_SUM_MAX_TRAVERSAL_FAILS];
72 inout[idx + SEARCH_METRIC_MAX_TRAVERSAL_STEPS] = PetscMax(inout[idx + SEARCH_METRIC_MAX_TRAVERSAL_STEPS],
73 in[idx + SEARCH_METRIC_MAX_TRAVERSAL_STEPS]);
74 inout[idx + SEARCH_METRIC_MAX_PASS_DEPTH] = PetscMax(inout[idx + SEARCH_METRIC_MAX_PASS_DEPTH],
75 in[idx + SEARCH_METRIC_MAX_PASS_DEPTH]);
76 }
77}
78
79// --------------------- Function Implementations ---------------------
80
81/**
82 * @brief Implementation of \ref get_log_level().
83 * @details Full API contract (arguments, ownership, side effects) is documented with
84 * the header declaration in `include/logging.h`.
85 * @see get_log_level()
86 */
87LogLevel get_log_level() {
88 if (current_log_level == -1) { // Log level not set yet
89 const char *env = getenv("LOG_LEVEL");
90 if (!env) {
91 current_log_level = LOG_ERROR; // Default level
92 }
93 else if (strcmp(env, "DEBUG") == 0) {
94 current_log_level = LOG_DEBUG;
95 }
96 else if (strcmp(env, "INFO") == 0) {
97 current_log_level = LOG_INFO;
98 }
99 else if (strcmp(env, "WARNING") == 0) {
100 current_log_level = LOG_WARNING;
101 }
102 else if (strcmp(env, "VERBOSE") == 0) {
103 current_log_level = LOG_VERBOSE;
104 }
105 else if (strcmp(env, "TRACE") == 0) {
106 current_log_level = LOG_TRACE;
107 }
108 else {
109 current_log_level = LOG_ERROR; // Default if unrecognized
110 }
111 }
112 return current_log_level;
113}
114
115/**
116 * @brief Internal helper implementation: `print_log_level()`.
117 * @details Local to this translation unit.
118 */
119PetscErrorCode print_log_level(void)
120{
121 PetscMPIInt rank;
122 PetscErrorCode ierr;
123 int level;
124 const char *level_name;
125
126 PetscFunctionBeginUser;
127 /* get MPI rank */
128 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRMPI(ierr);
129
130 /* decide level name */
131 level = get_log_level();
132 level_name = (level == LOG_ERROR) ? "ERROR" :
133 (level == LOG_WARNING) ? "WARNING" :
134 (level == LOG_INFO) ? "INFO" :
135 (level == LOG_DEBUG) ? "DEBUG" :
136 (level == LOG_VERBOSE) ? "VERBOSE" :
137 (level == LOG_TRACE) ? "TRACE" :
138 "UNKNOWN";
139
140 /* print it out */
141 ierr = PetscPrintf(PETSC_COMM_SELF,
142 "Current log level: %s (%d) | rank: %d\n",
143 level_name, level, (int)rank);
144 CHKERRMPI(ierr);
145
146 PetscFunctionReturn(PETSC_SUCCESS);
147}
148
149/**
150 * @brief Implementation of \ref set_allowed_functions().
151 * @details Full API contract (arguments, ownership, side effects) is documented with
152 * the header declaration in `include/logging.h`.
153 * @see set_allowed_functions()
154 */
155void set_allowed_functions(const char** functionList, int count)
156{
157 // 1. Free any existing entries
158 if (gAllowedFunctions) {
159 for (int i = 0; i < gNumAllowed; ++i) {
160 free(gAllowedFunctions[i]); // each was strdup'ed
161 }
162 free(gAllowedFunctions);
163 gAllowedFunctions = NULL;
164 gNumAllowed = 0;
165 }
166
167 // 2. Allocate new array
168 if (count > 0) {
169 gAllowedFunctions = (char**)malloc(sizeof(char*) * count);
170 }
171
172 // 3. Copy the new entries
173 for (int i = 0; i < count; ++i) {
174 // strdup is a POSIX function. If not available, implement your own string copy.
175 gAllowedFunctions[i] = strdup(functionList[i]);
176 }
177 gNumAllowed = count;
178}
179
180/**
181 * @brief Implementation of \ref is_function_allowed().
182 * @details Full API contract (arguments, ownership, side effects) is documented with
183 * the header declaration in `include/logging.h`.
184 * @see is_function_allowed()
185 */
186PetscBool is_function_allowed(const char* functionName)
187{
188 /* no list ⇒ allow all */
189 if (gNumAllowed == 0) {
190 return PETSC_TRUE;
191 }
192
193 /* otherwise only the listed functions are allowed */
194 for (int i = 0; i < gNumAllowed; ++i) {
195 if (strcmp(gAllowedFunctions[i], functionName) == 0) {
196 return PETSC_TRUE;
197 }
198 }
199 return PETSC_FALSE;
200}
201
202/**
203 * @brief Implementation of \ref LOG_CELL_VERTICES().
204 * @details Full API contract (arguments, ownership, side effects) is documented with
205 * the header declaration in `include/logging.h`.
206 * @see LOG_CELL_VERTICES()
207 */
208PetscErrorCode LOG_CELL_VERTICES(const Cell *cell, PetscMPIInt rank)
209{
210
211 // Validate input pointers
212 if (cell == NULL) {
213 LOG_ALLOW(LOCAL,LOG_ERROR, "'cell' is NULL.\n");
214 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "LOG_CELL_VERTICES - Input parameter 'cell' is NULL.");
215 }
216
217 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Rank %d, Cell Vertices:\n", rank);
218 for(int i = 0; i < 8; i++){
219 LOG_ALLOW(LOCAL,LOG_VERBOSE, " Vertex[%d]: (%.2f, %.2f, %.2f)\n",
220 i, cell->vertices[i].x, cell->vertices[i].y, cell->vertices[i].z);
221 }
222
223 return 0; // Indicate successful execution
224}
225
226
227/**
228 * @brief Implementation of \ref LOG_FACE_DISTANCES().
229 * @details Full API contract (arguments, ownership, side effects) is documented with
230 * the header declaration in `include/logging.h`.
231 * @see LOG_FACE_DISTANCES()
232 */
233PetscErrorCode LOG_FACE_DISTANCES(PetscReal* d)
234{
235
236 // Validate input array
237 if (d == NULL) {
238 LOG_ALLOW(LOCAL,LOG_ERROR, " 'd' is NULL.\n");
239 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, " Input array 'd' is NULL.");
240 }
241
242 PetscPrintf(PETSC_COMM_SELF, " Face Distances:\n");
243 PetscPrintf(PETSC_COMM_SELF, " LEFT(%d): %.15f\n", LEFT, d[LEFT]);
244 PetscPrintf(PETSC_COMM_SELF, " RIGHT(%d): %.15f\n", RIGHT, d[RIGHT]);
245 PetscPrintf(PETSC_COMM_SELF, " BOTTOM(%d): %.15f\n", BOTTOM, d[BOTTOM]);
246 PetscPrintf(PETSC_COMM_SELF, " TOP(%d): %.15f\n", TOP, d[TOP]);
247 PetscPrintf(PETSC_COMM_SELF, " FRONT(%d): %.15f\n", FRONT, d[FRONT]);
248 PetscPrintf(PETSC_COMM_SELF, " BACK(%d): %.15f\n", BACK, d[BACK]);
249
250 return 0; // Indicate successful execution
251}
252
253/*
254 * Helper function: Converts an integer (of type int) to a string.
255 */
256static void IntToStr(int value, char *buf, size_t bufsize)
257{
258 snprintf(buf, bufsize, "%d", value);
259}
260
261/*
262 * Helper function: Converts a 64‐bit integer to a string.
263 */
264static void Int64ToStr(PetscInt64 value, char *buf, size_t bufsize)
265{
266 snprintf(buf, bufsize, "%ld", value);
267}
268
269/*
270 * Helper function: Converts three integers into a formatted string "(i, j, k)".
271 */
272static void CellToStr(const PetscInt *cell, char *buf, size_t bufsize)
273{
274 snprintf(buf, bufsize, "(%d, %d, %d)", cell[0], cell[1], cell[2]);
275}
276
277/*
278 * Helper function: Converts three PetscReal values into a formatted string "(x, y, z)".
279 */
280static void TripleRealToStr(const PetscReal *arr, char *buf, size_t bufsize)
281{
282 snprintf(buf, bufsize, "(%.4f, %.4f, %.4f)", arr[0], arr[1], arr[2]);
283}
284
285/*
286 * Helper function: Computes the maximum string length for each column (across all particles).
287 *
288 * The function examines every particle (from 0 to nParticles-1) and converts the value to a
289 * string using the helper functions above. The maximum length is stored in the pointers provided.
290 *
291 * @param nParticles Number of particles.
292 * @param ranks Array of particle MPI ranks.
293 * @param pids Array of particle IDs.
294 * @param cellIDs Array of cell IDs (stored consecutively, 3 per particle).
295 * @param positions Array of positions (3 per particle).
296 * @param velocities Array of velocities (3 per particle).
297 * @param weights Array of weights (3 per particle).
298 * @param wRank [out] Maximum width for Rank column.
299 * @param wPID [out] Maximum width for PID column.
300 * @param wCell [out] Maximum width for Cell column.
301 * @param wPos [out] Maximum width for Position column.
302 * @param wVel [out] Maximum width for Velocity column.
303 * @param wWt [out] Maximum width for Weights column.
304 */
305static PetscErrorCode ComputeMaxColumnWidths(PetscInt nParticles,
306 const PetscMPIInt *ranks,
307 const PetscInt64 *pids,
308 const PetscInt *cellIDs,
309 const PetscReal *positions,
310 const PetscReal *velocities,
311 const PetscReal *weights,
312 int *wRank, int *wPID, int *wCell,
313 int *wPos, int *wVel, int *wWt)
314{
315 char tmp[TMP_BUF_SIZE];
316
317 *wRank = strlen("Rank"); /* Start with the header label lengths */
318 *wPID = strlen("PID");
319 *wCell = strlen("Cell (i,j,k)");
320 *wPos = strlen("Position (x,y,z)");
321 *wVel = strlen("Velocity (x,y,z)");
322 *wWt = strlen("Weights (a1,a2,a3)");
323
324 for (PetscInt i = 0; i < nParticles; i++) {
325 /* Rank */
326 IntToStr(ranks[i], tmp, TMP_BUF_SIZE);
327 if ((int)strlen(tmp) > *wRank) *wRank = (int)strlen(tmp);
328
329 /* PID */
330 Int64ToStr(pids[i], tmp, TMP_BUF_SIZE);
331 if ((int)strlen(tmp) > *wPID) *wPID = (int)strlen(tmp);
332
333 /* Cell: use the three consecutive values */
334 CellToStr(&cellIDs[3 * i], tmp, TMP_BUF_SIZE);
335 if ((int)strlen(tmp) > *wCell) *wCell = (int)strlen(tmp);
336
337 /* Position */
338 TripleRealToStr(&positions[3 * i], tmp, TMP_BUF_SIZE);
339 if ((int)strlen(tmp) > *wPos) *wPos = (int)strlen(tmp);
340
341 /* Velocity */
342 TripleRealToStr(&velocities[3 * i], tmp, TMP_BUF_SIZE);
343 if ((int)strlen(tmp) > *wVel) *wVel = (int)strlen(tmp);
344
345 /* Weights */
346 TripleRealToStr(&weights[3 * i], tmp, TMP_BUF_SIZE);
347 if ((int)strlen(tmp) > *wWt) *wWt = (int)strlen(tmp);
348 }
349 return 0;
350}
351
352/*
353 * Helper function: Builds a format string for a table row.
354 *
355 * The format string will include proper width specifiers for each column.
356 * For example, it might create something like:
357 *
358 * "| %-6s | %-8s | %-20s | %-25s | %-25s | %-25s |\n"
359 *
360 * @param wRank Maximum width for the Rank column.
361 * @param wPID Maximum width for the PID column.
362 * @param wCell Maximum width for the Cell column.
363 * @param wPos Maximum width for the Position column.
364 * @param wVel Maximum width for the Velocity column.
365 * @param wWt Maximum width for the Weights column.
366 * @param fmtStr Buffer in which to build the format string.
367 * @param bufSize Size of fmtStr.
368 */
369static void BuildRowFormatString(PetscMPIInt wRank, PetscInt wPID, PetscInt wCell, PetscInt wPos, PetscInt wVel, PetscInt wWt, char *fmtStr, size_t bufSize)
370{
371 // Build a format string using snprintf.
372 // We assume that the Rank is an int (%d), PID is a 64-bit int (%ld)
373 // and the remaining columns are strings (which have been formatted already).
374 snprintf(fmtStr, bufSize,
375 "| %%-%dd | %%-%dd | %%-%ds | %%-%ds | %%-%ds | %%-%ds |\n",
376 wRank, wPID, wCell, wPos, wVel, wWt);
377}
378
379/*
380 * Helper function: Builds a header string for the table using column titles.
381 */
382static void BuildHeaderString(char *headerStr, size_t bufSize, PetscMPIInt wRank, PetscInt wPID, PetscInt wCell, PetscInt wPos, PetscInt wVel, PetscInt wWt)
383{
384 snprintf(headerStr, bufSize,
385 "| %-*s | %-*s | %-*s | %-*s | %-*s | %-*s |\n",
386 (int)wRank, "Rank",
387 (int)wPID, "PID",
388 (int)wCell, "Cell (i,j,k)",
389 (int)wPos, "Position (x,y,z)",
390 (int)wVel, "Velocity (x,y,z)",
391 (int)wWt, "Weights (a1,a2,a3)");
392}
393
394/**
395 * @brief Implementation of \ref LOG_PARTICLE_FIELDS().
396 * @details Full API contract (arguments, ownership, side effects) is documented with
397 * the header declaration in `include/logging.h`.
398 * @see LOG_PARTICLE_FIELDS()
399 */
400PetscErrorCode LOG_PARTICLE_FIELDS(UserCtx* user, PetscInt printInterval)
401{
402 DM swarm = user->swarm;
403 PetscErrorCode ierr;
404 PetscInt localNumParticles;
405 PetscReal *positions = NULL;
406 PetscInt64 *particleIDs = NULL;
407 PetscMPIInt *particleRanks = NULL;
408 PetscInt *cellIDs = NULL;
409 PetscReal *weights = NULL;
410 PetscReal *velocities = NULL;
411 PetscMPIInt rank;
412
413 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
414 LOG_ALLOW(LOCAL,LOG_INFO, "Rank %d is retrieving particle data.\n", rank);
415
416 ierr = DMSwarmGetLocalSize(swarm, &localNumParticles); CHKERRQ(ierr);
417 LOG_ALLOW(LOCAL,LOG_DEBUG,"Rank %d has %d particles.\n", rank, localNumParticles);
418
419 ierr = DMSwarmGetField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_POSITION), NULL, NULL, (void**)&positions); CHKERRQ(ierr);
420 ierr = DMSwarmGetField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_PID), NULL, NULL, (void**)&particleIDs); CHKERRQ(ierr);
421 ierr = DMSwarmGetField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_RANK), NULL, NULL, (void**)&particleRanks); CHKERRQ(ierr);
422 ierr = DMSwarmGetField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_CELL_ID), NULL, NULL, (void**)&cellIDs); CHKERRQ(ierr);
423 ierr = DMSwarmGetField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_WEIGHT), NULL, NULL, (void**)&weights); CHKERRQ(ierr);
424 ierr = DMSwarmGetField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_VELOCITY), NULL, NULL, (void**)&velocities); CHKERRQ(ierr);
425
426 /* Compute maximum column widths. */
427 int wRank, wPID, wCell, wPos, wVel, wWt;
428 wRank = wPID = wCell = wPos = wVel = wWt = 0;
429 ierr = ComputeMaxColumnWidths(localNumParticles, particleRanks, particleIDs, cellIDs,
430 positions, velocities, weights,
431 &wRank, &wPID, &wCell, &wPos, &wVel, &wWt); CHKERRQ(ierr);
432
433 /* Build a header string and a row format string. */
434 char headerFmt[256];
435 char rowFmt[256];
436 BuildHeaderString(headerFmt, sizeof(headerFmt), wRank, wPID, wCell, wPos, wVel, wWt);
437 BuildRowFormatString(wRank, wPID, wCell, wPos, wVel, wWt, rowFmt, sizeof(rowFmt));
438
439 /* Print header (using synchronized printing for parallel output). */
440 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "--------------------------------------------------------------------------------------------------------------\n"); CHKERRQ(ierr);
441 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "%s", headerFmt); CHKERRQ(ierr);
442 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "--------------------------------------------------------------------------------------------------------------\n"); CHKERRQ(ierr);
443
444 /* Loop over particles and print every printInterval-th row. */
445 char rowStr[256];
446 for (PetscInt i = 0; i < localNumParticles; i++) {
447 if (i % printInterval == 0) {
448 // ------- DEBUG
449 //char cellStr[TMP_BUF_SIZE], posStr[TMP_BUF_SIZE], velStr[TMP_BUF_SIZE], wtStr[TMP_BUF_SIZE];
450 //CellToStr(&cellIDs[3*i], cellStr, TMP_BUF_SIZE);
451 //TripleRealToStr(&positions[3*i], posStr, TMP_BUF_SIZE);
452 //TripleRealToStr(&velocities[3*i], velStr, TMP_BUF_SIZE);
453 // TripleRealToStr(&weights[3*i], wtStr, TMP_BUF_SIZE);
454
455 // if (rank == 0) { // Or whatever rank is Rank 0
456 //PetscPrintf(PETSC_COMM_SELF, "[Rank 0 DEBUG LPF] Particle %lld: PID=%lld, Rank=%d\n", (long long)i, (long long)particleIDs[i], particleRanks[i]);
457 //PetscPrintf(PETSC_COMM_SELF, "[Rank 0 DEBUG LPF] Raw Pos: (%.10e, %.10e, %.10e)\n", positions[3*i+0], positions[3*i+1], positions[3*i+2]);
458 //PetscPrintf(PETSC_COMM_SELF, "[Rank 0 DEBUG LPF] Str Pos: %s\n", posStr);
459 //PetscPrintf(PETSC_COMM_SELF, "[Rank 0 DEBUG LPF] Raw Vel: (%.10e, %.10e, %.10e)\n", velocities[3*i+0], velocities[3*i+1], velocities[3*i+2]);
460 // PetscPrintf(PETSC_COMM_SELF, "[Rank 0 DEBUG LPF] Str Vel: %s\n", velStr);
461 // Add similar for cell, weights
462 // PetscPrintf(PETSC_COMM_SELF, "[Rank 0 DEBUG LPF] About to build rowStr for particle %lld\n", (long long)i);
463 // fflush(stdout);
464 // }
465
466 // snprintf(rowStr, sizeof(rowStr), rowFmt,
467 // particleRanks[i],
468 // particleIDs[i],
469 // cellStr,
470 // posStr,
471 // velStr,
472 // wtStr);
473
474
475 // ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "%s", rowStr); CHKERRQ(ierr);
476
477 // ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "%s", rowStr); CHKERRQ(ierr);
478
479 // -------- DEBUG
480 /* Format the row by converting each field to a string first.
481 * We use temporary buffers and then build the row string.
482 */
483
484 char cellStr[TMP_BUF_SIZE], posStr[TMP_BUF_SIZE], velStr[TMP_BUF_SIZE], wtStr[TMP_BUF_SIZE];
485 CellToStr(&cellIDs[3*i], cellStr, TMP_BUF_SIZE);
486 TripleRealToStr(&positions[3*i], posStr, TMP_BUF_SIZE);
487 TripleRealToStr(&velocities[3*i], velStr, TMP_BUF_SIZE);
488 TripleRealToStr(&weights[3*i], wtStr, TMP_BUF_SIZE);
489
490 /* Build the row string. Note that for the integer fields we can use the row format string. */
491 snprintf(rowStr, sizeof(rowStr), rowFmt,
492 particleRanks[i],
493 particleIDs[i],
494 cellStr,
495 posStr,
496 velStr,
497 wtStr);
498 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "%s", rowStr); CHKERRQ(ierr);
499 }
500 }
501
502
503 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "--------------------------------------------------------------------------------------------------------------\n"); CHKERRQ(ierr);
504 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "\n"); CHKERRQ(ierr);
505 ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT); CHKERRQ(ierr);
506
507 LOG_ALLOW_SYNC(GLOBAL,LOG_DEBUG,"Completed printing on Rank %d.\n", rank);
508
509 /* Restore fields */
510 ierr = DMSwarmRestoreField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_POSITION), NULL, NULL, (void**)&positions); CHKERRQ(ierr);
511 ierr = DMSwarmRestoreField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_PID), NULL, NULL, (void**)&particleIDs); CHKERRQ(ierr);
512 ierr = DMSwarmRestoreField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_RANK), NULL, NULL, (void**)&particleRanks); CHKERRQ(ierr);
513 ierr = DMSwarmRestoreField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_CELL_ID), NULL, NULL, (void**)&cellIDs); CHKERRQ(ierr);
514 ierr = DMSwarmRestoreField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_WEIGHT), NULL, NULL, (void**)&weights); CHKERRQ(ierr);
515 ierr = DMSwarmRestoreField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_VELOCITY), NULL, NULL, (void**)&velocities); CHKERRQ(ierr);
516
517 LOG_ALLOW(LOCAL,LOG_DEBUG, "Restored all particle fields.\n");
518 return 0;
519}
520
521/**
522 * @brief Implementation of \ref IsParticleConsoleSnapshotEnabled().
523 * @details Full API contract (arguments, ownership, side effects) is documented with
524 * the header declaration in `include/logging.h`.
525 * @see IsParticleConsoleSnapshotEnabled()
526 */
527
528PetscBool IsParticleConsoleSnapshotEnabled(const SimCtx *simCtx)
529{
530 if (!simCtx) {
531 return PETSC_FALSE;
532 }
533 return (PetscBool)(simCtx->np > 0 &&
534 simCtx->particleConsoleOutputFreq > 0 &&
535 get_log_level() >= LOG_INFO);
536}
537
538/**
539 * @brief Implementation of \ref ShouldEmitPeriodicParticleConsoleSnapshot().
540 * @details Full API contract (arguments, ownership, side effects) is documented with
541 * the header declaration in `include/logging.h`.
542 * @see ShouldEmitPeriodicParticleConsoleSnapshot()
543 */
544
545PetscBool ShouldEmitPeriodicParticleConsoleSnapshot(const SimCtx *simCtx, PetscInt completed_step)
546{
547 return (PetscBool)(IsParticleConsoleSnapshotEnabled(simCtx) &&
548 completed_step > 0 &&
549 completed_step % simCtx->particleConsoleOutputFreq == 0);
550}
551
552/**
553 * @brief Implementation of \ref EmitParticleConsoleSnapshot().
554 * @details Full API contract (arguments, ownership, side effects) is documented with
555 * the header declaration in `include/logging.h`.
556 * @see EmitParticleConsoleSnapshot()
557 */
558
559PetscErrorCode EmitParticleConsoleSnapshot(UserCtx *user, SimCtx *simCtx, PetscInt step)
560{
561 PetscErrorCode ierr;
562
563 PetscFunctionBeginUser;
564 LOG(GLOBAL, LOG_INFO, "Particle states at step %d:\n", step);
565 ierr = LOG_PARTICLE_FIELDS(user, simCtx->LoggingFrequency); CHKERRQ(ierr);
566 PetscFunctionReturn(0);
567}
568
569
570/**
571 * @brief Remove leading and trailing whitespace from a mutable configuration string.
572 */
573static void trim(char *s)
574{
575 if (!s) return;
576
577 /* ---- 1. strip leading blanks ----------------------------------- */
578 char *p = s;
579 while (*p && isspace((unsigned char)*p))
580 ++p;
581
582 if (p != s) /* move the trimmed text forward */
583 memmove(s, p, strlen(p) + 1); /* +1 to copy the final NUL */
584
585 /* ---- 2. strip trailing blanks ---------------------------------- */
586 size_t len = strlen(s);
587 while (len > 0 && isspace((unsigned char)s[len - 1]))
588 s[--len] = '\0';
589}
590
591/* ------------------------------------------------------------------------- */
592/**
593 * @brief Implementation of \ref LoadAllowedFunctionsFromFile().
594 * @details Full API contract (arguments, ownership, side effects) is documented with
595 * the header declaration in `include/logging.h`.
596 * @see LoadAllowedFunctionsFromFile()
597 */
598PetscErrorCode LoadAllowedFunctionsFromFile(const char filename[],
599 char ***funcsOut,
600 PetscInt *nOut)
601{
602 FILE *fp = NULL;
603 char **funcs = NULL;
604 size_t cap = 16; /* initial capacity */
605 size_t n = 0; /* number of names */
606 char line[PETSC_MAX_PATH_LEN];
607 PetscErrorCode ierr;
608
609 PetscFunctionBegin;
610
611 /* ---------------------------------------------------------------------- */
612 /* 1. Open file */
613 fp = fopen(filename, "r");
614 if (!fp) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN,
615 "Cannot open %s", filename);
616
617 /* 2. Allocate initial pointer array */
618 ierr = PetscMalloc1(cap, &funcs); CHKERRQ(ierr);
619
620 /* 3. Read file line by line */
621 while (fgets(line, sizeof line, fp)) {
622 /* Strip everything after a comment character '#'. */
623 char *hash = strchr(line, '#');
624 if (hash) *hash = '\0';
625
626 trim(line); /* remove leading/trailing blanks */
627 if (!*line) continue; /* skip if empty */
628
629 /* Grow the array if necessary */
630 if (n == cap) {
631 cap *= 2;
632 ierr = PetscRealloc(cap * sizeof(*funcs), (void **)&funcs); CHKERRQ(ierr);
633 }
634
635 /* Deep‑copy the cleaned identifier */
636 ierr = PetscStrallocpy(line, &funcs[n++]); CHKERRQ(ierr);
637 }
638 fclose(fp);
639
640 /* 4. Return results to caller */
641 *funcsOut = funcs;
642 *nOut = (PetscInt)n;
643
644 PetscFunctionReturn(0);
645}
646
647/* ------------------------------------------------------------------------- */
648/**
649 * @brief Internal helper implementation: `FreeAllowedFunctions()`.
650 * @details Local to this translation unit.
651 */
652PetscErrorCode FreeAllowedFunctions(char **funcs, PetscInt n)
653{
654 PetscErrorCode ierr;
655 PetscFunctionBegin;
656 if (funcs) {
657 for (PetscInt i = 0; i < n; ++i) {
658 ierr = PetscFree(funcs[i]); CHKERRQ(ierr);
659 }
660 ierr = PetscFree(funcs); CHKERRQ(ierr);
661 }
662 PetscFunctionReturn(0);
663}
664
665/**
666 * @brief Implementation of \ref BCFaceToString().
667 * @details Full API contract (arguments, ownership, side effects) is documented with
668 * the header declaration in `include/logging.h`.
669 * @see BCFaceToString()
670 */
671const char* BCFaceToString(BCFace face) {
672 switch (face) {
673 case BC_FACE_NEG_X: return "-Xi (I-Min)";
674 case BC_FACE_POS_X: return "+Xi (I-Max)";
675 case BC_FACE_NEG_Y: return "-Eta (J-Min)";
676 case BC_FACE_POS_Y: return "+Eta (J-Max)";
677 case BC_FACE_NEG_Z: return "-Zeta (K-Min)";
678 case BC_FACE_POS_Z: return "+Zeta (K-Max)";
679 default: return "Unknown Face";
680 }
681}
682
683/**
684 * @brief Implementation of \ref InitialConditionModeToString().
685 * @details Full API contract (arguments, ownership, side effects) is documented with
686 * the header declaration in `include/logging.h`.
687 * @see InitialConditionModeToString()
688 */
689const char* InitialConditionModeToString(InitialConditionMode mode)
690{
691 switch(mode){
692 case IC_MODE_ZERO: return "Zero";
693 case IC_MODE_CONSTANT_CARTESIAN: return "Cartesian Constant";
694 case IC_MODE_POISEUILLE: return "Poiseuille";
695 case IC_MODE_CONSTANT_STREAMWISE: return "Streamwise Constant";
696 case IC_MODE_FILE: return "File";
697 default: return "Unknown Initial Condition";
698 }
699}
700
701/*
702 * Converts a FlowDirection enum value to its token string. The public header
703 * owns the rendered API contract.
704 */
705const char* FlowDirectionToString(FlowDirection fd)
706{
707 switch ((int)fd) {
708 case FLOW_DIR_POS_XI: return "+Xi";
709 case FLOW_DIR_NEG_XI: return "-Xi";
710 case FLOW_DIR_POS_ETA: return "+Eta";
711 case FLOW_DIR_NEG_ETA: return "-Eta";
712 case FLOW_DIR_POS_ZETA: return "+Zeta";
713 case FLOW_DIR_NEG_ZETA: return "-Zeta";
714 default: return "from INLET";
715 }
716}
717
718/**
719 * @brief Implementation of \ref ParticleInitializationToString().
720 * @details Full API contract (arguments, ownership, side effects) is documented with
721 * the header declaration in `include/logging.h`.
722 * @see ParticleInitializationToString()
723 */
724const char* ParticleInitializationToString(ParticleInitializationType ParticleInitialization)
725{
726 switch(ParticleInitialization){
727 case PARTICLE_INIT_SURFACE_RANDOM: return "Surface: Random";
728 case PARTICLE_INIT_VOLUME: return "Volume";
729 case PARTICLE_INIT_POINT_SOURCE: return "Point Source";
730 case PARTICLE_INIT_SURFACE_EDGES: return "Surface: At edges";
731 default: return "Unknown Particle Initialization";
732 }
733}
734
735/**
736 * @brief Implementation of \ref LESModelToString().
737 * @details Full API contract (arguments, ownership, side effects) is documented with
738 * the header declaration in `include/logging.h`.
739 * @see LESModelToString()
740 */
741const char* LESModelToString(LESModelType LESFlag)
742{
743 switch(LESFlag){
744 case NO_LES_MODEL: return "No LES";
745 case CONSTANT_SMAGORINSKY: return "Constant Smagorinsky";
746 case DYNAMIC_SMAGORINSKY: return "Dynamic Smagorinsky";
747 default: return "Unknown LES Flag";
748 }
749}
750
751/**
752 * @brief Implementation of \ref MomentumSolverTypeToString().
753 * @details Full API contract (arguments, ownership, side effects) is documented with
754 * the header declaration in `include/logging.h`.
755 * @see MomentumSolverTypeToString()
756 */
757const char* MomentumSolverTypeToString(MomentumSolverType SolverFlag)
758{
759 switch(SolverFlag){
760 case MOMENTUM_SOLVER_EXPLICIT_RK: return "Explicit 4 stage Runge-Kutta ";
761 case MOMENTUM_SOLVER_DUALTIME_PICARD_JAMESON_RK: return "Dual Time Picard with 4-stage Jameson RK Smoothing";
762 case MOMENTUM_SOLVER_NEWTON_KRYLOV: return "Newton Krylov";
763 default: return "Unknown Momentum Solver Type";
764 }
765}
766
767/**
768 * @brief Implementation of \ref BCTypeToString().
769 * @details Full API contract (arguments, ownership, side effects) is documented with
770 * the header declaration in `include/logging.h`.
771 * @see BCTypeToString()
772 */
773const char* BCTypeToString(BCType type) {
774 switch (type) {
775 // case DIRICHLET: return "DIRICHLET";
776 // case NEUMANN: return "NEUMANN";
777 case WALL: return "WALL";
778 case INLET: return "INLET";
779 case OUTLET: return "OUTLET";
780 case FARFIELD: return "FARFIELD";
781 case PERIODIC: return "PERIODIC";
782 case INTERFACE: return "INTERFACE";
783
784 // case CUSTOM: return "CUSTOM";
785 default: return "Unknown BC Type";
786 }
787}
788
789/**
790 * @brief Internal helper implementation: `BCHandlerTypeToString()`.
791 * @details Local to this translation unit.
792 */
793const char* BCHandlerTypeToString(BCHandlerType handler_type) {
794 switch (handler_type) {
795 // Wall & Symmetry Handlers
796 case BC_HANDLER_WALL_NOSLIP: return "noslip";
797 case BC_HANDLER_WALL_MOVING: return "moving";
798 case BC_HANDLER_SYMMETRY_PLANE: return "symmetry_plane";
799
800 // Inlet Handlers
801 case BC_HANDLER_INLET_CONSTANT_VELOCITY: return "constant_velocity";
802 case BC_HANDLER_INLET_PULSATILE_FLUX: return "pulsatile_flux";
803 case BC_HANDLER_INLET_PARABOLIC: return "parabolic";
804 case BC_HANDLER_INLET_PROFILE_FROM_FILE: return "prescribed_flow";
805
806 // Outlet Handlers
807 case BC_HANDLER_OUTLET_CONSERVATION: return "conservation";
808 case BC_HANDLER_OUTLET_PRESSURE: return "pressure";
809
810 // Other Physical Handlers
811 case BC_HANDLER_FARFIELD_NONREFLECTING: return "nonreflecting";
812
813 // Multi-Block / Interface Handlers
814 case BC_HANDLER_PERIODIC_GEOMETRIC: return "geometric";
815 case BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX: return "constant flux";
816 case BC_HANDLER_PERIODIC_DRIVEN_INITIAL_FLUX: return "initial flux";
817 case BC_HANDLER_INTERFACE_OVERSET: return "overset";
818
819 // Default case
820 case BC_HANDLER_UNDEFINED:
821 default: return "UNKNOWN_HANDLER";
822 }
823}
824
825/**
826 * @brief Implementation of \ref DualMonitorDestroy().
827 * @details Full API contract (arguments, ownership, side effects) is documented with
828 * the header declaration in `include/logging.h`.
829 * @see DualMonitorDestroy()
830 */
831PetscErrorCode DualMonitorDestroy(void **ctx)
832{
833 DualMonitorCtx *monctx = (DualMonitorCtx*)*ctx;
834 PetscErrorCode ierr;
835 PetscMPIInt rank;
836
837 PetscFunctionBeginUser;
838 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank); CHKERRQ(ierr);
839 if(!rank && monctx->file_handle){
840 fclose(monctx->file_handle);
841 }
842
843 ierr = PetscFree(monctx); CHKERRQ(ierr);
844 *ctx = NULL;
845 PetscFunctionReturn(0);
846}
847
848/**
849 * @brief A custom KSP monitor that logs the true residual to a file and optionally to the console.
850 *
851 * This function replicates the behavior of KSPMonitorTrueResidualNorm by calculating
852 * the true residual norm ||b - Ax|| itself. It unconditionally logs to a file
853 * viewer and conditionally logs to the console based on a flag in the context.
854 *
855 * @param ksp The Krylov subspace context.
856 * @param it The current iteration number.
857 * @param rnorm The preconditioned residual norm (ignored, we compute our own).
858 * @param ctx A pointer to the DualMonitorCtx structure.
859 * @return PetscErrorCode 0 on success.
860 */
861#undef __FUNCT__
862#define __FUNCT__ "DualKSPMonitor"
863/**
864 * @brief Implementation of \ref DualKSPMonitor().
865 * @details Full API contract (arguments, ownership, side effects) is documented with
866 * the header declaration in `include/logging.h`.
867 * @see DualKSPMonitor()
868 */
869
870PetscErrorCode DualKSPMonitor(KSP ksp, PetscInt it, PetscReal rnorm, void *ctx)
871{
872 DualMonitorCtx *monctx = (DualMonitorCtx*)ctx;
873 PetscErrorCode ierr;
874 PetscReal trnorm, relnorm;
875 Vec r;
876 char norm_buf[256];
877 PetscMPIInt rank;
878
879 PetscFunctionBeginUser;
880 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank); CHKERRQ(ierr);
881
882 // 1. Calculate the true residual norm.
883 ierr = KSPBuildResidual(ksp, NULL, NULL, &r); CHKERRQ(ierr);
884 ierr = VecNorm(r, NORM_2, &trnorm); CHKERRQ(ierr);
885 ierr = VecDestroy(&r); CHKERRQ(ierr);
886
887 // 2. On the first iteration, compute and store the norm of the RHS vector `b`.
888 if (it == 0) {
889 Vec b;
890 ierr = KSPGetRhs(ksp, &b); CHKERRQ(ierr);
891 ierr = VecNorm(b, NORM_2, &monctx->bnorm); CHKERRQ(ierr);
892 }
893
894 if(!rank){
895 // 3. Compute the relative norm and format the output string.
896 if (monctx->bnorm > 1.e-15) {
897 relnorm = trnorm / monctx->bnorm;
898 sprintf(norm_buf, "ts: %-5d | block: %-2d | iter: %-3d | Unprecond Norm: %12.5e | True Norm: %12.5e | Rel Norm: %12.5e",(int)monctx->step, (int)monctx->block_id, (int)it, (double)rnorm, (double)trnorm, (double)relnorm);
899 } else {
900 sprintf(norm_buf,"ts: %-5d | block: %-2d | iter: %-3d | Unprecond Norm: %12.5e | True Norm: %12.5e",(int)monctx->step, (int)monctx->block_id, (int)it, (double)rnorm, (double)trnorm);
901 }
902
903 // 4. Log to the file viewer (unconditionally).
904 if(monctx->file_handle){
905 ierr = PetscFPrintf(PETSC_COMM_SELF,monctx->file_handle,"%s\n", norm_buf); CHKERRQ(ierr);
906 }
907 // 5. Log to the console (conditionally).
908 if (monctx->log_to_console) {
909 PetscFPrintf(PETSC_COMM_SELF,stdout, "%s\n", norm_buf); CHKERRQ(ierr);
910 }
911
912 } //rank
913
914 PetscFunctionReturn(0);
915}
916
917#define SOLUTION_CONVERGENCE_FLUID_THRESHOLD 0.1
918#define SOLUTION_CONVERGENCE_REL_EPS 1.0e-30
919
931
936
937/**
938 * @brief Forms a guarded relative metric for solution-convergence logging.
939 *
940 * This helper centralizes the divide-by-nearly-zero protection used by the
941 * solution-convergence logger when turning an absolute drift into a relative
942 * one. The denominator is clamped away from zero so warmup rows, quiescent
943 * fields, and statistically small observables do not generate infinities.
944 *
945 * @param[in] numerator Absolute quantity or drift magnitude.
946 * @param[in] denominator Reference magnitude used for normalization.
947 * @return Guarded relative value `numerator / max(|denominator|, eps)`.
948 */
949static PetscReal SolutionConvergenceSafeRelative(PetscReal numerator, PetscReal denominator)
950{
951 return numerator / PetscMax(PetscAbsReal(denominator), SOLUTION_CONVERGENCE_REL_EPS);
952}
953
954/**
955 * @brief Computes instantaneous global flow observables for statistical mode.
956 *
957 * The statistical solution-convergence path does not compare full Eulerian
958 * fields. Instead, it tracks a compact history of global observables derived
959 * from the completed Eulerian state. This helper computes the current
960 * volume-weighted fluid-domain mean speed and mean kinetic energy from `Ucat`.
961 *
962 * Only physical fluid cells contribute:
963 * - solid/immersed cells are excluded using `Nvert`
964 * - cells with near-zero metric Jacobian are ignored to avoid invalid volume
965 * weights
966 *
967 * Each MPI rank accumulates local partial sums and the routine reduces them to
968 * one global pair of observables.
969 *
970 * @param[in] simCtx Simulation context owning the finest-level flow
971 * fields.
972 * @param[out] mean_speed_out Volume-weighted domain mean of `|u|`.
973 * @param[out] mean_ke_out Volume-weighted domain mean of `0.5 |u|^2`.
974 * @return PetscErrorCode 0 on success.
975 */
976static PetscErrorCode ComputeCurrentFlowObservables(SimCtx *simCtx, PetscReal *mean_speed_out, PetscReal *mean_ke_out)
977{
978 PetscReal local[3] = {0.0, 0.0, 0.0};
979 PetscReal global[3] = {0.0, 0.0, 0.0};
980 UserCtx *user = NULL;
981
982 PetscFunctionBeginUser;
983 if (!simCtx || !mean_speed_out || !mean_ke_out) {
984 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "ComputeCurrentFlowObservables received a NULL argument.");
985 }
986
987 user = simCtx->usermg.mgctx[simCtx->usermg.mglevels - 1].user;
988
989 for (PetscInt bi = 0; bi < simCtx->block_number; ++bi) {
990 const DMDALocalInfo info = user[bi].info;
991 const PetscBool x_per = (PetscBool)(simCtx->i_periodic != 0);
992 const PetscBool y_per = (PetscBool)(simCtx->j_periodic != 0);
993 const PetscBool z_per = (PetscBool)(simCtx->k_periodic != 0);
994 const PetscInt i_end = (x_per && (info.xs + info.xm == info.mx)) ? info.mx - 1 : info.xs + info.xm;
995 const PetscInt j_end = (y_per && (info.ys + info.ym == info.my)) ? info.my - 1 : info.ys + info.ym;
996 const PetscInt k_end = (z_per && (info.zs + info.zm == info.mz)) ? info.mz - 1 : info.zs + info.zm;
997 Cmpnts ***ucat = NULL;
998 PetscReal ***aj = NULL;
999 PetscReal ***nvert = NULL;
1000
1001 PetscCall(DMDAVecGetArrayRead(user[bi].fda, user[bi].Ucat, &ucat));
1002 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].Aj, &aj));
1003 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].Nvert, &nvert));
1004
1005 for (PetscInt k = info.zs; k < k_end; ++k) {
1006 for (PetscInt j = info.ys; j < j_end; ++j) {
1007 for (PetscInt i = info.xs; i < i_end; ++i) {
1008 PetscReal jac = aj[k][j][i];
1009 PetscReal cell_volume = 0.0;
1010 PetscReal speed = 0.0;
1011 PetscReal ke = 0.0;
1012
1013 if (nvert[k][j][i] > SOLUTION_CONVERGENCE_FLUID_THRESHOLD) continue;
1014 if (PetscAbsReal(jac) <= 1.0e-14) continue;
1015
1016 cell_volume = 1.0 / jac;
1017 speed = PetscSqrtReal(ucat[k][j][i].x * ucat[k][j][i].x +
1018 ucat[k][j][i].y * ucat[k][j][i].y +
1019 ucat[k][j][i].z * ucat[k][j][i].z);
1020 ke = 0.5 * speed * speed;
1021
1022 local[0] += cell_volume;
1023 local[1] += speed * cell_volume;
1024 local[2] += ke * cell_volume;
1025 }
1026 }
1027 }
1028
1029 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].Nvert, &nvert));
1030 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].Aj, &aj));
1031 PetscCall(DMDAVecRestoreArrayRead(user[bi].fda, user[bi].Ucat, &ucat));
1032 }
1033
1034 PetscCallMPI(MPI_Allreduce(local, global, 3, MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD));
1035
1036 if (global[0] <= 0.0) {
1037 *mean_speed_out = 0.0;
1038 *mean_ke_out = 0.0;
1039 } else {
1040 *mean_speed_out = global[1] / global[0];
1041 *mean_ke_out = global[2] / global[0];
1042 }
1043
1044 PetscFunctionReturn(0);
1045}
1046
1047/**
1048 * @brief Computes deterministic solution-drift metrics for the current step.
1049 *
1050 * This helper powers both `steady_deterministic` and
1051 * `periodic_deterministic` solution-convergence modes. It compares the
1052 * completed Eulerian state against either:
1053 * - the previous physical timestep (`Ucat_o`, `P_o`) for steady/transient use
1054 * - the stored phase-aligned snapshot ring for periodic use
1055 *
1056 * The routine performs two passes over the fluid cells:
1057 * 1. velocity/observable pass
1058 * - computes current and reference mean speed / mean KE
1059 * - computes absolute/relative velocity L2 drift
1060 * - accumulates pressure means needed for gauge removal
1061 * 2. pressure-only pass
1062 * - subtracts the volume-weighted mean pressure from both states
1063 * - computes gauge-invariant pressure L2 drift
1064 *
1065 * Warmup behavior is handled here. If no valid reference exists yet, all drift
1066 * outputs are left at zero and `has_reference_out` is set to `PETSC_FALSE`,
1067 * while current observables are still reported.
1068 *
1069 * @param[in] simCtx Simulation context owning the current state.
1070 * @param[in] periodic_mode `PETSC_TRUE` when comparing against
1071 * phase-aligned periodic storage.
1072 * @param[in] phase_step Active phase slot for periodic mode, or `-1`
1073 * when unused.
1074 * @param[in] samples_before Number of solution-convergence samples
1075 * already recorded before this timestep.
1076 * @param[out] has_reference_out Whether a valid comparison state existed.
1077 * @param[out] u_abs_l2_out Absolute L2 drift of Cartesian velocity.
1078 * @param[out] u_rel_l2_out Relative L2 drift of Cartesian velocity.
1079 * @param[out] p_abs_l2_out Gauge-invariant absolute L2 pressure drift.
1080 * @param[out] p_rel_l2_out Gauge-invariant relative L2 pressure drift.
1081 * @param[out] mean_speed_out Current volume-weighted mean speed.
1082 * @param[out] mean_speed_ref_out Reference volume-weighted mean speed.
1083 * @param[out] mean_speed_abs_out Absolute drift of mean speed.
1084 * @param[out] mean_speed_rel_out Relative drift of mean speed.
1085 * @param[out] mean_ke_out Current volume-weighted mean kinetic energy.
1086 * @param[out] mean_ke_ref_out Reference volume-weighted mean kinetic
1087 * energy.
1088 * @param[out] mean_ke_abs_out Absolute drift of mean kinetic energy.
1089 * @param[out] mean_ke_rel_out Relative drift of mean kinetic energy.
1090 * @return PetscErrorCode 0 on success.
1091 */
1092static PetscErrorCode ComputeDeterministicSolutionMetrics(SimCtx *simCtx,
1093 PetscBool periodic_mode,
1094 PetscInt phase_step,
1095 PetscInt samples_before,
1096 PetscBool *has_reference_out,
1097 PetscReal *u_abs_l2_out,
1098 PetscReal *u_rel_l2_out,
1099 PetscReal *p_abs_l2_out,
1100 PetscReal *p_rel_l2_out,
1101 PetscReal *mean_speed_out,
1102 PetscReal *mean_speed_ref_out,
1103 PetscReal *mean_speed_abs_out,
1104 PetscReal *mean_speed_rel_out,
1105 PetscReal *mean_ke_out,
1106 PetscReal *mean_ke_ref_out,
1107 PetscReal *mean_ke_abs_out,
1108 PetscReal *mean_ke_rel_out)
1109{
1110 SolutionConvergenceDeterministicPass1 local_pass1 = {0};
1111 SolutionConvergenceDeterministicPass1 global_pass1 = {0};
1112 SolutionConvergenceDeterministicPass2 local_pass2 = {0};
1113 SolutionConvergenceDeterministicPass2 global_pass2 = {0};
1114 PetscReal current_pressure_mean = 0.0;
1115 PetscReal reference_pressure_mean = 0.0;
1116 UserCtx *user = NULL;
1117
1118 PetscFunctionBeginUser;
1119 if (!simCtx || !has_reference_out || !u_abs_l2_out || !u_rel_l2_out || !p_abs_l2_out || !p_rel_l2_out ||
1120 !mean_speed_out || !mean_speed_ref_out || !mean_speed_abs_out || !mean_speed_rel_out ||
1121 !mean_ke_out || !mean_ke_ref_out || !mean_ke_abs_out || !mean_ke_rel_out) {
1122 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "ComputeDeterministicSolutionMetrics received a NULL output pointer.");
1123 }
1124
1125 *has_reference_out = periodic_mode
1126 ? (PetscBool)(simCtx->solutionConvergencePeriodSteps > 0 &&
1127 phase_step >= 0 &&
1128 phase_step < simCtx->solutionConvergencePeriodSteps &&
1129 samples_before >= simCtx->solutionConvergencePeriodSteps)
1130 : (PetscBool)(samples_before > 0);
1131
1132 *u_abs_l2_out = 0.0;
1133 *u_rel_l2_out = 0.0;
1134 *p_abs_l2_out = 0.0;
1135 *p_rel_l2_out = 0.0;
1136 *mean_speed_out = 0.0;
1137 *mean_speed_ref_out = 0.0;
1138 *mean_speed_abs_out = 0.0;
1139 *mean_speed_rel_out = 0.0;
1140 *mean_ke_out = 0.0;
1141 *mean_ke_ref_out = 0.0;
1142 *mean_ke_abs_out = 0.0;
1143 *mean_ke_rel_out = 0.0;
1144
1145 user = simCtx->usermg.mgctx[simCtx->usermg.mglevels - 1].user;
1146
1147 for (PetscInt bi = 0; bi < simCtx->block_number; ++bi) {
1148 const DMDALocalInfo info = user[bi].info;
1149 const PetscBool x_per = (PetscBool)(simCtx->i_periodic != 0);
1150 const PetscBool y_per = (PetscBool)(simCtx->j_periodic != 0);
1151 const PetscBool z_per = (PetscBool)(simCtx->k_periodic != 0);
1152 const PetscInt i_end = (x_per && (info.xs + info.xm == info.mx)) ? info.mx - 1 : info.xs + info.xm;
1153 const PetscInt j_end = (y_per && (info.ys + info.ym == info.my)) ? info.my - 1 : info.ys + info.ym;
1154 const PetscInt k_end = (z_per && (info.zs + info.zm == info.mz)) ? info.mz - 1 : info.zs + info.zm;
1155 Cmpnts ***ucat = NULL;
1156 Cmpnts ***ucat_ref = NULL;
1157 PetscReal ***pressure = NULL;
1158 PetscReal ***pressure_ref = NULL;
1159 PetscReal ***aj = NULL;
1160 PetscReal ***nvert = NULL;
1161 Vec ucat_reference_vec = NULL;
1162 Vec pressure_reference_vec = NULL;
1163
1164 if (*has_reference_out) {
1165 if (periodic_mode) {
1166 ucat_reference_vec = user[bi].solutionConvergencePeriodicUcatRef[phase_step];
1167 pressure_reference_vec = user[bi].solutionConvergencePeriodicPRef[phase_step];
1168 } else {
1169 ucat_reference_vec = user[bi].Ucat_o;
1170 pressure_reference_vec = user[bi].P_o;
1171 }
1172 }
1173
1174 PetscCall(DMDAVecGetArrayRead(user[bi].fda, user[bi].Ucat, &ucat));
1175 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].P, &pressure));
1176 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].Aj, &aj));
1177 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].Nvert, &nvert));
1178 if (*has_reference_out) {
1179 PetscCall(DMDAVecGetArrayRead(user[bi].fda, ucat_reference_vec, &ucat_ref));
1180 PetscCall(DMDAVecGetArrayRead(user[bi].da, pressure_reference_vec, &pressure_ref));
1181 }
1182
1183 for (PetscInt k = info.zs; k < k_end; ++k) {
1184 for (PetscInt j = info.ys; j < j_end; ++j) {
1185 for (PetscInt i = info.xs; i < i_end; ++i) {
1186 PetscReal jac = aj[k][j][i];
1187 PetscReal cell_volume = 0.0;
1188 PetscReal speed = 0.0;
1189 PetscReal ke = 0.0;
1190
1191 if (nvert[k][j][i] > SOLUTION_CONVERGENCE_FLUID_THRESHOLD) continue;
1192 if (PetscAbsReal(jac) <= 1.0e-14) continue;
1193
1194 cell_volume = 1.0 / jac;
1195 speed = PetscSqrtReal(ucat[k][j][i].x * ucat[k][j][i].x +
1196 ucat[k][j][i].y * ucat[k][j][i].y +
1197 ucat[k][j][i].z * ucat[k][j][i].z);
1198 ke = 0.5 * speed * speed;
1199
1200 local_pass1.fluid_volume += cell_volume;
1201 local_pass1.current_speed_sum += speed * cell_volume;
1202 local_pass1.current_ke_sum += ke * cell_volume;
1203 local_pass1.current_u_norm_sq += (ucat[k][j][i].x * ucat[k][j][i].x +
1204 ucat[k][j][i].y * ucat[k][j][i].y +
1205 ucat[k][j][i].z * ucat[k][j][i].z) * cell_volume;
1206
1207 if (*has_reference_out) {
1208 PetscReal ref_speed = PetscSqrtReal(ucat_ref[k][j][i].x * ucat_ref[k][j][i].x +
1209 ucat_ref[k][j][i].y * ucat_ref[k][j][i].y +
1210 ucat_ref[k][j][i].z * ucat_ref[k][j][i].z);
1211 PetscReal ref_ke = 0.5 * ref_speed * ref_speed;
1212 PetscReal dux = ucat[k][j][i].x - ucat_ref[k][j][i].x;
1213 PetscReal duy = ucat[k][j][i].y - ucat_ref[k][j][i].y;
1214 PetscReal duz = ucat[k][j][i].z - ucat_ref[k][j][i].z;
1215
1216 local_pass1.reference_speed_sum += ref_speed * cell_volume;
1217 local_pass1.reference_ke_sum += ref_ke * cell_volume;
1218 local_pass1.delta_u_norm_sq += (dux * dux + duy * duy + duz * duz) * cell_volume;
1219 local_pass1.current_pressure_sum += pressure[k][j][i] * cell_volume;
1220 local_pass1.reference_pressure_sum += pressure_ref[k][j][i] * cell_volume;
1221 }
1222 }
1223 }
1224 }
1225
1226 if (*has_reference_out) {
1227 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, pressure_reference_vec, &pressure_ref));
1228 PetscCall(DMDAVecRestoreArrayRead(user[bi].fda, ucat_reference_vec, &ucat_ref));
1229 }
1230 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].Nvert, &nvert));
1231 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].Aj, &aj));
1232 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].P, &pressure));
1233 PetscCall(DMDAVecRestoreArrayRead(user[bi].fda, user[bi].Ucat, &ucat));
1234 }
1235
1236 PetscCallMPI(MPI_Allreduce(&local_pass1, &global_pass1,
1237 sizeof(SolutionConvergenceDeterministicPass1) / sizeof(PetscReal),
1238 MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD));
1239
1240 if (global_pass1.fluid_volume <= 0.0) PetscFunctionReturn(0);
1241
1242 *mean_speed_out = global_pass1.current_speed_sum / global_pass1.fluid_volume;
1243 *mean_ke_out = global_pass1.current_ke_sum / global_pass1.fluid_volume;
1244
1245 if (!*has_reference_out) PetscFunctionReturn(0);
1246
1247 *mean_speed_ref_out = global_pass1.reference_speed_sum / global_pass1.fluid_volume;
1248 *mean_speed_abs_out = PetscAbsReal(*mean_speed_out - *mean_speed_ref_out);
1249 *mean_speed_rel_out = SolutionConvergenceSafeRelative(*mean_speed_abs_out, *mean_speed_out);
1250 *mean_ke_ref_out = global_pass1.reference_ke_sum / global_pass1.fluid_volume;
1251 *mean_ke_abs_out = PetscAbsReal(*mean_ke_out - *mean_ke_ref_out);
1252 *mean_ke_rel_out = SolutionConvergenceSafeRelative(*mean_ke_abs_out, *mean_ke_out);
1253
1254 current_pressure_mean = global_pass1.current_pressure_sum / global_pass1.fluid_volume;
1255 reference_pressure_mean = global_pass1.reference_pressure_sum / global_pass1.fluid_volume;
1256
1257 for (PetscInt bi = 0; bi < simCtx->block_number; ++bi) {
1258 const DMDALocalInfo info = user[bi].info;
1259 const PetscBool x_per = (PetscBool)(simCtx->i_periodic != 0);
1260 const PetscBool y_per = (PetscBool)(simCtx->j_periodic != 0);
1261 const PetscBool z_per = (PetscBool)(simCtx->k_periodic != 0);
1262 const PetscInt i_end = (x_per && (info.xs + info.xm == info.mx)) ? info.mx - 1 : info.xs + info.xm;
1263 const PetscInt j_end = (y_per && (info.ys + info.ym == info.my)) ? info.my - 1 : info.ys + info.ym;
1264 const PetscInt k_end = (z_per && (info.zs + info.zm == info.mz)) ? info.mz - 1 : info.zs + info.zm;
1265 PetscReal ***pressure = NULL;
1266 PetscReal ***pressure_ref = NULL;
1267 PetscReal ***aj = NULL;
1268 PetscReal ***nvert = NULL;
1269 Vec pressure_reference_vec = periodic_mode ? user[bi].solutionConvergencePeriodicPRef[phase_step] : user[bi].P_o;
1270
1271 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].P, &pressure));
1272 PetscCall(DMDAVecGetArrayRead(user[bi].da, pressure_reference_vec, &pressure_ref));
1273 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].Aj, &aj));
1274 PetscCall(DMDAVecGetArrayRead(user[bi].da, user[bi].Nvert, &nvert));
1275
1276 for (PetscInt k = info.zs; k < k_end; ++k) {
1277 for (PetscInt j = info.ys; j < j_end; ++j) {
1278 for (PetscInt i = info.xs; i < i_end; ++i) {
1279 PetscReal jac = aj[k][j][i];
1280 PetscReal cell_volume = 0.0;
1281 PetscReal current_pressure = 0.0;
1282 PetscReal reference_pressure = 0.0;
1283 PetscReal delta_pressure = 0.0;
1284
1285 if (nvert[k][j][i] > SOLUTION_CONVERGENCE_FLUID_THRESHOLD) continue;
1286 if (PetscAbsReal(jac) <= 1.0e-14) continue;
1287
1288 cell_volume = 1.0 / jac;
1289 current_pressure = pressure[k][j][i] - current_pressure_mean;
1290 reference_pressure = pressure_ref[k][j][i] - reference_pressure_mean;
1291 delta_pressure = current_pressure - reference_pressure;
1292
1293 local_pass2.current_pressure_norm_sq += current_pressure * current_pressure * cell_volume;
1294 local_pass2.delta_pressure_norm_sq += delta_pressure * delta_pressure * cell_volume;
1295 }
1296 }
1297 }
1298
1299 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].Nvert, &nvert));
1300 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].Aj, &aj));
1301 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, pressure_reference_vec, &pressure_ref));
1302 PetscCall(DMDAVecRestoreArrayRead(user[bi].da, user[bi].P, &pressure));
1303 }
1304
1305 PetscCallMPI(MPI_Allreduce(&local_pass2, &global_pass2,
1306 sizeof(SolutionConvergenceDeterministicPass2) / sizeof(PetscReal),
1307 MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD));
1308
1309 *u_abs_l2_out = PetscSqrtReal(global_pass1.delta_u_norm_sq);
1310 *u_rel_l2_out = SolutionConvergenceSafeRelative(*u_abs_l2_out, PetscSqrtReal(global_pass1.current_u_norm_sq));
1311 *p_abs_l2_out = PetscSqrtReal(global_pass2.delta_pressure_norm_sq);
1312 *p_rel_l2_out = SolutionConvergenceSafeRelative(*p_abs_l2_out, PetscSqrtReal(global_pass2.current_pressure_norm_sq));
1313
1314 PetscFunctionReturn(0);
1315}
1316
1317/**
1318 * @brief Reads one sample from the statistical ring buffer by age.
1319 *
1320 * The statistical logger stores scalar observables in a compact circular
1321 * buffer. This helper interprets the buffer using `samples_available` as the
1322 * logical end of the history and returns the entry `offset_from_latest` steps
1323 * back from the newest stored sample.
1324 *
1325 * Out-of-range requests return zero so warmup handling can remain simple and
1326 * deterministic.
1327 *
1328 * @param[in] history Ring-buffer storage array.
1329 * @param[in] capacity Total ring-buffer capacity.
1330 * @param[in] samples_available Number of logical samples available to read.
1331 * @param[in] offset_from_latest `0` means newest sample, `1` previous sample,
1332 * and so on.
1333 * @return Requested historical sample, or `0.0` if unavailable.
1334 */
1335static PetscReal SolutionConvergenceHistoryGet(const PetscReal *history,
1336 PetscInt capacity,
1337 PetscInt samples_available,
1338 PetscInt offset_from_latest)
1339{
1340 PetscInt count = 0;
1341 PetscInt index = 0;
1342
1343 if (!history || capacity <= 0 || samples_available <= 0 || offset_from_latest < 0) {
1344 return 0.0;
1345 }
1346
1347 count = PetscMin(samples_available, capacity);
1348 if (offset_from_latest >= count) return 0.0;
1349
1350 index = (samples_available - 1 - offset_from_latest) % capacity;
1351 if (index < 0) index += capacity;
1352 return history[index];
1353}
1354
1355/**
1356 * @brief Appends one timestep's scalar observables to the statistical history.
1357 *
1358 * Statistical solution-convergence compares adjacent windows of scalar
1359 * observables rather than full fields. This helper writes the current
1360 * `mean_speed` and `mean_ke` into the rolling history arrays using the current
1361 * sample count to choose the circular-buffer slot.
1362 *
1363 * @param[in,out] simCtx Simulation context owning the history arrays.
1364 * @param[in] samples_before Number of samples present before appending the
1365 * current timestep.
1366 * @param[in] mean_speed Current timestep mean-speed observable.
1367 * @param[in] mean_ke Current timestep mean-KE observable.
1368 * @return PetscErrorCode 0 on success.
1369 */
1370static PetscErrorCode AppendStatisticalObservableSample(SimCtx *simCtx,
1371 PetscInt samples_before,
1372 PetscReal mean_speed,
1373 PetscReal mean_ke)
1374{
1375 PetscInt history_capacity = 0;
1376 PetscInt slot = 0;
1377
1378 PetscFunctionBeginUser;
1379 if (!simCtx || !simCtx->solutionConvergenceMeanSpeedHistory || !simCtx->solutionConvergenceMeanKEHistory) {
1380 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Statistical solution-convergence history is not allocated.");
1381 }
1382
1383 history_capacity = 2 * simCtx->solutionConvergenceWindowSteps;
1384 if (history_capacity <= 0) {
1385 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Statistical solution-convergence history capacity must be positive.");
1386 }
1387
1388 slot = samples_before % history_capacity;
1389 simCtx->solutionConvergenceMeanSpeedHistory[slot] = mean_speed;
1390 simCtx->solutionConvergenceMeanKEHistory[slot] = mean_ke;
1391
1392 PetscFunctionReturn(0);
1393}
1394
1395/**
1396 * @brief Computes adjacent-window drift metrics for statistical steady mode.
1397 *
1398 * Once enough samples have been accumulated, this helper forms:
1399 * - the current window over the most recent `window_steps` samples
1400 * - the previous adjacent window over the preceding `window_steps` samples
1401 *
1402 * From those windows it computes means, RMS values, and absolute/relative
1403 * drift for both tracked observables (`mean_speed` and `mean_ke`). When the
1404 * history is still warming up:
1405 * - fewer than `window_steps` samples: no window metrics are available
1406 * - between `window_steps` and `2*window_steps - 1` samples: current-window
1407 * metrics are available, but no reference window exists yet
1408 *
1409 * @param[in] simCtx Simulation context owning the
1410 * statistical history.
1411 * @param[in] samples_available Number of samples available after
1412 * appending the current timestep.
1413 * @param[out] has_reference_out Whether both adjacent windows
1414 * exist and drift metrics are
1415 * meaningful.
1416 * @param[out] mean_speed_window_out Mean speed over the current
1417 * window.
1418 * @param[out] mean_speed_window_prev_out Mean speed over the previous
1419 * window.
1420 * @param[out] mean_speed_window_abs_out Absolute drift between current
1421 * and previous window means.
1422 * @param[out] mean_speed_window_rel_out Relative drift between current
1423 * and previous window means.
1424 * @param[out] mean_speed_rms_window_out RMS of mean-speed samples in the
1425 * current window.
1426 * @param[out] mean_speed_rms_window_prev_out RMS of mean-speed samples in the
1427 * previous window.
1428 * @param[out] mean_speed_rms_window_abs_out Absolute drift between window RMS
1429 * values.
1430 * @param[out] mean_speed_rms_window_rel_out Relative drift between window RMS
1431 * values.
1432 * @param[out] mean_ke_window_out Mean kinetic energy over the
1433 * current window.
1434 * @param[out] mean_ke_window_prev_out Mean kinetic energy over the
1435 * previous window.
1436 * @param[out] mean_ke_window_abs_out Absolute drift between current
1437 * and previous KE-window means.
1438 * @param[out] mean_ke_window_rel_out Relative drift between current
1439 * and previous KE-window means.
1440 * @param[out] mean_ke_rms_window_out RMS of mean-KE samples in the
1441 * current window.
1442 * @param[out] mean_ke_rms_window_prev_out RMS of mean-KE samples in the
1443 * previous window.
1444 * @param[out] mean_ke_rms_window_abs_out Absolute drift between KE-window
1445 * RMS values.
1446 * @param[out] mean_ke_rms_window_rel_out Relative drift between KE-window
1447 * RMS values.
1448 * @return PetscErrorCode 0 on success.
1449 */
1450static PetscErrorCode ComputeStatisticalWindowMetrics(const SimCtx *simCtx,
1451 PetscInt samples_available,
1452 PetscBool *has_reference_out,
1453 PetscReal *mean_speed_window_out,
1454 PetscReal *mean_speed_window_prev_out,
1455 PetscReal *mean_speed_window_abs_out,
1456 PetscReal *mean_speed_window_rel_out,
1457 PetscReal *mean_speed_rms_window_out,
1458 PetscReal *mean_speed_rms_window_prev_out,
1459 PetscReal *mean_speed_rms_window_abs_out,
1460 PetscReal *mean_speed_rms_window_rel_out,
1461 PetscReal *mean_ke_window_out,
1462 PetscReal *mean_ke_window_prev_out,
1463 PetscReal *mean_ke_window_abs_out,
1464 PetscReal *mean_ke_window_rel_out,
1465 PetscReal *mean_ke_rms_window_out,
1466 PetscReal *mean_ke_rms_window_prev_out,
1467 PetscReal *mean_ke_rms_window_abs_out,
1468 PetscReal *mean_ke_rms_window_rel_out)
1469{
1470 PetscInt w = 0;
1471 PetscInt history_capacity = 0;
1472 PetscReal speed_sum = 0.0;
1473 PetscReal speed_sum_sq = 0.0;
1474 PetscReal speed_prev_sum = 0.0;
1475 PetscReal speed_prev_sum_sq = 0.0;
1476 PetscReal ke_sum = 0.0;
1477 PetscReal ke_sum_sq = 0.0;
1478 PetscReal ke_prev_sum = 0.0;
1479 PetscReal ke_prev_sum_sq = 0.0;
1480
1481 PetscFunctionBeginUser;
1482 if (!simCtx || !has_reference_out || !mean_speed_window_out || !mean_speed_window_prev_out ||
1483 !mean_speed_window_abs_out || !mean_speed_window_rel_out || !mean_speed_rms_window_out ||
1484 !mean_speed_rms_window_prev_out || !mean_speed_rms_window_abs_out || !mean_speed_rms_window_rel_out ||
1485 !mean_ke_window_out || !mean_ke_window_prev_out || !mean_ke_window_abs_out || !mean_ke_window_rel_out ||
1486 !mean_ke_rms_window_out || !mean_ke_rms_window_prev_out || !mean_ke_rms_window_abs_out ||
1487 !mean_ke_rms_window_rel_out) {
1488 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "ComputeStatisticalWindowMetrics received a NULL output pointer.");
1489 }
1490
1491 *has_reference_out = PETSC_FALSE;
1492 *mean_speed_window_out = 0.0;
1493 *mean_speed_window_prev_out = 0.0;
1494 *mean_speed_window_abs_out = 0.0;
1495 *mean_speed_window_rel_out = 0.0;
1496 *mean_speed_rms_window_out = 0.0;
1497 *mean_speed_rms_window_prev_out = 0.0;
1498 *mean_speed_rms_window_abs_out = 0.0;
1499 *mean_speed_rms_window_rel_out = 0.0;
1500 *mean_ke_window_out = 0.0;
1501 *mean_ke_window_prev_out = 0.0;
1502 *mean_ke_window_abs_out = 0.0;
1503 *mean_ke_window_rel_out = 0.0;
1504 *mean_ke_rms_window_out = 0.0;
1505 *mean_ke_rms_window_prev_out = 0.0;
1506 *mean_ke_rms_window_abs_out = 0.0;
1507 *mean_ke_rms_window_rel_out = 0.0;
1508
1509 w = simCtx->solutionConvergenceWindowSteps;
1510 history_capacity = 2 * w;
1511 if (w <= 0 || samples_available < w) PetscFunctionReturn(0);
1512
1513 for (PetscInt idx = 0; idx < w; ++idx) {
1514 PetscReal speed_value = SolutionConvergenceHistoryGet(simCtx->solutionConvergenceMeanSpeedHistory,
1515 history_capacity,
1516 samples_available,
1517 idx);
1518 PetscReal ke_value = SolutionConvergenceHistoryGet(simCtx->solutionConvergenceMeanKEHistory,
1519 history_capacity,
1520 samples_available,
1521 idx);
1522 speed_sum += speed_value;
1523 speed_sum_sq += speed_value * speed_value;
1524 ke_sum += ke_value;
1525 ke_sum_sq += ke_value * ke_value;
1526 }
1527
1528 *mean_speed_window_out = speed_sum / (PetscReal)w;
1529 *mean_speed_rms_window_out = PetscSqrtReal(PetscMax(0.0, speed_sum_sq / (PetscReal)w -
1530 (*mean_speed_window_out) * (*mean_speed_window_out)));
1531 *mean_ke_window_out = ke_sum / (PetscReal)w;
1532 *mean_ke_rms_window_out = PetscSqrtReal(PetscMax(0.0, ke_sum_sq / (PetscReal)w -
1533 (*mean_ke_window_out) * (*mean_ke_window_out)));
1534
1535 if (samples_available < 2 * w) PetscFunctionReturn(0);
1536
1537 for (PetscInt idx = w; idx < 2 * w; ++idx) {
1538 PetscReal speed_value = SolutionConvergenceHistoryGet(simCtx->solutionConvergenceMeanSpeedHistory,
1539 history_capacity,
1540 samples_available,
1541 idx);
1542 PetscReal ke_value = SolutionConvergenceHistoryGet(simCtx->solutionConvergenceMeanKEHistory,
1543 history_capacity,
1544 samples_available,
1545 idx);
1546 speed_prev_sum += speed_value;
1547 speed_prev_sum_sq += speed_value * speed_value;
1548 ke_prev_sum += ke_value;
1549 ke_prev_sum_sq += ke_value * ke_value;
1550 }
1551
1552 *has_reference_out = PETSC_TRUE;
1553 *mean_speed_window_prev_out = speed_prev_sum / (PetscReal)w;
1554 *mean_speed_window_abs_out = PetscAbsReal(*mean_speed_window_out - *mean_speed_window_prev_out);
1555 *mean_speed_window_rel_out = SolutionConvergenceSafeRelative(*mean_speed_window_abs_out, *mean_speed_window_out);
1556 *mean_speed_rms_window_prev_out = PetscSqrtReal(PetscMax(0.0, speed_prev_sum_sq / (PetscReal)w -
1557 (*mean_speed_window_prev_out) * (*mean_speed_window_prev_out)));
1558 *mean_speed_rms_window_abs_out = PetscAbsReal(*mean_speed_rms_window_out - *mean_speed_rms_window_prev_out);
1559 *mean_speed_rms_window_rel_out = SolutionConvergenceSafeRelative(*mean_speed_rms_window_abs_out, *mean_speed_rms_window_out);
1560
1561 *mean_ke_window_prev_out = ke_prev_sum / (PetscReal)w;
1562 *mean_ke_window_abs_out = PetscAbsReal(*mean_ke_window_out - *mean_ke_window_prev_out);
1563 *mean_ke_window_rel_out = SolutionConvergenceSafeRelative(*mean_ke_window_abs_out, *mean_ke_window_out);
1564 *mean_ke_rms_window_prev_out = PetscSqrtReal(PetscMax(0.0, ke_prev_sum_sq / (PetscReal)w -
1565 (*mean_ke_window_prev_out) * (*mean_ke_window_prev_out)));
1566 *mean_ke_rms_window_abs_out = PetscAbsReal(*mean_ke_rms_window_out - *mean_ke_rms_window_prev_out);
1567 *mean_ke_rms_window_rel_out = SolutionConvergenceSafeRelative(*mean_ke_rms_window_abs_out, *mean_ke_rms_window_out);
1568
1569 PetscFunctionReturn(0);
1570}
1571
1572/**
1573 * @brief Maps the internal solution-convergence mode enum to its log label.
1574 *
1575 * The logger writes a human-readable mode string into the
1576 * `solution_convergence.log` banner and `mode` column. This helper keeps the
1577 * formatting centralized so the file output stays consistent with the accepted
1578 * configuration names.
1579 *
1580 * @param[in] mode Internal solution-convergence mode selector.
1581 * @return Lowercase string label written to the log output.
1582 */
1583static const char *SolutionConvergenceModeToString(SolutionConvergenceMode mode)
1584{
1585 switch (mode) {
1586 case SOLUTION_CONVERGENCE_STEADY_DETERMINISTIC: return "steady_deterministic";
1587 case SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC: return "periodic_deterministic";
1588 case SOLUTION_CONVERGENCE_STATISTICAL_STEADY: return "statistical_steady";
1589 case SOLUTION_CONVERGENCE_TRANSIENT: return "transient";
1590 default: return "unknown";
1591 }
1592}
1593
1594/**
1595 * @brief Implementation of \ref LOG_SOLUTION_CONVERGENCE().
1596 * @details Full API contract (arguments, ownership, side effects) is documented with
1597 * the header declaration in `include/logging.h`.
1598 * @see LOG_SOLUTION_CONVERGENCE()
1599 */
1600PetscErrorCode LOG_SOLUTION_CONVERGENCE(SimCtx *simCtx)
1601{
1602 PetscMPIInt rank = 0;
1603 PetscBool has_reference = PETSC_FALSE;
1604 PetscInt phase_step = -1;
1605 PetscInt samples_before = 0;
1606 PetscReal u_abs_l2 = 0.0, u_rel_l2 = 0.0, p_abs_l2 = 0.0, p_rel_l2 = 0.0;
1607 PetscReal mean_speed = 0.0, mean_speed_reference = 0.0, mean_speed_abs_drift = 0.0, mean_speed_rel_drift = 0.0;
1608 PetscReal mean_ke = 0.0, mean_ke_reference = 0.0, mean_ke_abs_drift = 0.0, mean_ke_rel_drift = 0.0;
1609 PetscReal mean_speed_window = 0.0, mean_speed_window_prev = 0.0, mean_speed_window_abs_drift = 0.0, mean_speed_window_rel_drift = 0.0;
1610 PetscReal mean_speed_rms_window = 0.0, mean_speed_rms_window_prev = 0.0, mean_speed_rms_window_abs_drift = 0.0, mean_speed_rms_window_rel_drift = 0.0;
1611 PetscReal mean_ke_window = 0.0, mean_ke_window_prev = 0.0, mean_ke_window_abs_drift = 0.0, mean_ke_window_rel_drift = 0.0;
1612 PetscReal mean_ke_rms_window = 0.0, mean_ke_rms_window_prev = 0.0, mean_ke_rms_window_abs_drift = 0.0, mean_ke_rms_window_rel_drift = 0.0;
1613
1614 PetscFunctionBeginUser;
1615 if (!simCtx) PetscFunctionReturn(0);
1616 if (simCtx->exec_mode != EXEC_MODE_SOLVER) PetscFunctionReturn(0);
1617 if (!simCtx->solutionConvergenceEnabled) PetscFunctionReturn(0);
1618
1619 samples_before = simCtx->solutionConvergenceSamplesRecorded;
1620
1621 switch (simCtx->solutionConvergenceMode) {
1622 case SOLUTION_CONVERGENCE_STEADY_DETERMINISTIC:
1623 case SOLUTION_CONVERGENCE_TRANSIENT:
1624 PetscCall(ComputeDeterministicSolutionMetrics(simCtx, PETSC_FALSE, -1, samples_before,
1625 &has_reference,
1626 &u_abs_l2, &u_rel_l2,
1627 &p_abs_l2, &p_rel_l2,
1628 &mean_speed, &mean_speed_reference,
1629 &mean_speed_abs_drift, &mean_speed_rel_drift,
1630 &mean_ke, &mean_ke_reference,
1631 &mean_ke_abs_drift, &mean_ke_rel_drift));
1632 break;
1633 case SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC:
1634 phase_step = simCtx->solutionConvergencePeriodSteps > 0 ? (simCtx->step % simCtx->solutionConvergencePeriodSteps) : -1;
1635 PetscCall(ComputeDeterministicSolutionMetrics(simCtx, PETSC_TRUE, phase_step, samples_before,
1636 &has_reference,
1637 &u_abs_l2, &u_rel_l2,
1638 &p_abs_l2, &p_rel_l2,
1639 &mean_speed, &mean_speed_reference,
1640 &mean_speed_abs_drift, &mean_speed_rel_drift,
1641 &mean_ke, &mean_ke_reference,
1642 &mean_ke_abs_drift, &mean_ke_rel_drift));
1643 break;
1644 case SOLUTION_CONVERGENCE_STATISTICAL_STEADY:
1645 PetscCall(ComputeCurrentFlowObservables(simCtx, &mean_speed, &mean_ke));
1646 PetscCall(AppendStatisticalObservableSample(simCtx, samples_before, mean_speed, mean_ke));
1647 PetscCall(ComputeStatisticalWindowMetrics(simCtx, samples_before + 1,
1648 &has_reference,
1649 &mean_speed_window, &mean_speed_window_prev,
1650 &mean_speed_window_abs_drift, &mean_speed_window_rel_drift,
1651 &mean_speed_rms_window, &mean_speed_rms_window_prev,
1652 &mean_speed_rms_window_abs_drift, &mean_speed_rms_window_rel_drift,
1653 &mean_ke_window, &mean_ke_window_prev,
1654 &mean_ke_window_abs_drift, &mean_ke_window_rel_drift,
1655 &mean_ke_rms_window, &mean_ke_rms_window_prev,
1656 &mean_ke_rms_window_abs_drift, &mean_ke_rms_window_rel_drift));
1657 break;
1658 default:
1659 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown solution convergence mode %d.", (int)simCtx->solutionConvergenceMode);
1660 }
1661
1662 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
1663 if (rank == 0) {
1664 char log_path[PETSC_MAX_PATH_LEN + 32];
1665 FILE *f = NULL;
1666 const char *mode_str = SolutionConvergenceModeToString(simCtx->solutionConvergenceMode);
1667
1668 PetscCall(PetscSNPrintf(log_path, sizeof(log_path), "%s/solution_convergence.log", simCtx->log_dir));
1669 f = fopen(log_path, "a");
1670 if (!f) {
1671 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open solution convergence log: %s", log_path);
1672 }
1673
1674 if (ftell(f) == 0) {
1675 switch (simCtx->solutionConvergenceMode) {
1676 case SOLUTION_CONVERGENCE_STEADY_DETERMINISTIC:
1677 case SOLUTION_CONVERGENCE_TRANSIENT:
1678 fprintf(f, "==================== Solution Convergence Log [mode: %s] ====================\n", mode_str);
1679 /* 16 columns; header width = 314 chars */
1680 fprintf(f, "%-10s | %-18s | %-22s | %-3s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s\n",
1681 "step", "time", "mode", "ref",
1682 "u_abs_l2", "u_rel_l2", "p_abs_l2", "p_rel_l2",
1683 "mean_speed", "spd_ref", "spd_abs", "spd_rel",
1684 "mean_ke", "ke_ref", "ke_abs", "ke_rel");
1685 fprintf(f, "----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n");
1686 break;
1687 case SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC:
1688 fprintf(f, "==================== Solution Convergence Log [mode: %s | period_steps: %d] ====================\n",
1689 mode_str, (int)simCtx->solutionConvergencePeriodSteps);
1690 /* 18 columns; header width = 330 chars */
1691 fprintf(f, "%-10s | %-18s | %-22s | %-3s | %-5s | %-5s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s\n",
1692 "step", "time", "mode", "ref", "ph", "per",
1693 "u_abs_l2", "u_rel_l2", "p_abs_l2", "p_rel_l2",
1694 "mean_speed", "spd_ref", "spd_abs", "spd_rel",
1695 "mean_ke", "ke_ref", "ke_abs", "ke_rel");
1696 fprintf(f, "----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n");
1697 break;
1698 case SOLUTION_CONVERGENCE_STATISTICAL_STEADY:
1699 fprintf(f, "==================== Solution Convergence Log [mode: %s | window_steps: %d] ====================\n",
1700 mode_str, (int)simCtx->solutionConvergenceWindowSteps);
1701 /* 21 columns; header width = 406 chars */
1702 fprintf(f, "%-10s | %-18s | %-22s | %-3s | %-5s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s | %-18s\n",
1703 "step", "time", "mode", "ref", "win",
1704 "mean_speed", "mean_ke",
1705 "spd_win", "spd_win_prev", "spd_win_abs", "spd_win_rel",
1706 "spd_rms_win", "spd_rms_abs", "spd_rms_rel",
1707 "ke_win", "ke_win_prev", "ke_win_abs", "ke_win_rel",
1708 "ke_rms_win", "ke_rms_abs", "ke_rms_rel");
1709 fprintf(f, "------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n");
1710 break;
1711 default: break;
1712 }
1713 }
1714 if (simCtx->continueMode && simCtx->step == simCtx->StartStep + 1) {
1715 fprintf(f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
1716 }
1717
1718 switch (simCtx->solutionConvergenceMode) {
1719 case SOLUTION_CONVERGENCE_STEADY_DETERMINISTIC:
1720 case SOLUTION_CONVERGENCE_TRANSIENT:
1721 fprintf(f,
1722 "%-10d | %-18.10e | %-22s | %-3d | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e\n",
1723 (int)simCtx->step, (double)simCtx->ti, mode_str, has_reference ? 1 : 0,
1724 (double)u_abs_l2, (double)u_rel_l2, (double)p_abs_l2, (double)p_rel_l2,
1725 (double)mean_speed, (double)mean_speed_reference,
1726 (double)mean_speed_abs_drift, (double)mean_speed_rel_drift,
1727 (double)mean_ke, (double)mean_ke_reference,
1728 (double)mean_ke_abs_drift, (double)mean_ke_rel_drift);
1729 break;
1730 case SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC:
1731 fprintf(f,
1732 "%-10d | %-18.10e | %-22s | %-3d | %-5d | %-5d | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e\n",
1733 (int)simCtx->step, (double)simCtx->ti, mode_str, has_reference ? 1 : 0,
1734 (int)phase_step, (int)simCtx->solutionConvergencePeriodSteps,
1735 (double)u_abs_l2, (double)u_rel_l2, (double)p_abs_l2, (double)p_rel_l2,
1736 (double)mean_speed, (double)mean_speed_reference,
1737 (double)mean_speed_abs_drift, (double)mean_speed_rel_drift,
1738 (double)mean_ke, (double)mean_ke_reference,
1739 (double)mean_ke_abs_drift, (double)mean_ke_rel_drift);
1740 break;
1741 case SOLUTION_CONVERGENCE_STATISTICAL_STEADY:
1742 fprintf(f,
1743 "%-10d | %-18.10e | %-22s | %-3d | %-5d | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e | %-18.10e\n",
1744 (int)simCtx->step, (double)simCtx->ti, mode_str, has_reference ? 1 : 0,
1745 (int)simCtx->solutionConvergenceWindowSteps,
1746 (double)mean_speed, (double)mean_ke,
1747 (double)mean_speed_window, (double)mean_speed_window_prev,
1748 (double)mean_speed_window_abs_drift, (double)mean_speed_window_rel_drift,
1749 (double)mean_speed_rms_window,
1750 (double)mean_speed_rms_window_abs_drift, (double)mean_speed_rms_window_rel_drift,
1751 (double)mean_ke_window, (double)mean_ke_window_prev,
1752 (double)mean_ke_window_abs_drift, (double)mean_ke_window_rel_drift,
1753 (double)mean_ke_rms_window,
1754 (double)mean_ke_rms_window_abs_drift, (double)mean_ke_rms_window_rel_drift);
1755 break;
1756 default: break;
1757 }
1758 fclose(f);
1759 }
1760
1761 if (simCtx->solutionConvergenceMode == SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC &&
1762 phase_step >= 0 && phase_step < simCtx->solutionConvergencePeriodSteps) {
1763 UserCtx *user = simCtx->usermg.mgctx[simCtx->usermg.mglevels - 1].user;
1764 for (PetscInt bi = 0; bi < simCtx->block_number; ++bi) {
1765 PetscCall(VecCopy(user[bi].Ucat, user[bi].solutionConvergencePeriodicUcatRef[phase_step]));
1766 PetscCall(VecCopy(user[bi].P, user[bi].solutionConvergencePeriodicPRef[phase_step]));
1767 }
1768 }
1769
1770 simCtx->solutionConvergenceSamplesRecorded = samples_before + 1;
1771
1772 PetscFunctionReturn(0);
1773}
1774
1775/**
1776 * @brief Logs continuity metrics for a single block to a file.
1777 *
1778 * This function should be called for each block, once per timestep. It opens a
1779 * central log file in append mode. To ensure the header is written only once,
1780 * it checks if it is processing block 0 on the simulation's start step.
1781 *
1782 * @param user A pointer to the UserCtx for the specific block whose metrics
1783 * are to be logged. The function accesses both global (SimCtx)
1784 * and local (user->...) data.
1785 * @return PetscErrorCode 0 on success.
1786 */
1787#undef __FUNCT__
1788#define __FUNCT__ "LOG_CONTINUITY_METRICS"
1789/**
1790 * @brief Implementation of \ref LOG_CONTINUITY_METRICS().
1791 * @details Full API contract (arguments, ownership, side effects) is documented with
1792 * the header declaration in `include/logging.h`.
1793 * @see LOG_CONTINUITY_METRICS()
1794 */
1795
1796PetscErrorCode LOG_CONTINUITY_METRICS(UserCtx *user)
1797{
1798 PetscErrorCode ierr;
1799 PetscMPIInt rank;
1800 SimCtx *simCtx = user->simCtx; // Get the shared SimCtx
1801 const PetscInt bi = user->_this; // Get this block's specific ID
1802 const PetscInt ti = simCtx->step; // Get the current timestep
1803
1804 PetscFunctionBeginUser;
1805 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
1806
1807 // Only rank 0 performs file I/O.
1808 if (!rank) {
1809 FILE *f;
1810 char filen[PETSC_MAX_PATH_LEN + 64];
1811 ierr = PetscSNPrintf(filen, sizeof(filen), "%s/Continuity_Metrics.log", simCtx->log_dir); CHKERRQ(ierr);
1812
1813 // Open the log file in append mode.
1814 f = fopen(filen, "a");
1815 if (!f) {
1816 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open log file: %s", filen);
1817 }
1818
1819 // Write a header only when the file is empty and it's the first block (bi=0).
1820 // Using ftell() instead of step comparison ensures correctness across continuations.
1821 if (ftell(f) == 0 && bi == 0) {
1822 PetscFPrintf(PETSC_COMM_SELF, f, "%-10s | %-6s | %-18s | %-30s | %-18s | %-18s | %-18s | %-18s\n",
1823 "Timestep", "Block", "Max Divergence", "Max Divergence Location ([k][j][i]=idx)", "Sum(RHS)","Total Flux In", "Total Flux Out", "Net Flux");
1824 PetscFPrintf(PETSC_COMM_SELF, f, "------------------------------------------------------------------------------------------------------------------------------------------\n");
1825 }
1826 if (simCtx->continueMode && ti == simCtx->StartStep + 1 && bi == 0) {
1827 PetscFPrintf(PETSC_COMM_SELF, f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
1828 }
1829
1830 // Prepare the data strings and values for the current block.
1831 PetscReal net_flux = simCtx->FluxInSum - simCtx->FluxOutSum;
1832 char location_str[64];
1833 sprintf(location_str, "([%d][%d][%d] = %d)", (int)simCtx->MaxDivz, (int)simCtx->MaxDivy, (int)simCtx->MaxDivx, (int)simCtx->MaxDivFlatArg);
1834
1835 // Write the formatted line for the current block.
1836 PetscFPrintf(PETSC_COMM_SELF, f, "%-10d | %-6d | %-18.10e | %-39s | %-18.10e | %-18.10e | %-18.10e | %-18.10e\n",
1837 (int)ti,
1838 (int)bi,
1839 (double)simCtx->MaxDiv,
1840 location_str,
1841 (double)simCtx->summationRHS,
1842 (double)simCtx->FluxInSum,
1843 (double)simCtx->FluxOutSum,
1844 (double)net_flux);
1845
1846 fclose(f);
1847 }
1848
1849 PetscFunctionReturn(0);
1850}
1851
1852/**
1853 * @brief Implementation of \ref ParticleLocationStatusToString().
1854 * @details Full API contract (arguments, ownership, side effects) is documented with
1855 * the header declaration in `include/logging.h`.
1856 * @see ParticleLocationStatusToString()
1857 */
1858const char* ParticleLocationStatusToString(ParticleLocationStatus level)
1859{
1860 switch (level) {
1861 case NEEDS_LOCATION: return "NEEDS_LOCATION";
1862 case ACTIVE_AND_LOCATED: return "ACTIVE_AND_LOCATED";
1863 case MIGRATING_OUT: return "MIGRATING_OUT";
1864 case LOST: return "LOST";
1865 case UNINITIALIZED: return "UNINITIALIZED";
1866 default: return "UNKNOWN_LEVEL";
1867 }
1868}
1869
1870///////// Profiling System /////////
1871
1872// Data structure to hold profiling info for one function
1873typedef struct {
1874 const char *name;
1875 double total_time;
1876 double current_step_time;
1877 long long total_call_count;
1878 long long current_step_call_count;
1879 double start_time; // Timer for the current call
1880 PetscBool always_log;
1881} ProfiledFunction;
1882
1883// Global registry for all profiled functions
1884static ProfiledFunction *g_profiler_registry = NULL;
1885static PetscInt g_profiler_count = 0;
1886static PetscInt g_profiler_capacity = 0;
1887
1888// Internal helper to find a function in the registry or create it
1889/**
1890 * @brief Find a profiling record by name or allocate and register a new record.
1891 */
1892static PetscErrorCode _FindOrCreateEntry(const char *func_name, PetscInt *idx)
1893{
1894 PetscFunctionBeginUser;
1895 // Search for existing entry
1896 for (PetscInt i = 0; i < g_profiler_count; ++i) {
1897 if (strcmp(g_profiler_registry[i].name, func_name) == 0) {
1898 *idx = i;
1899 PetscFunctionReturn(0);
1900 }
1901 }
1902
1903 // Not found, create a new entry
1904 if (g_profiler_count >= g_profiler_capacity) {
1905 PetscInt new_capacity = g_profiler_capacity == 0 ? 16 : g_profiler_capacity * 2;
1906 PetscErrorCode ierr = PetscRealloc(sizeof(ProfiledFunction) * new_capacity, &g_profiler_registry); CHKERRQ(ierr);
1907 g_profiler_capacity = new_capacity;
1908 }
1909
1910 *idx = g_profiler_count;
1911 g_profiler_registry[*idx].name = func_name;
1912 g_profiler_registry[*idx].total_time = 0.0;
1913 g_profiler_registry[*idx].current_step_time = 0.0;
1914 g_profiler_registry[*idx].total_call_count = 0;
1915 g_profiler_registry[*idx].current_step_call_count = 0;
1916 g_profiler_registry[*idx].start_time = 0.0;
1917 g_profiler_registry[*idx].always_log = PETSC_FALSE;
1918 g_profiler_count++;
1919
1920 PetscFunctionReturn(0);
1921}
1922
1923// --- Public API Implementation ---
1924/**
1925 * @brief Internal helper implementation: `ProfilingInitialize()`.
1926 * @details Local to this translation unit.
1927 */
1928PetscErrorCode ProfilingInitialize(SimCtx *simCtx)
1929{
1930 PetscFunctionBeginUser;
1931 if (!simCtx) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "SimCtx cannot be null for ProfilingInitialize");
1932
1933 // Iterate through the list of critical functions provided in SimCtx
1934 for (PetscInt i = 0; i < simCtx->nProfilingSelectedFuncs; ++i) {
1935 PetscInt idx;
1936 const char *func_name = simCtx->profilingSelectedFuncs[i];
1937 PetscErrorCode ierr = _FindOrCreateEntry(func_name, &idx); CHKERRQ(ierr);
1938 g_profiler_registry[idx].always_log = PETSC_TRUE;
1939
1940 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Marked '%s' as a critical function for profiling.\n", func_name);
1941 }
1942 PetscFunctionReturn(0);
1943}
1944
1945/**
1946 * @brief Implementation of \ref _ProfilingStart().
1947 * @details Full API contract (arguments, ownership, side effects) is documented with
1948 * the header declaration in `include/logging.h`.
1949 * @see _ProfilingStart()
1950 */
1951
1952void _ProfilingStart(const char *func_name)
1953{
1954 PetscInt idx;
1955 if (_FindOrCreateEntry(func_name, &idx) != 0) return; // Fail silently
1956 PetscTime(&g_profiler_registry[idx].start_time);
1957}
1958
1959/**
1960 * @brief Implementation of \ref _ProfilingEnd().
1961 * @details Full API contract (arguments, ownership, side effects) is documented with
1962 * the header declaration in `include/logging.h`.
1963 * @see _ProfilingEnd()
1964 */
1965
1966void _ProfilingEnd(const char *func_name)
1967{
1968 double end_time;
1969 PetscTime(&end_time);
1970
1971 PetscInt idx;
1972 if (_FindOrCreateEntry(func_name, &idx) != 0) return; // Fail silently
1973
1974 double elapsed = end_time - g_profiler_registry[idx].start_time;
1975 g_profiler_registry[idx].total_time += elapsed;
1976 g_profiler_registry[idx].current_step_time += elapsed;
1977 g_profiler_registry[idx].total_call_count++;
1978 g_profiler_registry[idx].current_step_call_count++;
1979}
1980
1981/**
1982 * @brief Implementation of \ref ProfilingResetTimestepCounters().
1983 * @details Full API contract (arguments, ownership, side effects) is documented with
1984 * the header declaration in `include/logging.h`.
1985 * @see ProfilingResetTimestepCounters()
1986 */
1987
1988PetscErrorCode ProfilingResetTimestepCounters(void)
1989{
1990 PetscFunctionBeginUser;
1991 for (PetscInt i = 0; i < g_profiler_count; ++i) {
1992 g_profiler_registry[i].current_step_time = 0.0;
1993 g_profiler_registry[i].current_step_call_count = 0;
1994 }
1995 PetscFunctionReturn(0);
1996}
1997
1998/**
1999 * @brief Implementation of \ref ProfilingLogTimestepSummary().
2000 * @details Full API contract (arguments, ownership, side effects) is documented with
2001 * the header declaration in `include/logging.h`.
2002 * @see ProfilingLogTimestepSummary()
2003 */
2004
2005PetscErrorCode ProfilingLogTimestepSummary(SimCtx *simCtx, PetscInt step)
2006{
2007 PetscBool should_write = PETSC_FALSE;
2008 FILE *f = NULL;
2009 char filen[(2 * PETSC_MAX_PATH_LEN) + 16];
2010
2011 PetscFunctionBeginUser;
2012 if (!simCtx) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "SimCtx cannot be null for ProfilingLogTimestepSummary");
2013
2014 if (strcmp(simCtx->profilingTimestepMode, "off") == 0) {
2015 for (PetscInt i = 0; i < g_profiler_count; ++i) {
2016 g_profiler_registry[i].current_step_time = 0.0;
2017 g_profiler_registry[i].current_step_call_count = 0;
2018 }
2019 PetscFunctionReturn(0);
2020 }
2021
2022 for (PetscInt i = 0; i < g_profiler_count; ++i) {
2023 if (g_profiler_registry[i].current_step_call_count <= 0) {
2024 continue;
2025 }
2026 if (strcmp(simCtx->profilingTimestepMode, "all") == 0 || g_profiler_registry[i].always_log) {
2027 should_write = PETSC_TRUE;
2028 break;
2029 }
2030 }
2031
2032 if (should_write && simCtx->rank == 0) {
2033 snprintf(filen, sizeof(filen), "%s/%s", simCtx->log_dir, simCtx->profilingTimestepFile);
2034 if (step == simCtx->StartStep + 1 && !simCtx->continueMode) {
2035 f = fopen(filen, "w");
2036 if (!f) {
2037 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open profiling timestep log file: %s", filen);
2038 }
2039 PetscFPrintf(PETSC_COMM_SELF, f, "step,function,calls,step_time_s\n");
2040 } else {
2041 f = fopen(filen, "a");
2042 if (!f) {
2043 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open profiling timestep log file: %s", filen);
2044 }
2045 if (step == simCtx->StartStep + 1 && ftell(f) == 0) {
2046 PetscFPrintf(PETSC_COMM_SELF, f, "step,function,calls,step_time_s\n");
2047 }
2048 }
2049 if (simCtx->continueMode && step == simCtx->StartStep + 1) {
2050 PetscFPrintf(PETSC_COMM_SELF, f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
2051 }
2052
2053 for (PetscInt i = 0; i < g_profiler_count; ++i) {
2054 if (g_profiler_registry[i].current_step_call_count <= 0) {
2055 continue;
2056 }
2057 if (strcmp(simCtx->profilingTimestepMode, "all") == 0 || g_profiler_registry[i].always_log) {
2058 PetscFPrintf(
2059 PETSC_COMM_SELF,
2060 f,
2061 "%d,%s,%lld,%.6f\n",
2062 (int)step,
2063 g_profiler_registry[i].name,
2064 g_profiler_registry[i].current_step_call_count,
2065 g_profiler_registry[i].current_step_time
2066 );
2067 }
2068 }
2069 fclose(f);
2070 }
2071
2072 // Reset per-step counters for the next iteration
2073 for (PetscInt i = 0; i < g_profiler_count; ++i) {
2074 g_profiler_registry[i].current_step_time = 0.0;
2075 g_profiler_registry[i].current_step_call_count = 0;
2076 }
2077 PetscFunctionReturn(0);
2078}
2079
2080/**
2081 * @brief Implementation of \ref RuntimeMemoryLogSample().
2082 * @details Full API contract (arguments, ownership, side effects) is documented with
2083 * the header declaration in `include/logging.h`.
2084 * @see RuntimeMemoryLogSample()
2085 */
2086PetscErrorCode RuntimeMemoryLogSample(SimCtx *simCtx, PetscInt step, const char *event, const char *reason)
2087{
2088 PetscErrorCode ierr;
2089 PetscLogDouble process_current_bytes = 0.0;
2090 PetscLogDouble process_peak_bytes = 0.0;
2091 PetscLogDouble petsc_current_bytes = 0.0;
2092 PetscLogDouble petsc_peak_bytes = 0.0;
2093 PetscReal local_values[5];
2094 PetscReal global_values[5];
2095 PetscReal process_current_mb = 0.0;
2096 PetscReal process_peak_mb = 0.0;
2097 PetscReal petsc_current_mb = 0.0;
2098 PetscReal petsc_peak_mb = 0.0;
2099 PetscReal process_change_mb = 0.0;
2100 char path[(2 * PETSC_MAX_PATH_LEN) + 16];
2101 FILE *f = NULL;
2102
2103 PetscFunctionBeginUser;
2104 if (!simCtx) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "SimCtx cannot be null for RuntimeMemoryLogSample");
2105 if (!simCtx->runtimeMemoryLogEnabled) PetscFunctionReturn(0);
2106
2107 ierr = PetscMemoryGetCurrentUsage(&process_current_bytes); CHKERRQ(ierr);
2108 ierr = PetscMemoryGetMaximumUsage(&process_peak_bytes); CHKERRQ(ierr);
2109 ierr = PetscMallocGetCurrentUsage(&petsc_current_bytes); CHKERRQ(ierr);
2110 ierr = PetscMallocGetMaximumUsage(&petsc_peak_bytes); CHKERRQ(ierr);
2111
2112 process_current_mb = (PetscReal)(process_current_bytes / (1024.0 * 1024.0));
2113 process_peak_mb = (PetscReal)(process_peak_bytes / (1024.0 * 1024.0));
2114 petsc_current_mb = (PetscReal)(petsc_current_bytes / (1024.0 * 1024.0));
2115 petsc_peak_mb = (PetscReal)(petsc_peak_bytes / (1024.0 * 1024.0));
2116 if (simCtx->runtimeMemoryLogHasPrevious) {
2117 process_change_mb = process_current_mb - simCtx->runtimeMemoryLogPreviousProcessMB;
2118 }
2119
2120 local_values[0] = process_current_mb;
2121 local_values[1] = process_peak_mb;
2122 local_values[2] = petsc_current_mb;
2123 local_values[3] = petsc_peak_mb;
2124 local_values[4] = process_change_mb;
2125 ierr = MPI_Allreduce(local_values, global_values, 5, MPIU_REAL, MPI_MAX, PETSC_COMM_WORLD); CHKERRMPI(ierr);
2126
2127 simCtx->runtimeMemoryLogPreviousProcessMB = process_current_mb;
2128 simCtx->runtimeMemoryLogHasPrevious = PETSC_TRUE;
2129
2130 if (simCtx->rank == 0) {
2131 ierr = PetscSNPrintf(path, sizeof(path), "%s/%s", simCtx->log_dir, simCtx->runtimeMemoryLogFile); CHKERRQ(ierr);
2132 f = fopen(path, "a");
2133 if (!f) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open runtime memory log file: %s", path);
2134
2135 if (!simCtx->runtimeMemoryLogStarted) {
2136 fprintf(f, "# PICurv runtime memory log\n");
2137 if (simCtx->continueMode) {
2138 fprintf(f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
2139 }
2140 fprintf(
2141 f,
2142 "%-8s %-10s %22s %20s %22s %28s %22s %-18s\n",
2143 "Step",
2144 "Event",
2145 "Process Current MB Max",
2146 "Process Peak MB Max",
2147 "PETSc Allocated MB Max",
2148 "PETSc Peak Allocated MB Max",
2149 "Process Change MB Max",
2150 "Reason"
2151 );
2152 simCtx->runtimeMemoryLogStarted = PETSC_TRUE;
2153 }
2154
2155 fprintf(
2156 f,
2157 "%-8" PetscInt_FMT " %-10s %22.3f %20.3f %22.3f %28.3f %22.3f %-18s\n",
2158 step,
2159 event ? event : "-",
2160 (double)global_values[0],
2161 (double)global_values[1],
2162 (double)global_values[2],
2163 (double)global_values[3],
2164 (double)global_values[4],
2165 (reason && reason[0]) ? reason : "-"
2166 );
2167 if ((event && (strcmp(event, "Shutdown") == 0 || strcmp(event, "Final") == 0))) {
2168 fflush(f);
2169 }
2170 fclose(f);
2171 }
2172
2173 PetscFunctionReturn(0);
2174}
2175
2176// Comparison function for qsort to sort by total_time in descending order
2177/**
2178 * @brief Order profiling records by their accumulated execution time.
2179 */
2180static int _CompareProfiledFunctions(const void *a, const void *b)
2181{
2182 const ProfiledFunction *func_a = (const ProfiledFunction *)a;
2183 const ProfiledFunction *func_b = (const ProfiledFunction *)b;
2184
2185 if (func_a->total_time < func_b->total_time) return 1;
2186 if (func_a->total_time > func_b->total_time) return -1;
2187 return 0;
2188}
2189
2190/**
2191 * @brief Implementation of \ref ProfilingFinalize().
2192 * @details Full API contract (arguments, ownership, side effects) is documented with
2193 * the header declaration in `include/logging.h`.
2194 * @see ProfilingFinalize()
2195 */
2196PetscErrorCode ProfilingFinalize(SimCtx *simCtx)
2197{
2198 PetscErrorCode ierr;
2199 PetscInt rank = simCtx->rank;
2200 PetscFunctionBeginUser;
2201 if (!simCtx->profilingFinalSummary) PetscFunctionReturn(0);
2202 if (!rank) {
2203
2204 char exec_mode_modifier[32] = "Unknown";
2205 if(simCtx->exec_mode == EXEC_MODE_SOLVER) PetscCall(PetscStrncpy(exec_mode_modifier, "Solver", sizeof(exec_mode_modifier)));
2206 else if(simCtx->exec_mode == EXEC_MODE_POSTPROCESSOR) PetscCall(PetscStrncpy(exec_mode_modifier, "PostProcessor", sizeof(exec_mode_modifier)));
2207 //--- Step 0: Create a file viewer for log file
2208 FILE *f;
2209 char filen[PETSC_MAX_PATH_LEN + 128];
2210 ierr = PetscSNPrintf(filen, sizeof(filen), "%s/ProfilingSummary_%s.log",simCtx->log_dir,exec_mode_modifier); CHKERRQ(ierr);
2211
2212 // Open the log file: append with section label in continue mode, truncate otherwise.
2213 if (simCtx->continueMode) {
2214 f = fopen(filen, "a");
2215 if (!f) {
2216 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open log file: %s", filen);
2217 }
2218 fprintf(f, "\n=== Continuation from step %" PetscInt_FMT " ===\n", simCtx->StartStep);
2219 } else {
2220 f = fopen(filen, "w");
2221 if (!f) {
2222 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open log file: %s", filen);
2223 }
2224 }
2225
2226 // --- Step 1: Sort the data for readability ---
2227 qsort(g_profiler_registry, g_profiler_count, sizeof(ProfiledFunction), _CompareProfiledFunctions);
2228
2229 // --- Step 2: Dynamically determine the width for the function name column ---
2230 PetscInt max_name_len = strlen("Function"); // Start with the header's length
2231 for (PetscInt i = 0; i < g_profiler_count; ++i) {
2232 if (g_profiler_registry[i].total_call_count > 0) {
2233 PetscInt len = strlen(g_profiler_registry[i].name);
2234 if (len > max_name_len) {
2235 max_name_len = len;
2236 }
2237 }
2238 }
2239 // Add a little padding
2240 max_name_len += 2;
2241
2242 // --- Step 3: Define fixed widths for numeric columns for consistent alignment ---
2243 const int time_width = 18;
2244 const int count_width = 15;
2245 const int avg_width = 22;
2246
2247 // --- Step 4: Print the formatted table ---
2248 PetscFPrintf(PETSC_COMM_SELF, f, "=================================================================================================================\n");
2249 PetscFPrintf(PETSC_COMM_SELF, f, " FINAL PROFILING SUMMARY (Sorted by Total Time)\n");
2250 PetscFPrintf(PETSC_COMM_SELF, f, "=================================================================================================================\n");
2251
2252 // Header Row
2253 PetscFPrintf(PETSC_COMM_SELF, f, "%-*s | %-*s | %-*s | %-*s\n",
2254 max_name_len, "Function",
2255 time_width, "Total Time (s)",
2256 count_width, "Call Count",
2257 avg_width, "Avg. Time/Call (ms)");
2258
2259 // Separator Line (dynamically sized)
2260 for (int i = 0; i < max_name_len; i++) PetscFPrintf(PETSC_COMM_SELF, f, "-");
2261 PetscFPrintf(PETSC_COMM_SELF, f, "-|-");
2262 for (int i = 0; i < time_width; i++) PetscFPrintf(PETSC_COMM_SELF, f, "-");
2263 PetscFPrintf(PETSC_COMM_SELF, f, "-|-");
2264 for (int i = 0; i < count_width; i++) PetscFPrintf(PETSC_COMM_SELF, f, "-");
2265 PetscFPrintf(PETSC_COMM_SELF, f, "-|-");
2266 for (int i = 0; i < avg_width; i++) PetscFPrintf(PETSC_COMM_SELF, f, "-");
2267 PetscFPrintf(PETSC_COMM_SELF, f, "\n");
2268
2269 // Data Rows
2270 for (PetscInt i = 0; i < g_profiler_count; ++i) {
2271 if (g_profiler_registry[i].total_call_count > 0) {
2272 double avg_time_ms = (g_profiler_registry[i].total_time / g_profiler_registry[i].total_call_count) * 1000.0;
2273 PetscFPrintf(PETSC_COMM_SELF, f, "%-*s | %*.*f | %*lld | %*.*f\n",
2274 max_name_len, g_profiler_registry[i].name,
2275 time_width, 6, g_profiler_registry[i].total_time,
2276 count_width, g_profiler_registry[i].total_call_count,
2277 avg_width, 6, avg_time_ms);
2278 PetscFPrintf(PETSC_COMM_SELF, f, "------------------------------------------------------------------------------------------------------------------\n");
2279 }
2280 }
2281 PetscFPrintf(PETSC_COMM_SELF, f, "==================================================================================================================\n");
2282
2283 fclose(f);
2284 }
2285
2286 // --- Final Cleanup ---
2287 PetscFree(g_profiler_registry);
2288 g_profiler_registry = NULL;
2289 g_profiler_count = 0;
2290 g_profiler_capacity = 0;
2291 PetscFunctionReturn(0);
2292}
2293
2294/*================================================================================*
2295 * PROGRESS BAR UTILITY *
2296 *================================================================================*/
2297
2298/**
2299 * @brief Internal helper implementation: `PrintProgressBar()`.
2300 * @details Local to this translation unit.
2301 */
2302void PrintProgressBar(PetscInt step, PetscInt startStep, PetscInt totalSteps, PetscReal currentTime)
2303{
2304 if (totalSteps <= 0) return;
2305
2306 // --- Configuration ---
2307 const int barWidth = 50;
2308
2309 // --- Calculation ---
2310 // Calculate progress as a fraction from 0.0 to 1.0
2311 PetscReal progress = (PetscReal)(step - startStep + 1) / totalSteps;
2312 // Ensure progress doesn't exceed 1.0 due to floating point inaccuracies
2313 if (progress > 1.0) progress = 1.0;
2314
2315 int pos = (int)(barWidth * progress);
2316
2317 // --- Printing ---
2318 // Carriage return moves cursor to the beginning of the line
2319 PetscPrintf(PETSC_COMM_SELF, "\rProgress: [");
2320
2321 for (int i = 0; i < barWidth; ++i) {
2322 if (i < pos) {
2323 PetscPrintf(PETSC_COMM_SELF, "=");
2324 } else if (i == pos) {
2325 PetscPrintf(PETSC_COMM_SELF, ">");
2326 } else {
2327 PetscPrintf(PETSC_COMM_SELF, " ");
2328 }
2329 }
2330
2331 // Print percentage, step count, and current time
2332 PetscPrintf(PETSC_COMM_SELF, "] %3d%% (Step %" PetscInt_FMT "/%" PetscInt_FMT ", t=%.4f)",
2333 (int)(progress * 100.0),
2334 step + 1,
2335 startStep + totalSteps,
2336 currentTime);
2337
2338 // Flush the output buffer to ensure the bar is displayed immediately
2339 fflush(stdout);
2340}
2341
2342#undef __FUNCT__
2343#define __FUNCT__ "LOG_FIELD_MIN_MAX"
2344/**
2345 * @brief Implementation of \ref LOG_FIELD_MIN_MAX().
2346 * @details Full API contract is documented with the header declaration in `include/logging.h`.
2347 * @see LOG_FIELD_MIN_MAX()
2348 */
2349PetscErrorCode LOG_FIELD_MIN_MAX(UserCtx *user, FieldId field_id)
2350{
2351 PetscErrorCode ierr;
2352 PetscInt i, j, k;
2353 DMDALocalInfo info;
2354
2355 FieldView view;
2356 Vec fieldVec = NULL;
2357 DM dm = NULL;
2358 PetscInt dof;
2359 FieldLayout layout;
2360 const char *fieldName = NULL;
2361 const char *data_layout = NULL;
2362
2363 PetscFunctionBeginUser;
2364
2365 ierr = FieldGetView(user, field_id, &view); CHKERRQ(ierr);
2366 fieldName = view.descriptor->canonical_name;
2367 dm = view.dm;
2368 dof = view.descriptor->dof;
2369 layout = view.descriptor->layout;
2370 data_layout = FieldLayoutName(layout);
2371 fieldVec = (layout == FIELD_LAYOUT_COMPONENT_STAGGERED) ? view.local_vec : view.global_vec;
2372
2373 ierr = DMDAGetLocalInfo(dm, &info); CHKERRQ(ierr);
2374
2375 // --- 2. Define Architecture-Aware Loop Bounds ---
2376 PetscInt i_start, i_end, j_start, j_end, k_start, k_end;
2377
2378 if (layout == FIELD_LAYOUT_CELL_CENTERED) {
2379 // For cell-centered data, the physical values are stored from index 1 to N-1.
2380 // We find the intersection of the rank's owned range [xs, xe) with the
2381 // physical data range [1, IM-1).
2382 i_start = PetscMax(info.xs, 1); i_end = PetscMin(info.xs + info.xm, user->IM);
2383 j_start = PetscMax(info.ys, 1); j_end = PetscMin(info.ys + info.ym, user->JM);
2384 k_start = PetscMax(info.zs, 1); k_end = PetscMin(info.zs + info.zm, user->KM);
2385 } else { // For Node- or Face-Centered data
2386 // The physical values are stored from index 0 to N-1.
2387 // We find the intersection of the rank's owned range [xs, xe) with the
2388 // physical data range [0, IM-1].
2389 i_start = PetscMax(info.xs, 0); i_end = PetscMin(info.xs + info.xm, user->IM);
2390 j_start = PetscMax(info.ys, 0); j_end = PetscMin(info.ys + info.ym, user->JM);
2391 k_start = PetscMax(info.zs, 0); k_end = PetscMin(info.zs + info.zm, user->KM);
2392 }
2393
2394 // --- 3. Barrier for clean, grouped output ---
2395 ierr = MPI_Barrier(PETSC_COMM_WORLD); CHKERRQ(ierr);
2396 if (user->simCtx->rank == 0) {
2397 PetscPrintf(PETSC_COMM_SELF, "\n--- Field Ranges: [%s] (Layout: %s) ---\n", fieldName, data_layout);
2398 }
2399
2400 // --- 4. Branch on DoF and perform calculation with correct bounds ---
2401 if (dof == 1) {
2402 PetscReal localMin = PETSC_MAX_REAL, localMax = PETSC_MIN_REAL;
2403 PetscReal globalMin, globalMax;
2404 const PetscScalar ***array;
2405
2406 ierr = DMDAVecGetArrayRead(dm, fieldVec, &array); CHKERRQ(ierr);
2407 for (k = k_start; k < k_end; k++) {
2408 for (j = j_start; j < j_end; j++) {
2409 for (i = i_start; i < i_end; i++) {
2410 localMin = PetscMin(localMin, array[k][j][i]);
2411 localMax = PetscMax(localMax, array[k][j][i]);
2412 }
2413 }
2414 }
2415 ierr = DMDAVecRestoreArrayRead(dm, fieldVec, &array); CHKERRQ(ierr);
2416
2417 ierr = MPI_Allreduce(&localMin, &globalMin, 1, MPIU_REAL, MPI_MIN, PETSC_COMM_WORLD); CHKERRQ(ierr);
2418 ierr = MPI_Allreduce(&localMax, &globalMax, 1, MPIU_REAL, MPI_MAX, PETSC_COMM_WORLD); CHKERRQ(ierr);
2419
2420 PetscSynchronizedPrintf(PETSC_COMM_WORLD, " [Rank %d] Local Range: [ %11.4e , %11.4e ]\n", user->simCtx->rank, localMin, localMax);
2421 ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT); CHKERRQ(ierr);
2422 if (user->simCtx->rank == 0) {
2423 PetscPrintf(PETSC_COMM_SELF, " Global Range: [ %11.4e , %11.4e ]\n", globalMin, globalMax);
2424 }
2425
2426 } else if (dof == 3) {
2427 Cmpnts localMin = {PETSC_MAX_REAL, PETSC_MAX_REAL, PETSC_MAX_REAL};
2428 Cmpnts localMax = {PETSC_MIN_REAL, PETSC_MIN_REAL, PETSC_MIN_REAL};
2429 Cmpnts globalMin, globalMax;
2430 const Cmpnts ***array;
2431
2432 ierr = DMDAVecGetArrayRead(dm, fieldVec, &array); CHKERRQ(ierr);
2433 for (k = k_start; k < k_end; k++) {
2434 for (j = j_start; j < j_end; j++) {
2435 for (i = i_start; i < i_end; i++) {
2436 localMin.x = PetscMin(localMin.x, array[k][j][i].x);
2437 localMin.y = PetscMin(localMin.y, array[k][j][i].y);
2438 localMin.z = PetscMin(localMin.z, array[k][j][i].z);
2439 localMax.x = PetscMax(localMax.x, array[k][j][i].x);
2440 localMax.y = PetscMax(localMax.y, array[k][j][i].y);
2441 localMax.z = PetscMax(localMax.z, array[k][j][i].z);
2442 }
2443 }
2444 }
2445 ierr = DMDAVecRestoreArrayRead(dm, fieldVec, &array); CHKERRQ(ierr);
2446
2447 ierr = MPI_Allreduce(&localMin, &globalMin, 3, MPIU_REAL, MPI_MIN, PETSC_COMM_WORLD); CHKERRQ(ierr);
2448 ierr = MPI_Allreduce(&localMax, &globalMax, 3, MPIU_REAL, MPI_MAX, PETSC_COMM_WORLD); CHKERRQ(ierr);
2449
2450 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, " [Rank %d] Local X-Range: [ %11.4e , %11.4e ]\n", user->simCtx->rank, localMin.x, localMax.x);
2451 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, " [Rank %d] Local Y-Range: [ %11.4e , %11.4e ]\n", user->simCtx->rank, localMin.y, localMax.y);
2452 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, " [Rank %d] Local Z-Range: [ %11.4e , %11.4e ]\n", user->simCtx->rank, localMin.z, localMax.z);
2453 ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT); CHKERRQ(ierr);
2454
2455 if (user->simCtx->rank == 0) {
2456 PetscPrintf(PETSC_COMM_SELF, " [Global] X-Range: [ %11.4e , %11.4e ]\n", globalMin.x, globalMax.x);
2457 PetscPrintf(PETSC_COMM_SELF, " [Global] Y-Range: [ %11.4e , %11.4e ]\n", globalMin.y, globalMax.y);
2458 PetscPrintf(PETSC_COMM_SELF, " [Global] Z-Range: [ %11.4e , %11.4e ]\n", globalMin.z, globalMax.z);
2459 }
2460
2461 } else {
2462 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "LogFieldStatistics only supports fields with 1 or 3 components, but field '%s' has %" PetscInt_FMT ".", fieldName, dof);
2463 }
2464
2465 // --- 5. Final barrier for clean output ordering ---
2466 ierr = MPI_Barrier(PETSC_COMM_WORLD); CHKERRQ(ierr);
2467 if (user->simCtx->rank == 0) {
2468 PetscPrintf(PETSC_COMM_SELF, "--------------------------------------------\n\n");
2469 }
2470
2471 PetscFunctionReturn(0);
2472}
2473
2474#undef __FUNCT__
2475#define __FUNCT__ "LogFieldAnatomyView"
2476/**
2477 * @brief Shared architecture-aware anatomy logger for catalog and transient fields.
2478 */
2479static PetscErrorCode LogFieldAnatomyView(UserCtx *user, const char *field_name,
2480 const char *stage_name, DM dm,
2481 Vec vec_local, PetscInt dof,
2482 FieldLayout layout)
2483{
2484 PetscErrorCode ierr;
2485 DMDALocalInfo info;
2486 PetscMPIInt rank;
2487 const char *data_layout = FieldLayoutName(layout);
2488 char dominant_dir = '\0'; // 'x', 'y', 'z' for face-centered, 'm' for mixed (Ucont)
2489
2490 PetscFunctionBeginUser;
2491 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
2492
2493 PetscCheck(user != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "UserCtx cannot be NULL.");
2494 PetscCheck(field_name != NULL && stage_name != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
2495 "Field and stage labels cannot be NULL.");
2496 PetscCheck(dm != NULL && vec_local != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE,
2497 "Field '%s' has no available DM/local vector for anatomy logging.", field_name);
2498 PetscCheck(dof == 1 || dof == 3, PETSC_COMM_SELF, PETSC_ERR_SUP,
2499 "Field anatomy logging supports one- or three-component fields; '%s' has %d.",
2500 field_name, dof);
2501
2502 if (layout == FIELD_LAYOUT_I_FACE) dominant_dir = 'x';
2503 else if (layout == FIELD_LAYOUT_J_FACE) dominant_dir = 'y';
2504 else if (layout == FIELD_LAYOUT_K_FACE) dominant_dir = 'z';
2505 else if (layout == FIELD_LAYOUT_COMPONENT_STAGGERED) dominant_dir = 'm';
2506
2507 // --- 2. Get Grid Info and Array Pointers ---
2508 ierr = DMDAGetLocalInfo(dm, &info); CHKERRQ(ierr);
2509
2510 ierr = PetscBarrier(NULL);
2511 PetscPrintf(PETSC_COMM_WORLD, "\n--- Field Anatomy Log: [%s] | Stage: [%s] | Layout: [%s] ---\n", field_name, stage_name, data_layout);
2512
2513 // Global physical dimensions (number of cells)
2514 PetscInt im_phys = user->IM;
2515 PetscInt jm_phys = user->JM;
2516 PetscInt km_phys = user->KM;
2517
2518 // Slice through the center of the local domain
2519 PetscInt i_mid = (PetscInt)(info.xs + 0.5 * info.xm) - 1;
2520 PetscInt j_mid = (PetscInt)(info.ys + 0.5 * info.ym) - 1;
2521 PetscInt k_mid = (PetscInt)(info.zs + 0.5 * info.zm) - 1;
2522
2523 // --- 3. Print Boundary Information based on Data Layout ---
2524
2525 // ======================================================================
2526 // === CASE 1: Cell-Centered Fields (Ucat, P) - USES SHIFTED INDEX ===
2527 // ======================================================================
2528 if (layout == FIELD_LAYOUT_CELL_CENTERED) {
2529 const void *l_arr;
2530 ierr = DMDAVecGetArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2531
2532
2533 // --- I-Direction Boundaries ---
2534 if (info.xs == 0) { // Rank on -Xi boundary
2535 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Ghost for Cell[k][j][0]) = ", rank, 0);
2536 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][j_mid][0]);
2537 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][j_mid][0].x, ((const Cmpnts***)l_arr)[k_mid][j_mid][0].y, ((const Cmpnts***)l_arr)[k_mid][j_mid][0].z);
2538
2539 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Value for Cell[k][j][0]) = ", rank, 1);
2540 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][j_mid][1]);
2541 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][j_mid][1].x, ((const Cmpnts***)l_arr)[k_mid][j_mid][1].y, ((const Cmpnts***)l_arr)[k_mid][j_mid][1].z);
2542 }
2543 if (info.xs + info.xm == info.mx) { // Rank on +Xi boundary
2544 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Value for Cell[k][j][%d]) = ", rank, im_phys - 1, im_phys - 2);
2545 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][j_mid][im_phys - 1]);
2546 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][j_mid][im_phys - 1].x, ((const Cmpnts***)l_arr)[k_mid][j_mid][im_phys - 1].y, ((const Cmpnts***)l_arr)[k_mid][j_mid][im_phys - 1].z);
2547
2548 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Ghost for Cell[k][j][%d]) = ", rank, im_phys, im_phys - 2);
2549 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][j_mid][im_phys]);
2550 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][j_mid][im_phys].x, ((const Cmpnts***)l_arr)[k_mid][j_mid][im_phys].y, ((const Cmpnts***)l_arr)[k_mid][j_mid][im_phys].z);
2551 }
2552
2553 // --- J-Direction Boundaries ---
2554 if (info.ys == 0) { // Rank on -Eta boundary
2555 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Ghost for Cell[k][0][i]) = ", rank, 0);
2556 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][0][i_mid]);
2557 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][0][i_mid].x, ((const Cmpnts***)l_arr)[k_mid][0][i_mid].y, ((const Cmpnts***)l_arr)[k_mid][0][i_mid].z);
2558
2559 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Value for Cell[k][0][i]) = ", rank, 1);
2560 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][1][i_mid]);
2561 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][1][i_mid].x, ((const Cmpnts***)l_arr)[k_mid][1][i_mid].y, ((const Cmpnts***)l_arr)[k_mid][1][i_mid].z);
2562 }
2563
2564 if (info.ys + info.ym == info.my) { // Rank on +Eta boundary
2565 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Value for Cell[k][%d][i]) = ", rank, jm_phys - 1, jm_phys - 2);
2566 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][jm_phys - 1][i_mid]);
2567 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][jm_phys - 1][i_mid].x, ((const Cmpnts***)l_arr)[k_mid][jm_phys - 1][i_mid].y, ((const Cmpnts***)l_arr)[k_mid][jm_phys - 1][i_mid].z);
2568
2569 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Ghost for Cell[k][%d][i]) = ", rank, jm_phys, jm_phys - 2);
2570 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[k_mid][jm_phys][i_mid]);
2571 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[k_mid][jm_phys][i_mid].x, ((const Cmpnts***)l_arr)[k_mid][jm_phys][i_mid].y, ((const Cmpnts***)l_arr)[k_mid][jm_phys][i_mid].z);
2572 }
2573
2574 // --- K-Direction Boundaries ---
2575 if (info.zs == 0) { // Rank on -Zeta boundary
2576 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Ghost for Cell[0][j][i]) = ", rank, 0);
2577 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[0][j_mid][i_mid]);
2578 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[0][j_mid][i_mid].x, ((const Cmpnts***)l_arr)[0][j_mid][i_mid].y, ((const Cmpnts***)l_arr)[0][j_mid][i_mid].z);
2579 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Value for Cell[0][j][i]) = ", rank, 1);
2580 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[1][j_mid][i_mid]);
2581 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[1][j_mid][i_mid].x, ((const Cmpnts***)l_arr)[1][j_mid][i_mid].y, ((const Cmpnts***)l_arr)[1][j_mid][i_mid].z);
2582 }
2583 if (info.zs + info.zm == info.mz) { // Rank on +Zeta boundary
2584 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Value for Cell[%d][j][i]) = ", rank, km_phys - 1, km_phys - 2);
2585 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[km_phys - 1][j_mid][i_mid]);
2586 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[km_phys - 1][j_mid][i_mid].x, ((const Cmpnts***)l_arr)[km_phys - 1][j_mid][i_mid].y, ((const Cmpnts***)l_arr)[km_phys - 1][j_mid][i_mid].z);
2587 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Ghost for Cell[%d][j][i]) = ", rank, km_phys, km_phys - 2);
2588 if(dof==1) PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f)\n", ((const PetscReal***)l_arr)[km_phys][j_mid][i_mid]);
2589 else PetscSynchronizedPrintf(PETSC_COMM_WORLD, "(%.5f, %.5f, %.5f)\n", ((const Cmpnts***)l_arr)[km_phys][j_mid][i_mid].x, ((const Cmpnts***)l_arr)[km_phys][j_mid][i_mid].y, ((const Cmpnts***)l_arr)[km_phys][j_mid][i_mid].z);
2590 }
2591 ierr = DMDAVecRestoreArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2592 }
2593 // ======================================================================
2594 // === CASE 2: Face-Centered Fields - NUANCED DIRECTIONAL LOGIC ===
2595 // ======================================================================
2596 else if (layout == FIELD_LAYOUT_I_FACE ||
2597 layout == FIELD_LAYOUT_J_FACE ||
2598 layout == FIELD_LAYOUT_K_FACE ||
2599 layout == FIELD_LAYOUT_COMPONENT_STAGGERED) {
2600 const Cmpnts ***l_arr;
2601 ierr = DMDAVecGetArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2602
2603 // --- I-Direction Boundaries ---
2604 if (info.xs == 0) { // Rank on -Xi boundary
2605 if (dominant_dir == 'x') { // Node-like in I-dir
2606 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (First Phys. X-Face) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[k_mid][j_mid][0].x, l_arr[k_mid][j_mid][0].y, l_arr[k_mid][j_mid][0].z);
2607 } else if (dominant_dir == 'y' || dominant_dir == 'z') { // Cell-like in I-dir
2608 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Ghost for Cell[k][j][0]) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[k_mid][j_mid][0].x, l_arr[k_mid][j_mid][0].y, l_arr[k_mid][j_mid][0].z);
2609 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Value for Cell[k][j][0]) = (%.5f, %.5f, %.5f)\n", rank, 1, l_arr[k_mid][j_mid][1].x, l_arr[k_mid][j_mid][1].y, l_arr[k_mid][j_mid][1].z);
2610 } else if (dominant_dir == 'm') { // Ucont: Mixed
2611 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: u-comp @ Idx %2d (1st X-Face) = %.5f\n", rank, 0, l_arr[k_mid][j_mid][0].x);
2612 }
2613 }
2614 if (info.xs + info.xm == info.mx) { // Rank on +Xi boundary
2615 if (dominant_dir == 'x') { // Node-like in I-dir
2616 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Last Phys. X-Face) = (%.5f, %.5f, %.5f)\n", rank, im_phys - 1, l_arr[k_mid][j_mid][im_phys - 1].x, l_arr[k_mid][j_mid][im_phys-1].y, l_arr[k_mid][j_mid][im_phys - 1].z);
2617 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Ghost Location) = (%.5f, %.5f, %.5f)\n", rank, im_phys, l_arr[k_mid][j_mid][im_phys].x, l_arr[k_mid][j_mid][im_phys].y, l_arr[k_mid][j_mid][im_phys].z);
2618 } else if (dominant_dir == 'y' || dominant_dir == 'z') { // Cell-like in I-dir
2619 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Value for Cell[k][j][%d]) = (%.5f, %.5f, %.5f)\n", rank, im_phys - 1, im_phys - 2, l_arr[k_mid][j_mid][im_phys - 1].x, l_arr[k_mid][j_mid][im_phys - 1].y, l_arr[k_mid][j_mid][im_phys-1].z);
2620 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Ghost for Cell[k][j][%d]) = (%.5f, %.5f, %.5f)\n", rank, im_phys, im_phys - 2, l_arr[k_mid][j_mid][im_phys].x, l_arr[k_mid][j_mid][im_phys].y, l_arr[k_mid][j_mid][im_phys].z);
2621 } else if (dominant_dir == 'm') { // Ucont: Mixed
2622 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: u-comp @ Idx %2d (Last X-Face) = %.5f\n", rank, im_phys - 1, l_arr[k_mid][j_mid][im_phys - 1].x);
2623 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: u-comp @ Idx %2d (Ghost Location) = %.5f\n", rank, im_phys, l_arr[k_mid][j_mid][im_phys].x);
2624 }
2625 }
2626
2627 // --- J-Direction Boundaries ---
2628 if (info.ys == 0) { // Rank on -Eta boundary
2629 if (dominant_dir == 'y') { // Node-like in J-dir
2630 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (First Phys. Y-Face) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[k_mid][0][i_mid].x, l_arr[k_mid][0][i_mid].y, l_arr[k_mid][0][i_mid].z);
2631 } else if (dominant_dir == 'x' || dominant_dir == 'z') { // Cell-like in J-dir
2632 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Ghost for Cell[k][0][i]) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[k_mid][0][i_mid].x, l_arr[k_mid][0][i_mid].y, l_arr[k_mid][0][i_mid].z);
2633 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Value for Cell[k][0][i]) = (%.5f, %.5f, %.5f)\n", rank, 1, l_arr[k_mid][1][i_mid].x, l_arr[k_mid][1][i_mid].y, l_arr[k_mid][1][i_mid].z);
2634 } else if (dominant_dir == 'm') { // Ucont: Mixed
2635 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: v-comp @ Jdx %2d (1st Y-Face) = %.5f\n", rank, 0, l_arr[k_mid][0][i_mid].y);
2636 }
2637 }
2638 if (info.ys + info.ym == info.my) { // Rank on +Eta boundary
2639 if (dominant_dir == 'y') { // Node-like in J-dir
2640 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Last Phys. Y-Face) = (%.5f, %.5f, %.5f)\n", rank, jm_phys - 1, l_arr[k_mid][jm_phys - 1][i_mid].x, l_arr[k_mid][jm_phys - 1][i_mid].y, l_arr[k_mid][jm_phys - 1][i_mid].z);
2641 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Ghost Location) = (%.5f, %.5f, %.5f)\n", rank, jm_phys, l_arr[k_mid][jm_phys][i_mid].x, l_arr[k_mid][jm_phys][i_mid].y, l_arr[k_mid][jm_phys][i_mid].z);
2642 } else if (dominant_dir == 'x' || dominant_dir == 'z') { // Cell-like in J-dir
2643 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Value for Cell[k][%d][i]) = (%.5f, %.5f, %.5f)\n", rank, jm_phys-1, jm_phys-2, l_arr[k_mid][jm_phys - 1][i_mid].x, l_arr[k_mid][jm_phys - 1][i_mid].y, l_arr[k_mid][jm_phys - 1][i_mid].z);
2644 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Ghost for Cell[k][%d][i]) = (%.5f, %.5f, %.5f)\n", rank, jm_phys, jm_phys-2, l_arr[k_mid][jm_phys][i_mid].x, l_arr[k_mid][jm_phys][i_mid].y, l_arr[k_mid][jm_phys][i_mid].z);
2645 } else if (dominant_dir == 'm') { // Ucont: Mixed
2646 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: v-comp @ Jdx %2d (Last Y-Face) = %.5f\n", rank, jm_phys - 1, l_arr[k_mid][jm_phys - 1][i_mid].y);
2647 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: v-comp @ Jdx %2d (Ghost Location) = %.5f\n", rank, jm_phys, l_arr[k_mid][jm_phys][i_mid].y);
2648 }
2649 }
2650
2651 // --- K-Direction Boundaries ---
2652 if (info.zs == 0) { // Rank on -Zeta boundary
2653 if (dominant_dir == 'z') { // Node-like in K-dir
2654 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (First Phys. Z-Face) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[0][j_mid][i_mid].x, l_arr[0][j_mid][i_mid].y, l_arr[0][j_mid][i_mid].z);
2655 } else if (dominant_dir == 'x' || dominant_dir == 'y') { // Cell-like in K-dir
2656 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Ghost for Cell[0][j][i]) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[0][j_mid][i_mid].x, l_arr[0][j_mid][i_mid].y, l_arr[0][j_mid][i_mid].z);
2657 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Value for Cell[0][j][i]) = (%.5f, %.5f, %.5f)\n", rank, 1, l_arr[1][j_mid][i_mid].x, l_arr[1][j_mid][i_mid].y, l_arr[1][j_mid][i_mid].z);
2658 } else if (dominant_dir == 'm') { // Ucont: Mixed
2659 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: w-comp @ Idx %2d (1st Z-Face) = %.5f\n", rank, 0, l_arr[0][j_mid][i_mid].z);
2660 }
2661 }
2662 if (info.zs + info.zm == info.mz) { // Rank on +Zeta boundary
2663 if (dominant_dir == 'z') { // Node-like in K-dir
2664 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Idx %2d (Last Phys. Z-Face) = (%.5f, %.5f, %.5f)\n", rank, km_phys - 1, l_arr[km_phys - 1][j_mid][i_mid].x, l_arr[km_phys - 1][j_mid][i_mid].y, l_arr[km_phys - 1][j_mid][i_mid].z);
2665 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Idx %2d (Ghost Location) = (%.5f, %.5f, %.5f)\n", rank, km_phys, l_arr[km_phys][j_mid][i_mid].x, l_arr[km_phys][j_mid][i_mid].y, l_arr[km_phys][j_mid][i_mid].z);
2666 } else if (dominant_dir == 'x' || dominant_dir == 'y') { // Cell-like in K-dir
2667 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Idx %2d (Value for Cell[%d][j][i]) = (%.5f, %.5f, %.5f)\n", rank, km_phys-1, km_phys-2, l_arr[km_phys-1][j_mid][i_mid].x, l_arr[km_phys-1][j_mid][i_mid].y, l_arr[km_phys - 1][j_mid][i_mid].z);
2668 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Idx %2d (Ghost for Cell[%d][j][i]) = (%.5f, %.5f, %.5f)\n", rank, km_phys, km_phys-2, l_arr[km_phys][j_mid][i_mid].x, l_arr[km_phys][j_mid][i_mid].y, l_arr[km_phys][j_mid][i_mid].z);
2669 } else if (dominant_dir == 'm') { // Ucont: Mixed
2670 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: w-comp @ Idx %2d (Last Z-Face) = %.5f\n", rank, km_phys - 1, l_arr[km_phys - 1][j_mid][i_mid].z);
2671 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: w-comp @ Idx %2d (Ghost Loc.) = %.5f\n", rank, km_phys, l_arr[km_phys][j_mid][i_mid].z);
2672
2673 }
2674 }
2675 ierr = DMDAVecRestoreArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2676 }
2677 // ======================================================================
2678 // === CASE 3: Node-Centered Fields - USES DIRECT INDEX ===
2679 // ======================================================================
2680 else if (layout == FIELD_LAYOUT_NODE_CENTERED) {
2681 if (dof == 1) {
2682 const PetscReal ***l_arr;
2683 ierr = DMDAVecGetArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2684
2685 if (info.xs == 0)
2686 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (First Phys. Node) = %.5f\n", rank, 0, l_arr[k_mid][j_mid][0]);
2687 if (info.xs + info.xm == info.mx) {
2688 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Last Phys. Node) = %.5f\n", rank, im_phys - 1, l_arr[k_mid][j_mid][im_phys - 1]);
2689 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Unused/Ghost Loc) = %.5f\n", rank, im_phys, l_arr[k_mid][j_mid][im_phys]);
2690 }
2691 if (info.ys == 0)
2692 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (First Phys. Node) = %.5f\n", rank, 0, l_arr[k_mid][0][i_mid]);
2693 if (info.ys + info.ym == info.my) {
2694 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Last Phys. Node) = %.5f\n", rank, jm_phys - 1, l_arr[k_mid][jm_phys - 1][i_mid]);
2695 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Unused/Ghost Loc) = %.5f\n", rank, jm_phys, l_arr[k_mid][jm_phys][i_mid]);
2696 }
2697 if (info.zs == 0)
2698 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (First Phys. Node) = %.5f\n", rank, 0, l_arr[0][j_mid][i_mid]);
2699 if (info.zs + info.zm == info.mz) {
2700 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Last Phys. Node) = %.5f\n", rank, km_phys - 1, l_arr[km_phys - 1][j_mid][i_mid]);
2701 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Unused/Ghost Loc) = %.5f\n", rank, km_phys, l_arr[km_phys][j_mid][i_mid]);
2702 }
2703 ierr = DMDAVecRestoreArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2704 } else {
2705 const Cmpnts ***l_arr;
2706 ierr = DMDAVecGetArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2707
2708 // --- I-Direction Boundaries ---
2709 if (info.xs == 0) {
2710 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (First Phys. Node) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[k_mid][j_mid][0].x, l_arr[k_mid][j_mid][0].y, l_arr[k_mid][j_mid][0].z);
2711 }
2712 if (info.xs + info.xm == info.mx) {
2713 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Last Phys. Node) = (%.5f, %.5f, %.5f)\n", rank, im_phys - 1, l_arr[k_mid][j_mid][im_phys - 1].x, l_arr[k_mid][j_mid][im_phys - 1].y, l_arr[k_mid][j_mid][im_phys - 1].z);
2714 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, I-DIR]: Idx %2d (Unused/Ghost Loc) = (%.5f, %.5f, %.5f)\n", rank, im_phys, l_arr[k_mid][j_mid][im_phys].x, l_arr[k_mid][j_mid][im_phys].y, l_arr[k_mid][j_mid][im_phys].z);
2715 }
2716 // --- J-Direction Boundaries ---
2717 if (info.ys == 0) {
2718 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (First Phys. Node) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[k_mid][0][i_mid].x, l_arr[k_mid][0][i_mid].y, l_arr[k_mid][0][i_mid].z);
2719 }
2720 if (info.ys + info.ym == info.my) {
2721 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Last Phys. Node) = (%.5f, %.5f, %.5f)\n", rank, jm_phys - 1, l_arr[k_mid][jm_phys - 1][i_mid].x, l_arr[k_mid][jm_phys - 1][i_mid].y, l_arr[k_mid][jm_phys - 1][i_mid].z);
2722 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, J-DIR]: Jdx %2d (Unused/Ghost Loc) = (%.5f, %.5f, %.5f)\n", rank, jm_phys, l_arr[k_mid][jm_phys][i_mid].x, l_arr[k_mid][jm_phys][i_mid].y, l_arr[k_mid][jm_phys][i_mid].z);
2723 }
2724 // --- K-Direction Boundaries ---
2725 if (info.zs == 0) {
2726 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (First Phys. Node) = (%.5f, %.5f, %.5f)\n", rank, 0, l_arr[0][j_mid][i_mid].x, l_arr[0][j_mid][i_mid].y, l_arr[0][j_mid][i_mid].z);
2727 }
2728 if(info.zs + info.zm == info.mz) {
2729 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Last Phys. Node) = (%.5f, %.5f, %.5f)\n", rank, km_phys - 1, l_arr[km_phys - 1][j_mid][i_mid].x, l_arr[km_phys - 1][j_mid][i_mid].y, l_arr[km_phys - 1][j_mid][i_mid].z);
2730 PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[Rank %d, K-DIR]: Kdx %2d (Unused/Ghost Loc) = (%.5f, %.5f, %.5f)\n", rank, km_phys, l_arr[km_phys][j_mid][i_mid].x, l_arr[km_phys][j_mid][i_mid].y, l_arr[km_phys][j_mid][i_mid].z);
2731 }
2732 ierr = DMDAVecRestoreArrayRead(dm, vec_local, (void*)&l_arr); CHKERRQ(ierr);
2733 }
2734 }
2735 else {
2736 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG,
2737 "LOG_FIELD_ANATOMY encountered unsupported layout %d for field '%s'.",
2738 (int)layout, field_name);
2739 }
2740
2741 ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT); CHKERRQ(ierr);
2742 ierr = PetscBarrier(NULL);
2743 PetscFunctionReturn(0);
2744}
2745
2746#undef __FUNCT__
2747#define __FUNCT__ "LOG_FIELD_ANATOMY"
2748/**
2749 * @brief Resolves a persistent field view and emits its layout-aware boundary anatomy.
2750 * @see LOG_FIELD_ANATOMY()
2751 */
2752PetscErrorCode LOG_FIELD_ANATOMY(UserCtx *user, FieldId field_id, const char *stage_name)
2753{
2754 FieldView view;
2755
2756 PetscFunctionBeginUser;
2757 PetscCall(FieldGetView(user, field_id, &view));
2758 PetscCall(LogFieldAnatomyView(user, view.descriptor->canonical_name, stage_name,
2759 view.dm, view.local_vec, view.descriptor->dof,
2760 view.descriptor->layout));
2761 PetscFunctionReturn(0);
2762}
2763
2764/**
2765 * @brief Implementation of \ref IsStatisticsConsoleSnapshotEnabled().
2766 * @see IsStatisticsConsoleSnapshotEnabled()
2767 */
2768PetscBool IsStatisticsConsoleSnapshotEnabled(const SimCtx *simCtx)
2769{
2770 if (!FieldStatisticsIsActive(simCtx)) return PETSC_FALSE;
2771 return (PetscBool)(simCtx->statisticsConsoleOutputFreq > 0 && get_log_level() >= LOG_INFO);
2772}
2773
2774/**
2775 * @brief Implementation of \ref ShouldEmitPeriodicStatisticsConsoleSnapshot().
2776 * @see ShouldEmitPeriodicStatisticsConsoleSnapshot()
2777 */
2778PetscBool ShouldEmitPeriodicStatisticsConsoleSnapshot(const SimCtx *simCtx, PetscInt completed_step)
2779{
2780 if (!IsStatisticsConsoleSnapshotEnabled(simCtx)) return PETSC_FALSE;
2781 if (completed_step < 0) return PETSC_FALSE;
2782 return (PetscBool)((completed_step % simCtx->statisticsConsoleOutputFreq) == 0);
2783}
2784
2785/**
2786 * @brief Implementation of \ref EmitStatisticsConsoleSnapshot().
2787 * @see EmitStatisticsConsoleSnapshot()
2788 */
2789PetscErrorCode EmitStatisticsConsoleSnapshot(UserCtx *user, const SimCtx *simCtx, PetscInt step)
2790{
2791 PetscFunctionBeginUser;
2792 PetscCheck(simCtx != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "SimCtx cannot be NULL.");
2793 if (!IsStatisticsConsoleSnapshotEnabled(simCtx)) PetscFunctionReturn(0);
2794
2795 LOG(GLOBAL, LOG_INFO, "Statistics windows at step %d (%d window(s)):\n",
2796 step, simCtx->fieldStatisticsWindowCount);
2797 for (PetscInt w = 0; w < simCtx->fieldStatisticsWindowCount; ++w) {
2798 const PicurvWindow *window = &simCtx->fieldStatisticsWindows[w];
2799 char progress[16];
2800 char coverage[32] = "";
2801
2802 if (window->definition.bounded) {
2803 PetscCall(PetscSNPrintf(progress, sizeof(progress), "%5.1f%%",
2804 100.0 * (double)PicurvWindowProgress(window)));
2805 } else {
2806 PetscCall(PetscSNPrintf(progress, sizeof(progress), " open"));
2807 }
2808
2809 /* Mask health: with a moving body different points see different numbers of
2810 * states, and the spread of per-point valid fraction is what tells an
2811 * operator whether the window covers the domain evenly. Reported only on the
2812 * first block, since the reduction it performs is already collective. */
2813 if (user && user->fieldStatisticsStorage && window->sample_count > 0) {
2814 PetscReal lowest = 1.0, highest = 0.0;
2815
2816 PetscCall(PicurvWindowValidFractionRange(user, &window->definition,
2817 &user->fieldStatisticsStorage[w],
2818 window->sample_count, &lowest, &highest));
2819 PetscCall(PetscSNPrintf(coverage, sizeof(coverage), " valid=[%.2f,%.2f]",
2820 (double)lowest, (double)highest));
2821 }
2822
2823 LOG(GLOBAL, LOG_INFO,
2824 " %-24s %-8s samples=%-8d weight=%-12.6g represented=%-12.6g progress=%s%s\n",
2825 window->definition.name, PicurvWindowStateName(window->state),
2826 window->sample_count, (double)window->total_weight,
2827 (double)window->represented_time, progress, coverage);
2828 }
2829 PetscFunctionReturn(0);
2830}
2831
2832#undef __FUNCT__
2833#define __FUNCT__ "LOG_CORNER_FIELD_ANATOMY"
2834/**
2835 * @brief Emits anatomy for the transient scalar or vector corner-staging field.
2836 * @see LOG_CORNER_FIELD_ANATOMY()
2837 */
2838PetscErrorCode LOG_CORNER_FIELD_ANATOMY(UserCtx *user, FieldId corner_field_id, const char *stage_name)
2839{
2840 FieldView view;
2841
2842 PetscFunctionBeginUser;
2843 PetscCheck(user != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "UserCtx cannot be NULL.");
2844 PetscCheck(corner_field_id == FIELD_ID_CELL_SCALAR_AT_CORNER ||
2845 corner_field_id == FIELD_ID_CELL_VECTOR_AT_CORNER,
2846 PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE,
2847 "Corner anatomy logging expects a corner-staging field identity.");
2848 /* The caller states which workspace it used, so the degree of freedom and DM
2849 * come from the catalog instead of being inferred from a cached vector. */
2850 PetscCall(FieldGetView(user, corner_field_id, &view));
2851 PetscCall(LogFieldAnatomyView(user, view.descriptor->canonical_name, stage_name, view.dm,
2852 view.local_vec, view.descriptor->dof,
2853 view.descriptor->layout));
2854 PetscFunctionReturn(0);
2855}
2856
2857#undef __FUNCT__
2858#define __FUNCT__ "LOG_INTERPOLATION_ERROR"
2859/**
2860 * @brief Implementation of \ref LOG_INTERPOLATION_ERROR().
2861 * @details Full API contract (arguments, ownership, side effects) is documented with
2862 * the header declaration in `include/logging.h`.
2863 * @see LOG_INTERPOLATION_ERROR()
2864 */
2865PetscErrorCode LOG_INTERPOLATION_ERROR(UserCtx *user)
2866{
2867 SimCtx *simCtx = user->simCtx;
2868 PetscErrorCode ierr;
2869 DM swarm = user->swarm;
2870 Vec positionVec, analyticalvelocityVec, velocityVec, errorVec;
2871 PetscReal Interpolation_error = 0.0;
2872 PetscReal Maximum_Interpolation_error = 0.0;
2873 PetscReal AnalyticalSolution_magnitude = 0.0;
2874 PetscReal ErrorPercentage = 0.0;
2875
2876 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Creating global vectors.\n");
2877 ierr = DMSwarmCreateGlobalVectorFromField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_POSITION), &positionVec); CHKERRQ(ierr);
2878 ierr = DMSwarmCreateGlobalVectorFromField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_VELOCITY), &velocityVec); CHKERRQ(ierr);
2879
2880 ierr = VecDuplicate(positionVec, &analyticalvelocityVec); CHKERRQ(ierr);
2881 ierr = VecCopy(positionVec, analyticalvelocityVec); CHKERRQ(ierr);
2882
2883 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Computing analytical solution.\n");
2884 ierr = SetAnalyticalSolutionForParticles(analyticalvelocityVec, simCtx); CHKERRQ(ierr);
2885
2886 ierr = VecDuplicate(analyticalvelocityVec, &errorVec); CHKERRQ(ierr);
2887 ierr = VecCopy(analyticalvelocityVec, errorVec); CHKERRQ(ierr);
2888
2889 ierr = VecNorm(analyticalvelocityVec, NORM_2, &AnalyticalSolution_magnitude); CHKERRQ(ierr);
2890
2891 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Computing error.\n");
2892 ierr = VecAXPY(errorVec, -1.0, velocityVec); CHKERRQ(ierr);
2893 ierr = VecNorm(errorVec, NORM_2, &Interpolation_error); CHKERRQ(ierr);
2894 ierr = VecNorm(errorVec,NORM_INFINITY,&Maximum_Interpolation_error); CHKERRQ(ierr);
2895
2896 ErrorPercentage = (AnalyticalSolution_magnitude > 0) ?
2897 (Interpolation_error / AnalyticalSolution_magnitude * 100.0) : 0.0;
2898
2899 /* --- CSV output (always, rank 0 only) --- */
2900 if (simCtx->rank == 0) {
2901 char csv_path[PETSC_MAX_PATH_LEN + 32];
2902 ierr = PetscSNPrintf(csv_path, sizeof(csv_path), "%s/interpolation_error.csv", simCtx->log_dir); CHKERRQ(ierr);
2903 FILE *f = fopen(csv_path, "a");
2904 if (f) {
2905 if (ftell(f) == 0) {
2906 fprintf(f, "step,time,L2_error,Linf_error,L2_analytical,error_pct\n");
2907 }
2908 if (simCtx->continueMode && simCtx->step == simCtx->StartStep + 1) {
2909 fprintf(f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
2910 }
2911 PetscReal t = (PetscReal)simCtx->ti * simCtx->dt;
2912 fprintf(f, "%d,%.6e,%.6e,%.6e,%.6e,%.4f\n",
2913 (int)simCtx->step, t,
2914 Interpolation_error, Maximum_Interpolation_error,
2915 AnalyticalSolution_magnitude, ErrorPercentage);
2916 fclose(f);
2917 }
2918 }
2919
2920 /* --- Console output (only at INFO level or above) --- */
2921 if (get_log_level() >= LOG_INFO) {
2922 LOG_ALLOW(GLOBAL, LOG_INFO, "Interpolation error (%%): %g\n", ErrorPercentage);
2923 PetscPrintf(PETSC_COMM_WORLD, "Interpolation error (%%): %g\n", ErrorPercentage);
2924 LOG_ALLOW(GLOBAL, LOG_INFO, "Maximum Interpolation error: %g\n", Maximum_Interpolation_error);
2925 PetscPrintf(PETSC_COMM_WORLD, "Maximum Interpolation error: %g\n", Maximum_Interpolation_error);
2926 }
2927
2928 ierr = VecDestroy(&analyticalvelocityVec); CHKERRQ(ierr);
2929 ierr = VecDestroy(&errorVec); CHKERRQ(ierr);
2930 ierr = DMSwarmDestroyGlobalVectorFromField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_POSITION), &positionVec); CHKERRQ(ierr);
2931 ierr = DMSwarmDestroyGlobalVectorFromField(swarm, ParticleFieldName(PARTICLE_FIELD_ID_VELOCITY), &velocityVec); CHKERRQ(ierr);
2932
2933 return 0;
2934}
2935
2936#undef __FUNCT__
2937#define __FUNCT__ "LOG_SCATTER_METRICS"
2938/**
2939 * @brief Implementation of \ref LOG_SCATTER_METRICS().
2940 * @details Full API contract (arguments, ownership, side effects) is documented with
2941 * the header declaration in `include/logging.h`.
2942 * @see LOG_SCATTER_METRICS()
2943 */
2944PetscErrorCode LOG_SCATTER_METRICS(UserCtx *user)
2945{
2946 PetscErrorCode ierr;
2947 SimCtx *simCtx = NULL;
2948 DMDALocalInfo info;
2949 PetscInt xs, xe, ys, ye, zs, ze, mx, my, mz;
2950 PetscInt lxs, lxe, lys, lye, lzs, lze;
2951 Vec reference_vec = NULL;
2952 PetscReal ***psi = NULL;
2953 PetscReal ***psi_ref = NULL;
2954 PetscReal ***aj = NULL;
2955 PetscReal ***count = NULL;
2956 PetscReal *particle_psi = NULL;
2957 PetscInt nlocal = 0;
2958 PetscReal local_l1 = 0.0, global_l1 = 0.0;
2959 PetscReal local_l2_sq = 0.0, global_l2_sq = 0.0;
2960 PetscReal local_linf = 0.0, global_linf = 0.0;
2961 PetscReal local_ref_l2_sq = 0.0, global_ref_l2_sq = 0.0;
2962 PetscReal local_grid_integral = 0.0, global_grid_integral = 0.0;
2963 PetscReal local_domain_volume = 0.0, global_domain_volume = 0.0;
2964 PetscReal local_particle_sum = 0.0, global_particle_sum = 0.0;
2965 PetscInt64 local_particle_count = 0, global_particle_count = 0;
2966 PetscInt64 local_cell_count = 0, global_cell_count = 0;
2967 PetscInt64 local_occupied_count = 0, global_occupied_count = 0;
2968 PetscReal particle_integral = 0.0;
2969 PetscReal occupancy_fraction = 0.0;
2970 PetscReal mean_particles_per_occupied_cell = 0.0;
2971 PetscReal l2_error = 0.0;
2972 PetscReal relative_l2_error = 0.0;
2973
2974 PetscFunctionBeginUser;
2975 if (!user) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "UserCtx cannot be NULL.");
2976 simCtx = user->simCtx;
2977 if (!VerificationScalarOverrideActive(simCtx) || !user->swarm || !user->Psi || !user->ParticleCount) {
2978 PetscFunctionReturn(0);
2979 }
2980
2981 info = user->info;
2982 xs = info.xs; xe = info.xs + info.xm;
2983 ys = info.ys; ye = info.ys + info.ym;
2984 zs = info.zs; ze = info.zs + info.zm;
2985 mx = info.mx; my = info.my; mz = info.mz;
2986 lxs = (xs == 0) ? xs + 1 : xs; lxe = (xe == mx) ? xe - 1 : xe;
2987 lys = (ys == 0) ? ys + 1 : ys; lye = (ye == my) ? ye - 1 : ye;
2988 lzs = (zs == 0) ? zs + 1 : zs; lze = (ze == mz) ? ze - 1 : ze;
2989
2990 ierr = VecDuplicate(user->Psi, &reference_vec); CHKERRQ(ierr);
2991 ierr = SetAnalyticalScalarFieldAtCellCenters(user, reference_vec); CHKERRQ(ierr);
2992
2993 ierr = DMDAVecGetArrayRead(user->da, user->Psi, &psi); CHKERRQ(ierr);
2994 ierr = DMDAVecGetArrayRead(user->da, reference_vec, &psi_ref); CHKERRQ(ierr);
2995 ierr = DMDAVecGetArrayRead(user->da, user->Aj, &aj); CHKERRQ(ierr);
2996 ierr = DMDAVecGetArrayRead(user->da, user->ParticleCount, &count); CHKERRQ(ierr);
2997
2998 for (PetscInt k = lzs; k < lze; ++k) {
2999 for (PetscInt j = lys; j < lye; ++j) {
3000 for (PetscInt i = lxs; i < lxe; ++i) {
3001 const PetscReal cell_volume = (PetscAbsReal(aj[k][j][i]) > 1.0e-14) ? (1.0 / aj[k][j][i]) : 0.0;
3002 const PetscReal err = psi[k][j][i] - psi_ref[k][j][i];
3003 local_cell_count += 1;
3004 local_domain_volume += cell_volume;
3005 local_grid_integral += psi[k][j][i] * cell_volume;
3006 local_l1 += PetscAbsReal(err) * cell_volume;
3007 local_l2_sq += err * err * cell_volume;
3008 local_ref_l2_sq += psi_ref[k][j][i] * psi_ref[k][j][i] * cell_volume;
3009 local_linf = PetscMax(local_linf, PetscAbsReal(err));
3010 if (count[k][j][i] > 0.0) local_occupied_count += 1;
3011 }
3012 }
3013 }
3014
3015 ierr = DMDAVecRestoreArrayRead(user->da, user->ParticleCount, &count); CHKERRQ(ierr);
3016 ierr = DMDAVecRestoreArrayRead(user->da, user->Aj, &aj); CHKERRQ(ierr);
3017 ierr = DMDAVecRestoreArrayRead(user->da, reference_vec, &psi_ref); CHKERRQ(ierr);
3018 ierr = DMDAVecRestoreArrayRead(user->da, user->Psi, &psi); CHKERRQ(ierr);
3019 ierr = VecDestroy(&reference_vec); CHKERRQ(ierr);
3020
3021 ierr = DMSwarmGetLocalSize(user->swarm, &nlocal); CHKERRQ(ierr);
3022 local_particle_count = (PetscInt64)nlocal;
3023 if (nlocal > 0) {
3024 ierr = DMSwarmGetField(user->swarm, ParticleFieldName(PARTICLE_FIELD_ID_PSI), NULL, NULL, (void **)&particle_psi); CHKERRQ(ierr);
3025 for (PetscInt p = 0; p < nlocal; ++p) local_particle_sum += particle_psi[p];
3026 ierr = DMSwarmRestoreField(user->swarm, ParticleFieldName(PARTICLE_FIELD_ID_PSI), NULL, NULL, (void **)&particle_psi); CHKERRQ(ierr);
3027 }
3028
3029 ierr = MPI_Allreduce(&local_l1, &global_l1, 1, MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3030 ierr = MPI_Allreduce(&local_l2_sq, &global_l2_sq, 1, MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3031 ierr = MPI_Allreduce(&local_linf, &global_linf, 1, MPIU_REAL, MPI_MAX, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3032 ierr = MPI_Allreduce(&local_ref_l2_sq, &global_ref_l2_sq, 1, MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3033 ierr = MPI_Allreduce(&local_grid_integral, &global_grid_integral, 1, MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3034 ierr = MPI_Allreduce(&local_domain_volume, &global_domain_volume, 1, MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3035 ierr = MPI_Allreduce(&local_particle_sum, &global_particle_sum, 1, MPIU_REAL, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3036 ierr = MPI_Allreduce(&local_particle_count, &global_particle_count, 1, MPIU_INT64, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3037 ierr = MPI_Allreduce(&local_cell_count, &global_cell_count, 1, MPIU_INT64, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3038 ierr = MPI_Allreduce(&local_occupied_count, &global_occupied_count, 1, MPIU_INT64, MPI_SUM, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3039
3040 l2_error = PetscSqrtReal(global_l2_sq);
3041 relative_l2_error = (global_ref_l2_sq > 0.0) ? (l2_error / PetscSqrtReal(global_ref_l2_sq)) : 0.0;
3042 occupancy_fraction = (global_cell_count > 0) ? ((PetscReal)global_occupied_count / (PetscReal)global_cell_count) : 0.0;
3043 mean_particles_per_occupied_cell =
3044 (global_occupied_count > 0) ? ((PetscReal)global_particle_count / (PetscReal)global_occupied_count) : 0.0;
3045 particle_integral =
3046 (global_particle_count > 0) ? (global_domain_volume * global_particle_sum / (PetscReal)global_particle_count) : 0.0;
3047
3048 if (simCtx->rank == 0) {
3049 char csv_path[PETSC_MAX_PATH_LEN + 32];
3050 FILE *f = NULL;
3051 ierr = PetscSNPrintf(csv_path, sizeof(csv_path), "%s/scatter_metrics.csv", simCtx->log_dir); CHKERRQ(ierr);
3052 f = fopen(csv_path, "a");
3053 if (f) {
3054 if (ftell(f) == 0) {
3055 fprintf(f,
3056 "step,time,total_particles,total_cells,occupied_cells,occupancy_fraction,"
3057 "mean_particles_per_occupied_cell,particle_integral,grid_integral,"
3058 "conservation_error_abs,L1_error,L2_error,Linf_error,relative_L2_error\n");
3059 }
3060 if (simCtx->continueMode && simCtx->step == simCtx->StartStep + 1) {
3061 fprintf(f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
3062 }
3063 fprintf(f, "%d,%.6e,%lld,%lld,%lld,%.6e,%.6e,%.6e,%.6e,%.6e,%.6e,%.6e,%.6e,%.6e\n",
3064 (int)simCtx->step,
3065 (double)simCtx->ti,
3066 (long long)global_particle_count,
3067 (long long)global_cell_count,
3068 (long long)global_occupied_count,
3069 (double)occupancy_fraction,
3070 (double)mean_particles_per_occupied_cell,
3071 (double)particle_integral,
3072 (double)global_grid_integral,
3073 (double)PetscAbsReal(global_grid_integral - particle_integral),
3074 (double)global_l1,
3075 (double)l2_error,
3076 (double)global_linf,
3077 (double)relative_l2_error);
3078 fclose(f);
3079 }
3080 }
3081
3082 if (get_log_level() >= LOG_INFO) {
3083 LOG_ALLOW(GLOBAL, LOG_INFO, "Scatter relative L2 error: %.6e\n", (double)relative_l2_error);
3084 LOG_ALLOW(GLOBAL, LOG_INFO, "Scatter occupancy fraction: %.6e\n", (double)occupancy_fraction);
3085 }
3086
3087 PetscFunctionReturn(0);
3088}
3089
3090#undef __FUNCT__
3091#define __FUNCT__ "ResetSearchMetrics"
3092/**
3093 * @brief Implementation of \ref ResetSearchMetrics().
3094 * @details Full API contract (arguments, ownership, side effects) is documented with
3095 * the header declaration in `include/logging.h`.
3096 * @see ResetSearchMetrics()
3097 */
3098PetscErrorCode ResetSearchMetrics(SimCtx *simCtx)
3099{
3100 PetscFunctionBeginUser;
3101 if (!simCtx) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "SimCtx cannot be NULL for ResetSearchMetrics.");
3102
3103 simCtx->searchMetrics.searchAttempts = 0;
3104 simCtx->searchMetrics.searchPopulation = 0;
3105 simCtx->searchMetrics.searchLocatedCount = 0;
3106 simCtx->searchMetrics.searchLostCount = 0;
3107 simCtx->searchMetrics.traversalStepsSum = 0;
3108 simCtx->searchMetrics.reSearchCount = 0;
3109 simCtx->searchMetrics.maxTraversalSteps = 0;
3110 simCtx->searchMetrics.maxTraversalFailCount = 0;
3111 simCtx->searchMetrics.tieBreakCount = 0;
3112 simCtx->searchMetrics.boundaryClampCount = 0;
3113 simCtx->searchMetrics.bboxGuessSuccessCount = 0;
3114 simCtx->searchMetrics.bboxGuessFallbackCount = 0;
3115 simCtx->searchMetrics.maxParticlePassDepth = 0;
3116 simCtx->searchMetrics.currentSettlementPass = 0;
3117
3118 PetscFunctionReturn(0);
3119}
3120
3121#undef __FUNCT__
3122#define __FUNCT__ "LOG_SEARCH_METRICS"
3123/**
3124 * @brief Implementation of \ref LOG_SEARCH_METRICS().
3125 * @details Full API contract (arguments, ownership, side effects) is documented with
3126 * the header declaration in `include/logging.h`.
3127 * @see LOG_SEARCH_METRICS()
3128 */
3129PetscErrorCode LOG_SEARCH_METRICS(UserCtx *user)
3130{
3131 PetscErrorCode ierr;
3132 SimCtx *simCtx = NULL;
3133 PetscInt totalParticles = 0;
3134 PetscReal local_metrics[SEARCH_METRIC_REDUCTION_LEN] = {0.0};
3135 PetscReal global_metrics[SEARCH_METRIC_REDUCTION_LEN] = {0.0};
3136 PetscReal meanTraversalSteps = 0.0;
3137 PetscReal searchFailureFraction = 0.0;
3138 PetscReal searchWorkIndex = 0.0;
3139 PetscReal reSearchFraction = 0.0;
3140 long long searchAttempts = 0;
3141 long long searchPopulation = 0;
3142 long long searchLocatedCount = 0;
3143 long long searchLostCount = 0;
3144 long long traversalStepsSum = 0;
3145 long long reSearchCount = 0;
3146 long long tieBreakCount = 0;
3147 long long boundaryClampCount = 0;
3148 long long bboxGuessSuccessCount = 0;
3149 long long bboxGuessFallbackCount = 0;
3150 long long maxTraversalFailCount = 0;
3151 long long maxTraversalSteps = 0;
3152 long long maxPassDepth = 0;
3153 MPI_Op reduction_op = MPI_OP_NULL;
3154
3155 PetscFunctionBeginUser;
3156 if (!user || !user->simCtx) {
3157 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "UserCtx and SimCtx are required for LOG_SEARCH_METRICS.");
3158 }
3159 simCtx = user->simCtx;
3160
3161 if (simCtx->np <= 0) {
3162 PetscFunctionReturn(0);
3163 }
3164
3165 ierr = DMSwarmGetSize(user->swarm, &totalParticles); CHKERRQ(ierr);
3166
3167 local_metrics[SEARCH_METRIC_SUM_SEARCH_ATTEMPTS] = (PetscReal)simCtx->searchMetrics.searchAttempts;
3168 local_metrics[SEARCH_METRIC_SUM_SEARCH_POPULATION] = (PetscReal)simCtx->searchMetrics.searchPopulation;
3169 local_metrics[SEARCH_METRIC_SUM_SEARCH_LOCATED] = (PetscReal)simCtx->searchMetrics.searchLocatedCount;
3170 local_metrics[SEARCH_METRIC_SUM_SEARCH_LOST] = (PetscReal)simCtx->searchMetrics.searchLostCount;
3171 local_metrics[SEARCH_METRIC_SUM_TRAVERSAL_STEPS] = (PetscReal)simCtx->searchMetrics.traversalStepsSum;
3172 local_metrics[SEARCH_METRIC_SUM_RESEARCH] = (PetscReal)simCtx->searchMetrics.reSearchCount;
3173 local_metrics[SEARCH_METRIC_SUM_TIE_BREAKS] = (PetscReal)simCtx->searchMetrics.tieBreakCount;
3174 local_metrics[SEARCH_METRIC_SUM_BOUNDARY_CLAMPS] = (PetscReal)simCtx->searchMetrics.boundaryClampCount;
3175 local_metrics[SEARCH_METRIC_SUM_BBOX_GUESS_SUCCESS] = (PetscReal)simCtx->searchMetrics.bboxGuessSuccessCount;
3176 local_metrics[SEARCH_METRIC_SUM_BBOX_GUESS_FALLBACK] = (PetscReal)simCtx->searchMetrics.bboxGuessFallbackCount;
3177 local_metrics[SEARCH_METRIC_SUM_MAX_TRAVERSAL_FAILS] = (PetscReal)simCtx->searchMetrics.maxTraversalFailCount;
3178 local_metrics[SEARCH_METRIC_MAX_TRAVERSAL_STEPS] = (PetscReal)simCtx->searchMetrics.maxTraversalSteps;
3179 local_metrics[SEARCH_METRIC_MAX_PASS_DEPTH] = (PetscReal)simCtx->searchMetrics.maxParticlePassDepth;
3180
3181 ierr = MPI_Op_create(SearchMetricsReduceOp, PETSC_TRUE, &reduction_op); CHKERRMPI(ierr);
3182 ierr = MPI_Allreduce(local_metrics, global_metrics, SEARCH_METRIC_REDUCTION_LEN, MPIU_REAL, reduction_op, PETSC_COMM_WORLD); CHKERRMPI(ierr);
3183 ierr = MPI_Op_free(&reduction_op); CHKERRMPI(ierr);
3184 reduction_op = MPI_OP_NULL;
3185
3186 searchAttempts = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_SEARCH_ATTEMPTS] + 0.5);
3187 searchPopulation = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_SEARCH_POPULATION] + 0.5);
3188 searchLocatedCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_SEARCH_LOCATED] + 0.5);
3189 searchLostCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_SEARCH_LOST] + 0.5);
3190 traversalStepsSum = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_TRAVERSAL_STEPS] + 0.5);
3191 reSearchCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_RESEARCH] + 0.5);
3192 tieBreakCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_TIE_BREAKS] + 0.5);
3193 boundaryClampCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_BOUNDARY_CLAMPS] + 0.5);
3194 bboxGuessSuccessCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_BBOX_GUESS_SUCCESS] + 0.5);
3195 bboxGuessFallbackCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_BBOX_GUESS_FALLBACK] + 0.5);
3196 maxTraversalFailCount = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_SUM_MAX_TRAVERSAL_FAILS] + 0.5);
3197 maxTraversalSteps = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_MAX_TRAVERSAL_STEPS] + 0.5);
3198 maxPassDepth = (long long)PetscFloorReal(global_metrics[SEARCH_METRIC_MAX_PASS_DEPTH] + 0.5);
3199
3200 if (searchAttempts > 0) {
3201 meanTraversalSteps = (PetscReal)traversalStepsSum / (PetscReal)searchAttempts;
3202 }
3203 if (searchPopulation > 0) {
3204 searchFailureFraction = (PetscReal)searchLostCount / (PetscReal)searchPopulation;
3205 searchWorkIndex = (PetscReal)traversalStepsSum / (PetscReal)searchPopulation;
3206 reSearchFraction = (PetscReal)reSearchCount / (PetscReal)searchPopulation;
3207 }
3208
3209 if (simCtx->rank == 0) {
3210 char csv_path[PETSC_MAX_PATH_LEN + 32];
3211 FILE *f = NULL;
3212
3213 ierr = PetscSNPrintf(csv_path, sizeof(csv_path), "%s/search_metrics.csv", simCtx->log_dir); CHKERRQ(ierr);
3214 f = fopen(csv_path, "a");
3215 if (!f) {
3216 LOG_ALLOW(GLOBAL, LOG_WARNING, "LOG_SEARCH_METRICS: could not open '%s' for writing.\n", csv_path);
3217 } else {
3218 if (ftell(f) == 0) {
3219 fprintf(f,
3220 "step,time,total_particles,lost,lost_cumulative,migrated,migration_passes,search_attempts,"
3221 "mean_traversal_steps,max_traversal_steps,tie_break_count,boundary_clamp_count,"
3222 "bbox_guess_success_count,bbox_guess_fallback_count,max_particle_pass_depth,load_imbalance,"
3223 "search_population,search_located_count,search_lost_count,traversal_steps_sum,re_search_count,"
3224 "max_traversal_fail_count,search_failure_fraction,search_work_index,re_search_fraction\n");
3225 }
3226 if (simCtx->continueMode && simCtx->step == simCtx->StartStep + 1) {
3227 fprintf(f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
3228 }
3229 fprintf(f,
3230 "%d,%.6e,%d,%d,%d,%d,%d,%lld,%.6e,%lld,%lld,%lld,%lld,%lld,%lld,%.6e,%lld,%lld,%lld,%lld,%lld,%lld,%.6e,%.6e,%.6e\n",
3231 (int)simCtx->step,
3232 (double)simCtx->ti,
3233 (int)totalParticles,
3234 (int)simCtx->particlesLostLastStep,
3235 (int)simCtx->particlesLostCumulative,
3236 (int)simCtx->particlesMigratedLastStep,
3237 (int)simCtx->migrationPassesLastStep,
3238 searchAttempts,
3239 (double)meanTraversalSteps,
3240 maxTraversalSteps,
3241 tieBreakCount,
3242 boundaryClampCount,
3243 bboxGuessSuccessCount,
3244 bboxGuessFallbackCount,
3245 maxPassDepth,
3246 (double)simCtx->particleLoadImbalance,
3247 searchPopulation,
3248 searchLocatedCount,
3249 searchLostCount,
3250 traversalStepsSum,
3251 reSearchCount,
3252 maxTraversalFailCount,
3253 (double)searchFailureFraction,
3254 (double)searchWorkIndex,
3255 (double)reSearchFraction);
3256 fclose(f);
3257 }
3258 }
3259
3260 LOG_ALLOW(GLOBAL, LOG_DEBUG,
3261 "Search metrics: sff=%.3e swi=%.3e re_search=%.3e lost(step/total)=%d/%d migrated=%d passes=%d traversal(mean/max)=%.2f/%lld tie_breaks=%lld max_pass_depth=%lld\n",
3262 (double)searchFailureFraction,
3263 (double)searchWorkIndex,
3264 (double)reSearchFraction,
3265 (int)simCtx->particlesLostLastStep,
3266 (int)simCtx->particlesLostCumulative,
3267 (int)simCtx->particlesMigratedLastStep,
3268 (int)simCtx->migrationPassesLastStep,
3269 (double)meanTraversalSteps,
3270 maxTraversalSteps,
3271 tieBreakCount,
3272 maxPassDepth);
3273
3274 PetscFunctionReturn(0);
3275}
3276
3277#undef __FUNCT__
3278#define __FUNCT__ "CalculateAdvancedParticleMetrics"
3279/**
3280 * @brief Internal helper implementation: `CalculateAdvancedParticleMetrics()`.
3281 * @details Local to this translation unit.
3282 */
3283PetscErrorCode CalculateAdvancedParticleMetrics(UserCtx *user)
3284{
3285 PetscErrorCode ierr;
3286 SimCtx *simCtx = user->simCtx;
3287 PetscMPIInt size, rank;
3288
3289 PetscFunctionBeginUser;
3290 ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size); CHKERRQ(ierr);
3291 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
3292
3293 // --- 1. Particle Load Imbalance ---
3294 PetscInt nLocal, nGlobal, nLocalMax;
3295 ierr = DMSwarmGetLocalSize(user->swarm, &nLocal); CHKERRQ(ierr);
3296 ierr = DMSwarmGetSize(user->swarm, &nGlobal); CHKERRQ(ierr);
3297 ierr = MPI_Allreduce(&nLocal, &nLocalMax, 1, MPIU_INT, MPI_MAX, PETSC_COMM_WORLD); CHKERRQ(ierr);
3298
3299 PetscReal avg_per_rank = (size > 0) ? ((PetscReal)nGlobal / size) : 0.0;
3300 // Handle division by zero if there are no particles
3301 simCtx->particleLoadImbalance = (avg_per_rank > 1e-9) ? (nLocalMax / avg_per_rank) : 1.0;
3302
3303
3304 // --- 2. Number of Occupied Cells ---
3305 // This part requires access to the user->ParticleCount vector.
3306 PetscInt local_occupied_cells = 0;
3307 PetscInt global_occupied_cells;
3308 const PetscScalar *count_array;
3309 PetscInt vec_local_size;
3310
3311 ierr = VecGetLocalSize(user->ParticleCount, &vec_local_size); CHKERRQ(ierr);
3312 ierr = VecGetArrayRead(user->ParticleCount, &count_array); CHKERRQ(ierr);
3313
3314 for (PetscInt i = 0; i < vec_local_size; ++i) {
3315 if (count_array[i] > 0.5) { // Use 0.5 to be safe with floating point
3316 local_occupied_cells++;
3317 }
3318 }
3319 ierr = VecRestoreArrayRead(user->ParticleCount, &count_array); CHKERRQ(ierr);
3320
3321 ierr = MPI_Allreduce(&local_occupied_cells, &global_occupied_cells, 1, MPIU_INT, MPI_SUM, PETSC_COMM_WORLD); CHKERRQ(ierr);
3322 simCtx->occupiedCellCount = global_occupied_cells;
3323
3324 LOG_ALLOW_SYNC(GLOBAL, LOG_INFO, "[Rank %d] Advanced Metrics: Imbalance=%.2f, OccupiedCells=%d\n", rank, simCtx->particleLoadImbalance, simCtx->occupiedCellCount);
3325
3326 PetscFunctionReturn(0);
3327}
3328
3329#undef __FUNCT__
3330#define __FUNCT__ "LOG_PARTICLE_METRICS"
3331/**
3332 * @brief Implementation of \ref LOG_PARTICLE_METRICS().
3333 * @details Full API contract (arguments, ownership, side effects) is documented with
3334 * the header declaration in `include/logging.h`.
3335 * @see LOG_PARTICLE_METRICS()
3336 */
3337PetscErrorCode LOG_PARTICLE_METRICS(UserCtx *user, const char *stageName)
3338{
3339 PetscErrorCode ierr;
3340 PetscMPIInt rank;
3341 SimCtx *simCtx = user->simCtx;
3342 const char *stage_label = (stageName && stageName[0] != '\0') ? stageName : "N/A";
3343
3344 PetscFunctionBeginUser;
3345 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
3346
3347 PetscInt totalParticles;
3348 ierr = DMSwarmGetSize(user->swarm, &totalParticles); CHKERRQ(ierr);
3349
3350 if (!rank) {
3351 FILE *f;
3352 char filen[PETSC_MAX_PATH_LEN + 64];
3353 ierr = PetscSNPrintf(filen, sizeof(filen), "%s/Particle_Metrics.log", simCtx->log_dir); CHKERRQ(ierr);
3354 f = fopen(filen, "a");
3355 if (!f) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open particle log file: %s", filen);
3356
3357 if (ftell(f) == 0) {
3358 PetscFPrintf(PETSC_COMM_SELF, f, "%-18s | %-10s | %-12s | %-10s | %-10s | %-10s | %-15s | %-10s | %-10s\n",
3359 "Stage", "Timestep", "Total Ptls", "Lost", "Lost Total", "Migrated", "Occupied Cells", "Imbalance", "Mig Passes");
3360 PetscFPrintf(PETSC_COMM_SELF, f, "-------------------------------------------------------------------------------------------------------------------------------------------\n");
3361 }
3362 if (simCtx->continueMode && simCtx->step == simCtx->StartStep + 1) {
3363 PetscFPrintf(PETSC_COMM_SELF, f, "# Continuation from step %" PetscInt_FMT "\n", simCtx->StartStep);
3364 }
3365
3366 PetscFPrintf(PETSC_COMM_SELF, f, "%-18s | %-10d | %-12d | %-10d | %-10d | %-10d | %-15d | %-10.2f | %-10d\n",
3367 stage_label, (int)simCtx->step, (int)totalParticles, (int)simCtx->particlesLostLastStep,
3368 (int)simCtx->particlesLostCumulative, (int)simCtx->particlesMigratedLastStep, (int)simCtx->occupiedCellCount,
3369 (double)simCtx->particleLoadImbalance, (int)simCtx->migrationPassesLastStep);
3370 fclose(f);
3371 }
3372 PetscFunctionReturn(0);
3373}
SetAnalyticalScalarFieldAtCellCenters
PetscErrorCode SetAnalyticalScalarFieldAtCellCenters(UserCtx *user, Vec targetVec)
Writes the configured verification scalar profile at physical cell centers into a scalar Vec.
Definition AnalyticalSolutions.c:675
SetAnalyticalSolutionForParticles
PetscErrorCode SetAnalyticalSolutionForParticles(Vec tempVec, SimCtx *simCtx)
Applies the analytical solution to particle velocity vector.
Definition AnalyticalSolutions.c:536
FieldDescriptor::layout
FieldLayout layout
Definition field_catalog.h:126
FieldView::descriptor
const FieldDescriptor * descriptor
Definition field_catalog.h:136
FieldGetView
PetscErrorCode FieldGetView(UserCtx *user, FieldId field_id, FieldView *view)
Resolve the existing DM and global/local vectors for one field.
Definition field_catalog.c:259
FieldView::local_vec
Vec local_vec
Definition field_catalog.h:139
FieldView::dm
DM dm
Definition field_catalog.h:137
FieldDescriptor::dof
PetscInt dof
Definition field_catalog.h:124
FieldView::global_vec
Vec global_vec
Definition field_catalog.h:138
FieldLayout
FieldLayout
Logical storage topology of a field.
Definition field_catalog.h:71
FIELD_LAYOUT_K_FACE
@ FIELD_LAYOUT_K_FACE
Definition field_catalog.h:76
FIELD_LAYOUT_I_FACE
@ FIELD_LAYOUT_I_FACE
Definition field_catalog.h:74
FIELD_LAYOUT_CELL_CENTERED
@ FIELD_LAYOUT_CELL_CENTERED
Definition field_catalog.h:73
FIELD_LAYOUT_COMPONENT_STAGGERED
@ FIELD_LAYOUT_COMPONENT_STAGGERED
Definition field_catalog.h:77
FIELD_LAYOUT_NODE_CENTERED
@ FIELD_LAYOUT_NODE_CENTERED
Definition field_catalog.h:72
FIELD_LAYOUT_J_FACE
@ FIELD_LAYOUT_J_FACE
Definition field_catalog.h:75
FieldDescriptor::canonical_name
const char * canonical_name
Definition field_catalog.h:121
FieldLayoutName
const char * FieldLayoutName(FieldLayout layout)
Return a stable printable label for a field layout.
Definition field_catalog.c:240
FieldId
FieldId
Compile-time identity for a catalogued Eulerian field.
Definition field_catalog.h:23
FIELD_ID_CELL_SCALAR_AT_CORNER
@ FIELD_ID_CELL_SCALAR_AT_CORNER
Definition field_catalog.h:63
FIELD_ID_CELL_VECTOR_AT_CORNER
@ FIELD_ID_CELL_VECTOR_AT_CORNER
Definition field_catalog.h:64
FieldView
Non-owning runtime objects resolved for one field and UserCtx.
Definition field_catalog.h:135
LOG_FIELD_MIN_MAX
PetscErrorCode LOG_FIELD_MIN_MAX(UserCtx *user, FieldId field_id)
Implementation of LOG_FIELD_MIN_MAX().
Definition logging.c:2349
SEARCH_METRIC_SUM_SEARCH_LOCATED
@ SEARCH_METRIC_SUM_SEARCH_LOCATED
Definition logging.c:36
SEARCH_METRIC_SUM_MAX_TRAVERSAL_FAILS
@ SEARCH_METRIC_SUM_MAX_TRAVERSAL_FAILS
Definition logging.c:44
SEARCH_METRIC_REDUCTION_LEN
@ SEARCH_METRIC_REDUCTION_LEN
Definition logging.c:47
SEARCH_METRIC_SUM_BBOX_GUESS_FALLBACK
@ SEARCH_METRIC_SUM_BBOX_GUESS_FALLBACK
Definition logging.c:43
SEARCH_METRIC_SUM_RESEARCH
@ SEARCH_METRIC_SUM_RESEARCH
Definition logging.c:39
SEARCH_METRIC_SUM_TRAVERSAL_STEPS
@ SEARCH_METRIC_SUM_TRAVERSAL_STEPS
Definition logging.c:38
SEARCH_METRIC_MAX_TRAVERSAL_STEPS
@ SEARCH_METRIC_MAX_TRAVERSAL_STEPS
Definition logging.c:45
SEARCH_METRIC_SUM_BOUNDARY_CLAMPS
@ SEARCH_METRIC_SUM_BOUNDARY_CLAMPS
Definition logging.c:41
SEARCH_METRIC_SUM_SEARCH_POPULATION
@ SEARCH_METRIC_SUM_SEARCH_POPULATION
Definition logging.c:35
SEARCH_METRIC_MAX_PASS_DEPTH
@ SEARCH_METRIC_MAX_PASS_DEPTH
Definition logging.c:46
SEARCH_METRIC_SUM_TIE_BREAKS
@ SEARCH_METRIC_SUM_TIE_BREAKS
Definition logging.c:40
SEARCH_METRIC_SUM_SEARCH_LOST
@ SEARCH_METRIC_SUM_SEARCH_LOST
Definition logging.c:37
SEARCH_METRIC_SUM_BBOX_GUESS_SUCCESS
@ SEARCH_METRIC_SUM_BBOX_GUESS_SUCCESS
Definition logging.c:42
SEARCH_METRIC_SUM_SEARCH_ATTEMPTS
@ SEARCH_METRIC_SUM_SEARCH_ATTEMPTS
Definition logging.c:34
EmitStatisticsConsoleSnapshot
PetscErrorCode EmitStatisticsConsoleSnapshot(UserCtx *user, const SimCtx *simCtx, PetscInt step)
Implementation of EmitStatisticsConsoleSnapshot().
Definition logging.c:2789
set_allowed_functions
void set_allowed_functions(const char **functionList, int count)
Implementation of set_allowed_functions().
Definition logging.c:155
ProfiledFunction::always_log
PetscBool always_log
Definition logging.c:1880
AppendStatisticalObservableSample
static PetscErrorCode AppendStatisticalObservableSample(SimCtx *simCtx, PetscInt samples_before, PetscReal mean_speed, PetscReal mean_ke)
Appends one timestep's scalar observables to the statistical history.
Definition logging.c:1370
LogFieldAnatomyView
static PetscErrorCode LogFieldAnatomyView(UserCtx *user, const char *field_name, const char *stage_name, DM dm, Vec vec_local, PetscInt dof, FieldLayout layout)
Shared architecture-aware anatomy logger for catalog and transient fields.
Definition logging.c:2479
LOG_PARTICLE_METRICS
PetscErrorCode LOG_PARTICLE_METRICS(UserCtx *user, const char *stageName)
Implementation of LOG_PARTICLE_METRICS().
Definition logging.c:3337
ShouldEmitPeriodicStatisticsConsoleSnapshot
PetscBool ShouldEmitPeriodicStatisticsConsoleSnapshot(const SimCtx *simCtx, PetscInt completed_step)
Implementation of ShouldEmitPeriodicStatisticsConsoleSnapshot().
Definition logging.c:2778
BCHandlerTypeToString
const char * BCHandlerTypeToString(BCHandlerType handler_type)
Internal helper implementation: BCHandlerTypeToString().
Definition logging.c:793
is_function_allowed
PetscBool is_function_allowed(const char *functionName)
Implementation of is_function_allowed().
Definition logging.c:186
g_profiler_count
static PetscInt g_profiler_count
Definition logging.c:1885
DualMonitorDestroy
PetscErrorCode DualMonitorDestroy(void **ctx)
Implementation of DualMonitorDestroy().
Definition logging.c:831
TMP_BUF_SIZE
#define TMP_BUF_SIZE
Definition logging.c:10
gAllowedFunctions
static char ** gAllowedFunctions
Global/static array of function names allowed to log.
Definition logging.c:26
current_log_level
static LogLevel current_log_level
Static variable to cache the current logging level.
Definition logging.c:19
LOG_INTERPOLATION_ERROR
PetscErrorCode LOG_INTERPOLATION_ERROR(UserCtx *user)
Implementation of LOG_INTERPOLATION_ERROR().
Definition logging.c:2865
SolutionConvergenceHistoryGet
static PetscReal SolutionConvergenceHistoryGet(const PetscReal *history, PetscInt capacity, PetscInt samples_available, PetscInt offset_from_latest)
Reads one sample from the statistical ring buffer by age.
Definition logging.c:1335
ShouldEmitPeriodicParticleConsoleSnapshot
PetscBool ShouldEmitPeriodicParticleConsoleSnapshot(const SimCtx *simCtx, PetscInt completed_step)
Implementation of ShouldEmitPeriodicParticleConsoleSnapshot().
Definition logging.c:545
BCFaceToString
const char * BCFaceToString(BCFace face)
Implementation of BCFaceToString().
Definition logging.c:671
FreeAllowedFunctions
PetscErrorCode FreeAllowedFunctions(char **funcs, PetscInt n)
Internal helper implementation: FreeAllowedFunctions().
Definition logging.c:652
IsParticleConsoleSnapshotEnabled
PetscBool IsParticleConsoleSnapshotEnabled(const SimCtx *simCtx)
Implementation of IsParticleConsoleSnapshotEnabled().
Definition logging.c:528
print_log_level
PetscErrorCode print_log_level(void)
Internal helper implementation: print_log_level().
Definition logging.c:119
ProfiledFunction::total_call_count
long long total_call_count
Definition logging.c:1877
SolutionConvergenceSafeRelative
static PetscReal SolutionConvergenceSafeRelative(PetscReal numerator, PetscReal denominator)
Forms a guarded relative metric for solution-convergence logging.
Definition logging.c:949
EmitParticleConsoleSnapshot
PetscErrorCode EmitParticleConsoleSnapshot(UserCtx *user, SimCtx *simCtx, PetscInt step)
Implementation of EmitParticleConsoleSnapshot().
Definition logging.c:559
g_profiler_capacity
static PetscInt g_profiler_capacity
Definition logging.c:1886
ProfilingFinalize
PetscErrorCode ProfilingFinalize(SimCtx *simCtx)
Implementation of ProfilingFinalize().
Definition logging.c:2196
BuildRowFormatString
static void BuildRowFormatString(PetscMPIInt wRank, PetscInt wPID, PetscInt wCell, PetscInt wPos, PetscInt wVel, PetscInt wWt, char *fmtStr, size_t bufSize)
Definition logging.c:369
trim
static void trim(char *s)
Remove leading and trailing whitespace from a mutable configuration string.
Definition logging.c:573
LoadAllowedFunctionsFromFile
PetscErrorCode LoadAllowedFunctionsFromFile(const char filename[], char ***funcsOut, PetscInt *nOut)
Implementation of LoadAllowedFunctionsFromFile().
Definition logging.c:598
SOLUTION_CONVERGENCE_REL_EPS
#define SOLUTION_CONVERGENCE_REL_EPS
Definition logging.c:918
gNumAllowed
static int gNumAllowed
Number of entries in the gAllowedFunctions array.
Definition logging.c:31
ProfiledFunction::total_time
double total_time
Definition logging.c:1875
BuildHeaderString
static void BuildHeaderString(char *headerStr, size_t bufSize, PetscMPIInt wRank, PetscInt wPID, PetscInt wCell, PetscInt wPos, PetscInt wVel, PetscInt wWt)
Definition logging.c:382
PrintProgressBar
void PrintProgressBar(PetscInt step, PetscInt startStep, PetscInt totalSteps, PetscReal currentTime)
Internal helper implementation: PrintProgressBar().
Definition logging.c:2302
SolutionConvergenceModeToString
static const char * SolutionConvergenceModeToString(SolutionConvergenceMode mode)
Maps the internal solution-convergence mode enum to its log label.
Definition logging.c:1583
_FindOrCreateEntry
static PetscErrorCode _FindOrCreateEntry(const char *func_name, PetscInt *idx)
Find a profiling record by name or allocate and register a new record.
Definition logging.c:1892
CellToStr
static void CellToStr(const PetscInt *cell, char *buf, size_t bufsize)
Definition logging.c:272
RuntimeMemoryLogSample
PetscErrorCode RuntimeMemoryLogSample(SimCtx *simCtx, PetscInt step, const char *event, const char *reason)
Implementation of RuntimeMemoryLogSample().
Definition logging.c:2086
get_log_level
LogLevel get_log_level()
Implementation of get_log_level().
Definition logging.c:87
ProfilingLogTimestepSummary
PetscErrorCode ProfilingLogTimestepSummary(SimCtx *simCtx, PetscInt step)
Implementation of ProfilingLogTimestepSummary().
Definition logging.c:2005
LOG_FACE_DISTANCES
PetscErrorCode LOG_FACE_DISTANCES(PetscReal *d)
Implementation of LOG_FACE_DISTANCES().
Definition logging.c:233
LOG_PARTICLE_FIELDS
PetscErrorCode LOG_PARTICLE_FIELDS(UserCtx *user, PetscInt printInterval)
Implementation of LOG_PARTICLE_FIELDS().
Definition logging.c:400
_ProfilingEnd
void _ProfilingEnd(const char *func_name)
Implementation of _ProfilingEnd().
Definition logging.c:1966
TripleRealToStr
static void TripleRealToStr(const PetscReal *arr, char *buf, size_t bufsize)
Definition logging.c:280
Int64ToStr
static void Int64ToStr(PetscInt64 value, char *buf, size_t bufsize)
Definition logging.c:264
BCTypeToString
const char * BCTypeToString(BCType type)
Implementation of BCTypeToString().
Definition logging.c:773
CalculateAdvancedParticleMetrics
PetscErrorCode CalculateAdvancedParticleMetrics(UserCtx *user)
Internal helper implementation: CalculateAdvancedParticleMetrics().
Definition logging.c:3283
ParticleLocationStatusToString
const char * ParticleLocationStatusToString(ParticleLocationStatus level)
Implementation of ParticleLocationStatusToString().
Definition logging.c:1858
LOG_SCATTER_METRICS
PetscErrorCode LOG_SCATTER_METRICS(UserCtx *user)
Implementation of LOG_SCATTER_METRICS().
Definition logging.c:2944
_CompareProfiledFunctions
static int _CompareProfiledFunctions(const void *a, const void *b)
Order profiling records by their accumulated execution time.
Definition logging.c:2180
LOG_SOLUTION_CONVERGENCE
PetscErrorCode LOG_SOLUTION_CONVERGENCE(SimCtx *simCtx)
Implementation of LOG_SOLUTION_CONVERGENCE().
Definition logging.c:1600
FlowDirectionToString
const char * FlowDirectionToString(FlowDirection fd)
Convert a FlowDirection enum value to its YAML token string.
Definition logging.c:705
DualKSPMonitor
PetscErrorCode DualKSPMonitor(KSP ksp, PetscInt it, PetscReal rnorm, void *ctx)
Implementation of DualKSPMonitor().
Definition logging.c:870
LOG_CONTINUITY_METRICS
PetscErrorCode LOG_CONTINUITY_METRICS(UserCtx *user)
Implementation of LOG_CONTINUITY_METRICS().
Definition logging.c:1796
LOG_CORNER_FIELD_ANATOMY
PetscErrorCode LOG_CORNER_FIELD_ANATOMY(UserCtx *user, FieldId corner_field_id, const char *stage_name)
Emits anatomy for the transient scalar or vector corner-staging field.
Definition logging.c:2838
ProfiledFunction::current_step_time
double current_step_time
Definition logging.c:1876
LOG_FIELD_ANATOMY
PetscErrorCode LOG_FIELD_ANATOMY(UserCtx *user, FieldId field_id, const char *stage_name)
Resolves a persistent field view and emits its layout-aware boundary anatomy.
Definition logging.c:2752
ProfiledFunction::current_step_call_count
long long current_step_call_count
Definition logging.c:1878
ComputeStatisticalWindowMetrics
static PetscErrorCode ComputeStatisticalWindowMetrics(const SimCtx *simCtx, PetscInt samples_available, PetscBool *has_reference_out, PetscReal *mean_speed_window_out, PetscReal *mean_speed_window_prev_out, PetscReal *mean_speed_window_abs_out, PetscReal *mean_speed_window_rel_out, PetscReal *mean_speed_rms_window_out, PetscReal *mean_speed_rms_window_prev_out, PetscReal *mean_speed_rms_window_abs_out, PetscReal *mean_speed_rms_window_rel_out, PetscReal *mean_ke_window_out, PetscReal *mean_ke_window_prev_out, PetscReal *mean_ke_window_abs_out, PetscReal *mean_ke_window_rel_out, PetscReal *mean_ke_rms_window_out, PetscReal *mean_ke_rms_window_prev_out, PetscReal *mean_ke_rms_window_abs_out, PetscReal *mean_ke_rms_window_rel_out)
Computes adjacent-window drift metrics for statistical steady mode.
Definition logging.c:1450
LOG_SEARCH_METRICS
PetscErrorCode LOG_SEARCH_METRICS(UserCtx *user)
Implementation of LOG_SEARCH_METRICS().
Definition logging.c:3129
g_profiler_registry
static ProfiledFunction * g_profiler_registry
Definition logging.c:1884
InitialConditionModeToString
const char * InitialConditionModeToString(InitialConditionMode mode)
Implementation of InitialConditionModeToString().
Definition logging.c:689
ProfilingInitialize
PetscErrorCode ProfilingInitialize(SimCtx *simCtx)
Internal helper implementation: ProfilingInitialize().
Definition logging.c:1928
LESModelToString
const char * LESModelToString(LESModelType LESFlag)
Implementation of LESModelToString().
Definition logging.c:741
LOG_CELL_VERTICES
PetscErrorCode LOG_CELL_VERTICES(const Cell *cell, PetscMPIInt rank)
Implementation of LOG_CELL_VERTICES().
Definition logging.c:208
IntToStr
static void IntToStr(int value, char *buf, size_t bufsize)
Definition logging.c:256
ComputeDeterministicSolutionMetrics
static PetscErrorCode ComputeDeterministicSolutionMetrics(SimCtx *simCtx, PetscBool periodic_mode, PetscInt phase_step, PetscInt samples_before, PetscBool *has_reference_out, PetscReal *u_abs_l2_out, PetscReal *u_rel_l2_out, PetscReal *p_abs_l2_out, PetscReal *p_rel_l2_out, PetscReal *mean_speed_out, PetscReal *mean_speed_ref_out, PetscReal *mean_speed_abs_out, PetscReal *mean_speed_rel_out, PetscReal *mean_ke_out, PetscReal *mean_ke_ref_out, PetscReal *mean_ke_abs_out, PetscReal *mean_ke_rel_out)
Computes deterministic solution-drift metrics for the current step.
Definition logging.c:1092
ProfilingResetTimestepCounters
PetscErrorCode ProfilingResetTimestepCounters(void)
Implementation of ProfilingResetTimestepCounters().
Definition logging.c:1988
IsStatisticsConsoleSnapshotEnabled
PetscBool IsStatisticsConsoleSnapshotEnabled(const SimCtx *simCtx)
Implementation of IsStatisticsConsoleSnapshotEnabled().
Definition logging.c:2768
ResetSearchMetrics
PetscErrorCode ResetSearchMetrics(SimCtx *simCtx)
Implementation of ResetSearchMetrics().
Definition logging.c:3098
MomentumSolverTypeToString
const char * MomentumSolverTypeToString(MomentumSolverType SolverFlag)
Implementation of MomentumSolverTypeToString().
Definition logging.c:757
ComputeMaxColumnWidths
static PetscErrorCode ComputeMaxColumnWidths(PetscInt nParticles, const PetscMPIInt *ranks, const PetscInt64 *pids, const PetscInt *cellIDs, const PetscReal *positions, const PetscReal *velocities, const PetscReal *weights, int *wRank, int *wPID, int *wCell, int *wPos, int *wVel, int *wWt)
Definition logging.c:305
ComputeCurrentFlowObservables
static PetscErrorCode ComputeCurrentFlowObservables(SimCtx *simCtx, PetscReal *mean_speed_out, PetscReal *mean_ke_out)
Computes instantaneous global flow observables for statistical mode.
Definition logging.c:976
SOLUTION_CONVERGENCE_FLUID_THRESHOLD
#define SOLUTION_CONVERGENCE_FLUID_THRESHOLD
Definition logging.c:917
ProfiledFunction::start_time
double start_time
Definition logging.c:1879
ParticleInitializationToString
const char * ParticleInitializationToString(ParticleInitializationType ParticleInitialization)
Implementation of ParticleInitializationToString().
Definition logging.c:724
_ProfilingStart
void _ProfilingStart(const char *func_name)
Implementation of _ProfilingStart().
Definition logging.c:1952
SearchMetricsReduceOp
static void SearchMetricsReduceOp(void *invec, void *inoutvec, int *len, MPI_Datatype *datatype)
Internal reduction callback for packed search metrics.
Definition logging.c:54
ProfiledFunction::name
const char * name
Definition logging.c:1874
logging.h
Logging utilities and macros for PETSc-based applications.
LOG_ALLOW_SYNC
#define LOG_ALLOW_SYNC(scope, level, fmt,...)
Synchronized logging macro that checks both the log level and whether the calling function is in the ...
Definition logging.h:253
DualMonitorCtx::log_to_console
PetscBool log_to_console
Definition logging.h:58
LOCAL
#define LOCAL
Logging scope definitions for controlling message output.
Definition logging.h:45
GLOBAL
#define GLOBAL
Scope for global logging across all processes.
Definition logging.h:46
LOG_ALLOW
#define LOG_ALLOW(scope, level, fmt,...)
Logging macro that checks both the log level and whether the calling function is in the allowed-funct...
Definition logging.h:200
DualMonitorCtx::bnorm
PetscReal bnorm
Definition logging.h:59
DualMonitorCtx::step
PetscInt step
Definition logging.h:60
LOG
#define LOG(scope, level, fmt,...)
Logging macro for PETSc-based applications with scope control.
Definition logging.h:84
LogLevel
LogLevel
Enumeration of logging levels.
Definition logging.h:28
LOG_ERROR
@ LOG_ERROR
Critical errors that may halt the program.
Definition logging.h:29
LOG_TRACE
@ LOG_TRACE
Very fine-grained tracing information for in-depth debugging.
Definition logging.h:33
LOG_INFO
@ LOG_INFO
Informational messages about program execution.
Definition logging.h:31
LOG_WARNING
@ LOG_WARNING
Non-critical issues that warrant attention.
Definition logging.h:30
LOG_DEBUG
@ LOG_DEBUG
Detailed debugging information.
Definition logging.h:32
LOG_VERBOSE
@ LOG_VERBOSE
Extremely detailed logs, typically for development use only.
Definition logging.h:34
DualMonitorCtx::file_handle
FILE * file_handle
Definition logging.h:57
DualMonitorCtx::block_id
PetscInt block_id
Definition logging.h:61
DualMonitorCtx
Context for a dual-purpose KSP monitor.
Definition logging.h:56
particle_field_catalog.h
Typed identities and metadata for persistent solver-particle fields.
ParticleFieldName
const char * ParticleFieldName(ParticleFieldId field_id)
Return the canonical PETSc DMSwarm name for an ID.
Definition particle_field_catalog.c:111
PARTICLE_FIELD_ID_WEIGHT
@ PARTICLE_FIELD_ID_WEIGHT
Definition particle_field_catalog.h:28
PARTICLE_FIELD_ID_POSITION
@ PARTICLE_FIELD_ID_POSITION
Definition particle_field_catalog.h:25
PARTICLE_FIELD_ID_PID
@ PARTICLE_FIELD_ID_PID
Definition particle_field_catalog.h:33
PARTICLE_FIELD_ID_CELL_ID
@ PARTICLE_FIELD_ID_CELL_ID
Definition particle_field_catalog.h:27
PARTICLE_FIELD_ID_PSI
@ PARTICLE_FIELD_ID_PSI
Definition particle_field_catalog.h:31
PARTICLE_FIELD_ID_RANK
@ PARTICLE_FIELD_ID_RANK
Definition particle_field_catalog.h:34
PARTICLE_FIELD_ID_VELOCITY
@ PARTICLE_FIELD_ID_VELOCITY
Definition particle_field_catalog.h:26
statistics_accumulator.h
Per-window PETSc accumulator storage and pointwise application.
PicurvWindowValidFractionRange
PetscErrorCode PicurvWindowValidFractionRange(UserCtx *user, const PicurvWindowDefinition *definition, const PicurvWindowStorage *storage, PetscInt sample_count, PetscReal *minimum, PetscReal *maximum)
Reports the range of per-point valid fraction across a window's domain.
Definition statistics_accumulator.c:780
statistics_window.h
Window lifecycle, scheduling, and weighting for the field-statistics pipeline.
PicurvWindow::sample_count
PetscInt sample_count
Definition statistics_window.h:81
PicurvWindow::state
PicurvWindowState state
Definition statistics_window.h:77
PicurvWindowDefinition::name
char name[64]
Definition statistics_window.h:60
PicurvWindow::total_weight
PetscReal total_weight
Definition statistics_window.h:82
PicurvWindowStateName
const char * PicurvWindowStateName(PicurvWindowState state)
Returns a stable human-readable name for a window state.
Definition statistics_window.c:23
PicurvWindowDefinition::bounded
PetscBool bounded
False for an open-ended window.
Definition statistics_window.h:63
PicurvWindow::definition
PicurvWindowDefinition definition
Definition statistics_window.h:76
FieldStatisticsIsActive
PetscBool FieldStatisticsIsActive(const struct SimCtx *simCtx)
Reports whether this run has live field-statistics state.
PicurvWindowProgress
PetscReal PicurvWindowProgress(const PicurvWindow *window)
Reports the fraction of a bounded window's span that has been represented.
Definition statistics_window.c:414
PicurvWindow::represented_time
PetscReal represented_time
Physical time the window covers.
Definition statistics_window.h:83
PicurvWindow
Runtime state of one window.
Definition statistics_window.h:75
LESModelType
LESModelType
Identifies the six logical faces of a structured computational block.
Definition variables.h:520
DYNAMIC_SMAGORINSKY
@ DYNAMIC_SMAGORINSKY
Definition variables.h:523
NO_LES_MODEL
@ NO_LES_MODEL
Definition variables.h:521
CONSTANT_SMAGORINSKY
@ CONSTANT_SMAGORINSKY
Definition variables.h:522
SimCtx::fieldStatisticsWindowCount
PetscInt fieldStatisticsWindowCount
Definition variables.h:770
BCType
BCType
Defines the general mathematical/physical Category of a boundary.
Definition variables.h:283
INLET
@ INLET
Definition variables.h:290
INTERFACE
@ INTERFACE
Definition variables.h:285
FARFIELD
@ FARFIELD
Definition variables.h:291
OUTLET
@ OUTLET
Definition variables.h:289
PERIODIC
@ PERIODIC
Definition variables.h:292
WALL
@ WALL
Definition variables.h:286
SimCtx::continueMode
PetscBool continueMode
Definition variables.h:712
MGCtx::user
UserCtx * user
Definition variables.h:571
SimCtx::profilingFinalSummary
PetscBool profilingFinalSummary
Definition variables.h:868
SimCtx::rank
PetscMPIInt rank
Definition variables.h:698
SimCtx::profilingTimestepFile
char profilingTimestepFile[PETSC_MAX_PATH_LEN]
Definition variables.h:867
SearchMetricsState::searchLocatedCount
PetscInt64 searchLocatedCount
Definition variables.h:241
SimCtx::statisticsConsoleOutputFreq
PetscInt statisticsConsoleOutputFreq
Definition variables.h:772
SimCtx::block_number
PetscInt block_number
Definition variables.h:790
SearchMetricsState::searchLostCount
PetscInt64 searchLostCount
Definition variables.h:242
UserCtx::solutionConvergencePeriodicPRef
Vec * solutionConvergencePeriodicPRef
Definition variables.h:949
UserCtx::simCtx
SimCtx * simCtx
Back-pointer to the master simulation context.
Definition variables.h:909
ParticleInitializationType
ParticleInitializationType
Enumerator to identify the particle initialization strategy.
Definition variables.h:551
PARTICLE_INIT_SURFACE_RANDOM
@ PARTICLE_INIT_SURFACE_RANDOM
Random placement on the inlet face.
Definition variables.h:552
PARTICLE_INIT_SURFACE_EDGES
@ PARTICLE_INIT_SURFACE_EDGES
Deterministic placement at inlet face edges.
Definition variables.h:555
PARTICLE_INIT_POINT_SOURCE
@ PARTICLE_INIT_POINT_SOURCE
All particles at a fixed (psrc_x,psrc_y,psrc_z) — for validation.
Definition variables.h:554
PARTICLE_INIT_VOLUME
@ PARTICLE_INIT_VOLUME
Random volumetric distribution across the domain.
Definition variables.h:553
ParticleLocationStatus
ParticleLocationStatus
Defines the state of a particle with respect to its location and migration status during the iterativ...
Definition variables.h:137
LOST
@ LOST
Definition variables.h:141
NEEDS_LOCATION
@ NEEDS_LOCATION
Definition variables.h:138
ACTIVE_AND_LOCATED
@ ACTIVE_AND_LOCATED
Definition variables.h:139
UNINITIALIZED
@ UNINITIALIZED
Definition variables.h:142
MIGRATING_OUT
@ MIGRATING_OUT
Definition variables.h:140
SimCtx::solutionConvergenceMeanSpeedHistory
PetscReal * solutionConvergenceMeanSpeedHistory
Definition variables.h:766
SimCtx::FluxOutSum
PetscReal FluxOutSum
Definition variables.h:799
SimCtx::runtimeMemoryLogEnabled
PetscBool runtimeMemoryLogEnabled
Enable the rank-reduced runtime memory log.
Definition variables.h:883
SearchMetricsState::boundaryClampCount
PetscInt64 boundaryClampCount
Definition variables.h:248
SimCtx::particlesLostLastStep
PetscInt particlesLostLastStep
Definition variables.h:834
UserCtx::KM
PetscInt KM
Definition variables.h:920
SimCtx::usermg
UserMG usermg
Definition variables.h:852
UserCtx::solutionConvergencePeriodicUcatRef
Vec * solutionConvergencePeriodicUcatRef
Definition variables.h:948
SearchMetricsState::traversalStepsSum
PetscInt64 traversalStepsSum
Definition variables.h:243
BCHandlerType
BCHandlerType
Defines the specific computational "strategy" for a boundary handler.
Definition variables.h:303
BC_HANDLER_INLET_PULSATILE_FLUX
@ BC_HANDLER_INLET_PULSATILE_FLUX
Definition variables.h:312
BC_HANDLER_PERIODIC_GEOMETRIC
@ BC_HANDLER_PERIODIC_GEOMETRIC
Definition variables.h:316
BC_HANDLER_INLET_PARABOLIC
@ BC_HANDLER_INLET_PARABOLIC
Definition variables.h:309
BC_HANDLER_INLET_CONSTANT_VELOCITY
@ BC_HANDLER_INLET_CONSTANT_VELOCITY
Definition variables.h:308
BC_HANDLER_PERIODIC_DRIVEN_INITIAL_FLUX
@ BC_HANDLER_PERIODIC_DRIVEN_INITIAL_FLUX
Definition variables.h:319
BC_HANDLER_INTERFACE_OVERSET
@ BC_HANDLER_INTERFACE_OVERSET
Definition variables.h:317
BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX
@ BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX
Definition variables.h:318
BC_HANDLER_WALL_MOVING
@ BC_HANDLER_WALL_MOVING
Definition variables.h:306
BC_HANDLER_INLET_PROFILE_FROM_FILE
@ BC_HANDLER_INLET_PROFILE_FROM_FILE
Definition variables.h:310
BC_HANDLER_WALL_NOSLIP
@ BC_HANDLER_WALL_NOSLIP
Definition variables.h:305
BC_HANDLER_OUTLET_CONSERVATION
@ BC_HANDLER_OUTLET_CONSERVATION
Definition variables.h:314
BC_HANDLER_FARFIELD_NONREFLECTING
@ BC_HANDLER_FARFIELD_NONREFLECTING
Definition variables.h:313
BC_HANDLER_OUTLET_PRESSURE
@ BC_HANDLER_OUTLET_PRESSURE
Definition variables.h:315
BC_HANDLER_SYMMETRY_PLANE
@ BC_HANDLER_SYMMETRY_PLANE
Definition variables.h:307
BC_HANDLER_UNDEFINED
@ BC_HANDLER_UNDEFINED
Definition variables.h:304
SimCtx::solutionConvergenceSamplesRecorded
PetscInt solutionConvergenceSamplesRecorded
Definition variables.h:765
UserCtx::_this
PetscInt _this
Definition variables.h:924
SearchMetricsState::searchPopulation
PetscInt64 searchPopulation
Definition variables.h:240
SimCtx::solutionConvergenceEnabled
PetscBool solutionConvergenceEnabled
Definition variables.h:761
SimCtx::solutionConvergenceMeanKEHistory
PetscReal * solutionConvergenceMeanKEHistory
Definition variables.h:767
SimCtx::dt
PetscReal dt
Definition variables.h:710
SimCtx::runtimeMemoryLogFile
char runtimeMemoryLogFile[PETSC_MAX_PATH_LEN]
File name written under log_dir.
Definition variables.h:884
SimCtx::runtimeMemoryLogStarted
PetscBool runtimeMemoryLogStarted
True after rank 0 writes the log header.
Definition variables.h:885
SimCtx::occupiedCellCount
PetscInt occupiedCellCount
Definition variables.h:838
SimCtx::profilingTimestepMode
char profilingTimestepMode[32]
Definition variables.h:866
SimCtx::k_periodic
PetscInt k_periodic
Definition variables.h:791
SearchMetricsState::currentSettlementPass
PetscInt currentSettlementPass
Definition variables.h:252
SimCtx::np
PetscInt np
Definition variables.h:827
SimCtx::StartStep
PetscInt StartStep
Definition variables.h:705
MomentumSolverType
MomentumSolverType
Enumerator to identify the implemented momentum solver strategies.
Definition variables.h:534
MOMENTUM_SOLVER_DUALTIME_PICARD_JAMESON_RK
@ MOMENTUM_SOLVER_DUALTIME_PICARD_JAMESON_RK
Definition variables.h:536
MOMENTUM_SOLVER_EXPLICIT_RK
@ MOMENTUM_SOLVER_EXPLICIT_RK
Definition variables.h:535
MOMENTUM_SOLVER_NEWTON_KRYLOV
@ MOMENTUM_SOLVER_NEWTON_KRYLOV
Definition variables.h:537
SimCtx::solutionConvergencePeriodSteps
PetscInt solutionConvergencePeriodSteps
Definition variables.h:763
Cmpnts::x
PetscScalar x
Definition variables.h:103
SearchMetricsState::reSearchCount
PetscInt64 reSearchCount
Definition variables.h:244
SimCtx::MaxDiv
PetscReal MaxDiv
Definition variables.h:859
SearchMetricsState::bboxGuessFallbackCount
PetscInt64 bboxGuessFallbackCount
Definition variables.h:250
UserCtx::Ucat_o
Vec Ucat_o
Definition variables.h:946
SimCtx::MaxDivx
PetscInt MaxDivx
Definition variables.h:860
SimCtx::MaxDivy
PetscInt MaxDivy
Definition variables.h:860
SearchMetricsState::bboxGuessSuccessCount
PetscInt64 bboxGuessSuccessCount
Definition variables.h:249
SimCtx::MaxDivz
PetscInt MaxDivz
Definition variables.h:860
SimCtx::fieldStatisticsWindows
struct PicurvWindow * fieldStatisticsWindows
Definition variables.h:771
SimCtx::log_dir
char log_dir[PETSC_MAX_PATH_LEN]
Definition variables.h:718
SimCtx::MaxDivFlatArg
PetscInt MaxDivFlatArg
Definition variables.h:860
SimCtx::FluxInSum
PetscReal FluxInSum
Definition variables.h:799
SearchMetricsState::maxParticlePassDepth
PetscInt64 maxParticlePassDepth
Definition variables.h:251
SearchMetricsState::maxTraversalSteps
PetscInt64 maxTraversalSteps
Definition variables.h:245
Cmpnts::z
PetscScalar z
Definition variables.h:103
UserCtx::ParticleCount
Vec ParticleCount
Definition variables.h:996
UserCtx::JM
PetscInt JM
Definition variables.h:920
SimCtx::runtimeMemoryLogHasPrevious
PetscBool runtimeMemoryLogHasPrevious
True after the first process-memory sample.
Definition variables.h:886
UserMG::mglevels
PetscInt mglevels
Definition variables.h:578
SimCtx::profilingSelectedFuncs
char ** profilingSelectedFuncs
Definition variables.h:864
SimCtx::solutionConvergenceWindowSteps
PetscInt solutionConvergenceWindowSteps
Definition variables.h:764
FlowDirection
FlowDirection
Primary flow direction for streamwise IC and Poiseuille modes.
Definition variables.h:272
FLOW_DIR_NEG_ZETA
@ FLOW_DIR_NEG_ZETA
Definition variables.h:278
FLOW_DIR_NEG_ETA
@ FLOW_DIR_NEG_ETA
Definition variables.h:276
FLOW_DIR_POS_ZETA
@ FLOW_DIR_POS_ZETA
Definition variables.h:277
FLOW_DIR_POS_XI
@ FLOW_DIR_POS_XI
Definition variables.h:273
FLOW_DIR_NEG_XI
@ FLOW_DIR_NEG_XI
Definition variables.h:274
FLOW_DIR_POS_ETA
@ FLOW_DIR_POS_ETA
Definition variables.h:275
SimCtx::particlesLostCumulative
PetscInt particlesLostCumulative
Definition variables.h:835
SimCtx::nProfilingSelectedFuncs
PetscInt nProfilingSelectedFuncs
Definition variables.h:865
SimCtx::particlesMigratedLastStep
PetscInt particlesMigratedLastStep
Definition variables.h:837
UserCtx::fieldStatisticsStorage
struct PicurvWindowStorage * fieldStatisticsStorage
Definition variables.h:962
InitialConditionMode
InitialConditionMode
Selects the algorithm used to populate a fresh Eulerian velocity field.
Definition variables.h:151
IC_MODE_CONSTANT_CARTESIAN
@ IC_MODE_CONSTANT_CARTESIAN
Definition variables.h:153
IC_MODE_POISEUILLE
@ IC_MODE_POISEUILLE
Definition variables.h:154
IC_MODE_CONSTANT_STREAMWISE
@ IC_MODE_CONSTANT_STREAMWISE
Definition variables.h:155
IC_MODE_FILE
@ IC_MODE_FILE
Definition variables.h:156
IC_MODE_ZERO
@ IC_MODE_ZERO
Definition variables.h:152
SimCtx::particleConsoleOutputFreq
PetscInt particleConsoleOutputFreq
Definition variables.h:708
SimCtx::searchMetrics
SearchMetricsState searchMetrics
Definition variables.h:840
SimCtx::i_periodic
PetscInt i_periodic
Definition variables.h:791
SimCtx::runtimeMemoryLogPreviousProcessMB
PetscReal runtimeMemoryLogPreviousProcessMB
Previous local process memory sample in MB.
Definition variables.h:887
SimCtx::step
PetscInt step
Definition variables.h:703
UserCtx::info
DMDALocalInfo info
Definition variables.h:918
SimCtx::migrationPassesLastStep
PetscInt migrationPassesLastStep
Definition variables.h:836
Cmpnts::y
PetscScalar y
Definition variables.h:103
EXEC_MODE_SOLVER
@ EXEC_MODE_SOLVER
Definition variables.h:668
EXEC_MODE_POSTPROCESSOR
@ EXEC_MODE_POSTPROCESSOR
Definition variables.h:669
UserCtx::IM
PetscInt IM
Definition variables.h:920
TOP
@ TOP
Definition variables.h:147
FRONT
@ FRONT
Definition variables.h:147
BOTTOM
@ BOTTOM
Definition variables.h:147
BACK
@ BACK
Definition variables.h:147
LEFT
@ LEFT
Definition variables.h:147
RIGHT
@ RIGHT
Definition variables.h:147
UserMG::mgctx
MGCtx * mgctx
Definition variables.h:581
SolutionConvergenceMode
SolutionConvergenceMode
Selects the runtime solution-convergence diagnostics mode.
Definition variables.h:543
SOLUTION_CONVERGENCE_TRANSIENT
@ SOLUTION_CONVERGENCE_TRANSIENT
Definition variables.h:547
SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC
@ SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC
Definition variables.h:545
SOLUTION_CONVERGENCE_STATISTICAL_STEADY
@ SOLUTION_CONVERGENCE_STATISTICAL_STEADY
Definition variables.h:546
SOLUTION_CONVERGENCE_STEADY_DETERMINISTIC
@ SOLUTION_CONVERGENCE_STEADY_DETERMINISTIC
Definition variables.h:544
SimCtx::solutionConvergenceMode
SolutionConvergenceMode solutionConvergenceMode
Definition variables.h:762
SearchMetricsState::searchAttempts
PetscInt64 searchAttempts
Definition variables.h:239
SimCtx::exec_mode
ExecutionMode exec_mode
Definition variables.h:714
SearchMetricsState::tieBreakCount
PetscInt64 tieBreakCount
Definition variables.h:247
SimCtx::ti
PetscReal ti
Definition variables.h:704
SimCtx::summationRHS
PetscReal summationRHS
Definition variables.h:858
SearchMetricsState::maxTraversalFailCount
PetscInt64 maxTraversalFailCount
Definition variables.h:246
SimCtx::LoggingFrequency
PetscInt LoggingFrequency
Definition variables.h:857
Cell::vertices
Cmpnts vertices[8]
Coordinates of the eight vertices of the cell.
Definition variables.h:178
UserCtx::Psi
Vec Psi
Definition variables.h:997
SimCtx::particleLoadImbalance
PetscReal particleLoadImbalance
Definition variables.h:839
UserCtx::P_o
Vec P_o
Definition variables.h:946
BCFace
BCFace
Identifies the six logical faces of a structured computational block.
Definition variables.h:261
BC_FACE_NEG_X
@ BC_FACE_NEG_X
Definition variables.h:262
BC_FACE_POS_Z
@ BC_FACE_POS_Z
Definition variables.h:264
BC_FACE_POS_Y
@ BC_FACE_POS_Y
Definition variables.h:263
BC_FACE_NEG_Z
@ BC_FACE_NEG_Z
Definition variables.h:264
BC_FACE_POS_X
@ BC_FACE_POS_X
Definition variables.h:262
BC_FACE_NEG_Y
@ BC_FACE_NEG_Y
Definition variables.h:263
SimCtx::j_periodic
PetscInt j_periodic
Definition variables.h:791
Cell
Defines the vertices of a single hexahedral grid cell.
Definition variables.h:177
Cmpnts
A 3D point or vector with PetscScalar components.
Definition variables.h:102
SimCtx
The master context for the entire simulation.
Definition variables.h:695
UserCtx
User-defined context containing data specific to a single computational grid level.
Definition variables.h:906
verification_sources.h
VerificationScalarOverrideActive
PetscBool VerificationScalarOverrideActive(const SimCtx *simCtx)
Reports whether a verification-only scalar override is active.
Definition verification_sources.c:21
Generated by doxygen 1.9.8