PICurv 0.1.0
A Parallel Particle-In-Cell Solver for Curvilinear LES
Loading...
Searching...
No Matches
Macros | Functions | Variables
io.c File Reference

Implementation of data input/output routines, focusing on grid configuration. More...

#include "io.h"
Include dependency graph for io.c:

Go to the source code of this file.

Macros

#define __FUNCT__   "ReadGridGenerationInputs"
 
#define __FUNCT__   "ReadGridFile"
 
#define __FUNCT__   "ParseAllBoundaryConditions"
 
#define __FUNCT__   "DeterminePeriodicity"
 
#define __FUNCT__   "ReadFieldData"
 
#define __FUNCT__   "WriteFieldData"
 
#define __FUNCT__   "ParsePostProcessingSettings"
 
#define __FUNCT__   "ParseScalingInformation"
 

Functions

void TrimWhitespace (char *str)
 Implementation of TrimWhitespace().
 
PetscBool ShouldWriteDataOutput (const SimCtx *simCtx, PetscInt completed_step)
 Implementation of ShouldWriteDataOutput().
 
PetscErrorCode ReadGridGenerationInputs (UserCtx *user)
 Internal helper implementation: ReadGridGenerationInputs().
 
PetscErrorCode PopulateFinestUserGridResolutionFromOptions (UserCtx *finest_users, PetscInt nblk)
 Internal helper implementation: PopulateFinestUserGridResolutionFromOptions().
 
PetscErrorCode ReadGridFile (UserCtx *user)
 Internal helper implementation: ReadGridFile().
 
void FreeBC_ParamList (BC_Param *head)
 Implementation of FreeBC_ParamList().
 
PetscErrorCode StringToBCFace (const char *str, BCFace *face_out)
 Internal helper implementation: StringToBCFace().
 
PetscErrorCode StringToBCType (const char *str, BCType *type_out)
 Internal helper implementation: StringToBCType().
 
PetscErrorCode StringToBCHandlerType (const char *str, BCHandlerType *handler_out)
 Internal helper implementation: StringToBCHandlerType().
 
PetscErrorCode ValidateBCHandlerForBCType (BCType type, BCHandlerType handler)
 Internal helper implementation: ValidateBCHandlerForBCType().
 
PetscErrorCode GetBCParamReal (BC_Param *params, const char *key, PetscReal *value_out, PetscBool *found)
 Internal helper implementation: GetBCParamReal().
 
PetscErrorCode GetBCParamBool (BC_Param *params, const char *key, PetscBool *value_out, PetscBool *found)
 Internal helper implementation: GetBCParamBool().
 
PetscErrorCode ParseAllBoundaryConditions (UserCtx *user, const char *bcs_input_filename)
 Internal helper implementation: ParseAllBoundaryConditions().
 
PetscErrorCode DeterminePeriodicity (SimCtx *simCtx)
 Internal helper implementation: DeterminePeriodicity().
 
PetscErrorCode VerifyPathExistence (const char *path, PetscBool is_dir, PetscBool is_optional, const char *description, PetscBool *exists)
 Internal helper implementation: VerifyPathExistence().
 
static PetscErrorCode CheckDataFile (UserCtx *user, PetscInt ti, const char *fieldName, const char *ext, PetscBool *fileExists)
 Internal helper implementation: CheckDataFile().
 
static PetscErrorCode ReadOptionalField (UserCtx *user, const char *field_name, const char *field_label, Vec field_vec, PetscInt ti, const char *ext)
 Internal helper implementation: ReadOptionalField().
 
static PetscErrorCode ReadOptionalSwarmField (UserCtx *user, const char *field_name, const char *field_label, PetscInt ti, const char *ext)
 Internal helper implementation: ReadOptionalSwarmField().
 
PetscErrorCode ReadFieldData (UserCtx *user, const char *field_name, Vec field_vec, PetscInt ti, const char *ext)
 Internal helper implementation: ReadFieldData().
 
PetscErrorCode ReadSimulationFields (UserCtx *user, PetscInt ti)
 Internal helper implementation: ReadSimulationFields().
 
PetscErrorCode ReadStatisticalFields (UserCtx *user, PetscInt ti)
 Implementation of ReadStatisticalFields().
 
PetscErrorCode ReadLESFields (UserCtx *user, PetscInt ti)
 Internal helper implementation: ReadLESFields().
 
PetscErrorCode ReadRANSFields (UserCtx *user, PetscInt ti)
 Internal helper implementation: ReadRANSFields().
 
PetscErrorCode ReadSwarmField (UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
 Internal helper implementation: ReadSwarmField().
 
PetscErrorCode ReadSwarmIntField (UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
 Internal helper implementation: ReadSwarmIntField().
 
PetscErrorCode ReadAllSwarmFields (UserCtx *user, PetscInt ti)
 Internal helper implementation: ReadAllSwarmFields().
 
PetscErrorCode WriteFieldData (UserCtx *user, const char *field_name, Vec field_vec, PetscInt ti, const char *ext)
 Internal helper implementation: WriteFieldData().
 
PetscErrorCode WriteSimulationFields (UserCtx *user)
 Implementation of WriteSimulationFields().
 
PetscErrorCode WriteStatisticalFields (UserCtx *user)
 Implementation of WriteStatisticalFields().
 
PetscErrorCode WriteLESFields (UserCtx *user)
 Internal helper implementation: WriteLESFields().
 
PetscErrorCode WriteRANSFields (UserCtx *user)
 Implementation of WriteRANSFields().
 
PetscErrorCode WriteSwarmField (UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
 Implementation of WriteSwarmField().
 
PetscErrorCode WriteSwarmIntField (UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
 Internal helper implementation: WriteSwarmIntField().
 
PetscErrorCode WriteAllSwarmFields (UserCtx *user)
 Internal helper implementation: WriteAllSwarmFields().
 
PetscErrorCode VecToArrayOnRank0 (Vec inVec, PetscInt *N, double **arrayOut)
 Internal helper implementation: VecToArrayOnRank0().
 
PetscErrorCode SwarmFieldToArrayOnRank0 (DM swarm, const char *field_name, PetscInt *n_total_particles, PetscInt *n_components, void **gathered_array)
 Internal helper implementation: SwarmFieldToArrayOnRank0().
 
PetscErrorCode DisplayBanner (SimCtx *simCtx)
 Internal helper implementation: DisplayBanner().
 
PetscErrorCode ParsePostProcessingSettings (SimCtx *simCtx)
 Internal helper implementation: ParsePostProcessingSettings().
 
PetscErrorCode ParseScalingInformation (SimCtx *simCtx)
 Implementation of ParseScalingInformation().
 
PetscInt ReadDataFileToArray (const char *filename, double **data_out, PetscInt *Nout, MPI_Comm comm)
 Implementation of ReadDataFileToArray().
 
PetscErrorCode ReadPositionsFromFile (PetscInt timeIndex, UserCtx *user, double **coordsArray, PetscInt *Ncoords)
 Internal helper implementation: ReadPositionsFromFile().
 
PetscErrorCode ReadFieldDataToRank0 (PetscInt timeIndex, const char *fieldName, UserCtx *user, double **scalarArray, PetscInt *Nscalars)
 Internal helper implementation: ReadFieldDataToRank0().
 

Variables

static PetscInt g_nblk_from_file = 0
 Stores the number of blocks read from the grid file.
 
static PetscInt * g_IMs_from_file = NULL
 Caches the IM dimensions for all blocks read from the grid file.
 
static PetscInt * g_JMs_from_file = NULL
 Caches the JM dimensions for all blocks read from the grid file.
 
static PetscInt * g_KMs_from_file = NULL
 Caches the KM dimensions for all blocks read from the grid file.
 
static PetscBool g_file_has_been_read = PETSC_FALSE
 A flag to ensure the grid file is read only once.
 

Detailed Description

Implementation of data input/output routines, focusing on grid configuration.

This module provides functions to parse grid geometry information, either from command-line options for programmatically generated grids or by reading the header of a grid definition file.

Definition in file io.c.

Macro Definition Documentation

◆ __FUNCT__ [1/8]

#define __FUNCT__   "ReadGridGenerationInputs"

Definition at line 80 of file io.c.

◆ __FUNCT__ [2/8]

#define __FUNCT__   "ReadGridFile"

Definition at line 80 of file io.c.

◆ __FUNCT__ [3/8]

#define __FUNCT__   "ParseAllBoundaryConditions"

Definition at line 80 of file io.c.

◆ __FUNCT__ [4/8]

#define __FUNCT__   "DeterminePeriodicity"

Definition at line 80 of file io.c.

◆ __FUNCT__ [5/8]

#define __FUNCT__   "ReadFieldData"

Definition at line 80 of file io.c.

◆ __FUNCT__ [6/8]

#define __FUNCT__   "WriteFieldData"

Definition at line 80 of file io.c.

◆ __FUNCT__ [7/8]

#define __FUNCT__   "ParsePostProcessingSettings"

Definition at line 80 of file io.c.

◆ __FUNCT__ [8/8]

#define __FUNCT__   "ParseScalingInformation"

Definition at line 80 of file io.c.

Function Documentation

◆ TrimWhitespace()

void TrimWhitespace ( char *  str)

Implementation of TrimWhitespace().

Helper function to trim leading/trailing whitespace from a string.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
TrimWhitespace()

Definition at line 38 of file io.c.

38 {
39 if (!str) return;
40
41 char *start = str;
42 // Find the first non-whitespace character
43 while (isspace((unsigned char)*start)) {
44 start++;
45 }
46
47 // Find the end of the string
48 char *end = str + strlen(str) - 1;
49 // Move backwards from the end to find the last non-whitespace character
50 while (end > start && isspace((unsigned char)*end)) {
51 end--;
52 }
53
54 // Null-terminate after the last non-whitespace character
55 *(end + 1) = '\0';
56
57 // If there was leading whitespace, shift the string to the left
58 if (str != start) {
59 memmove(str, start, (end - start) + 2); // +2 to include the new null terminator
60 }
61}
Here is the caller graph for this function:

◆ ShouldWriteDataOutput()

PetscBool ShouldWriteDataOutput ( const SimCtx *  simCtx,
PetscInt  completed_step 
)

Implementation of ShouldWriteDataOutput().

Returns whether full field/restart output should be written for the.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
ShouldWriteDataOutput()

Definition at line 70 of file io.c.

71{
72 if (!simCtx) {
73 return PETSC_FALSE;
74 }
75 return (PetscBool)(simCtx->tiout > 0 && completed_step > 0 && completed_step % simCtx->tiout == 0);
76}
SimCtx::tiout
PetscInt tiout
Definition variables.h:696
Here is the caller graph for this function:

◆ ReadGridGenerationInputs()

PetscErrorCode ReadGridGenerationInputs ( UserCtx *  user)

Internal helper implementation: ReadGridGenerationInputs().

Parses command-line options for a programmatically generated grid for a SINGLE block.

Local to this translation unit.

Definition at line 85 of file io.c.

86{
87 PetscErrorCode ierr;
88 SimCtx *simCtx = user->simCtx;
89 PetscInt nblk = simCtx->block_number;
90 PetscInt block_index = user->_this;
91 PetscBool found;
92
93 // Temporary arrays to hold the parsed values for ALL blocks
94 PetscInt *IMs = NULL, *JMs = NULL, *KMs = NULL, *cgrids = NULL;
95 PetscReal *xMins = NULL, *xMaxs = NULL, *rxs = NULL;
96 PetscReal *yMins = NULL, *yMaxs = NULL, *rys = NULL;
97 PetscReal *zMins = NULL, *zMaxs = NULL, *rzs = NULL;
98
99 PetscFunctionBeginUser;
100 PROFILE_FUNCTION_BEGIN;
101
102 LOG_ALLOW(LOCAL, LOG_DEBUG, "Rank %d: Reading generated grid inputs for block %d.\n", simCtx->rank, block_index);
103
104 if (block_index >= nblk) {
105 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Block index %d is out of range for nblk=%d", block_index, nblk);
106 }
107
108 // --- Allocate temporary storage for all array options ---
109 ierr = PetscMalloc4(nblk, &IMs, nblk, &JMs, nblk, &KMs, nblk, &cgrids); CHKERRQ(ierr);
110 ierr = PetscMalloc6(nblk, &xMins, nblk, &xMaxs, nblk, &rxs, nblk, &yMins, nblk, &yMaxs, nblk, &rys); CHKERRQ(ierr);
111 ierr = PetscMalloc3(nblk, &zMins, nblk, &zMaxs, nblk, &rzs); CHKERRQ(ierr);
112
113 // --- Set default values for the temporary arrays ---
114 for (PetscInt i = 0; i < nblk; ++i) {
115 IMs[i] = 10; JMs[i] = 10; KMs[i] = 10; cgrids[i] = 0;
116 xMins[i] = 0.0; xMaxs[i] = 1.0; rxs[i] = 1.0;
117 yMins[i] = 0.0; yMaxs[i] = 1.0; rys[i] = 1.0;
118 zMins[i] = 0.0; zMaxs[i] = 1.0; rzs[i] = 1.0;
119 }
120
121 // --- Parse the array options from the command line / control file ---
122 PetscInt count;
123 count = nblk; ierr = PetscOptionsGetIntArray(NULL, NULL, "-im", IMs, &count, &found); CHKERRQ(ierr);
124 count = nblk; ierr = PetscOptionsGetIntArray(NULL, NULL, "-jm", JMs, &count, &found); CHKERRQ(ierr);
125 count = nblk; ierr = PetscOptionsGetIntArray(NULL, NULL, "-km", KMs, &count, &found); CHKERRQ(ierr);
126 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-xMins", xMins, &count, &found); CHKERRQ(ierr);
127 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-xMaxs", xMaxs, &count, &found); CHKERRQ(ierr);
128 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-rxs", rxs, &count, &found); CHKERRQ(ierr);
129 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-yMins", yMins, &count, &found); CHKERRQ(ierr);
130 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-yMaxs", yMaxs, &count, &found); CHKERRQ(ierr);
131 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-rys", rys, &count, &found); CHKERRQ(ierr);
132 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-zMins", zMins, &count, &found); CHKERRQ(ierr);
133 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-zMaxs", zMaxs, &count, &found); CHKERRQ(ierr);
134 count = nblk; ierr = PetscOptionsGetRealArray(NULL, NULL, "-rzs", rzs, &count, &found); CHKERRQ(ierr);
135 count = nblk; ierr = PetscOptionsGetIntArray(NULL, NULL, "-cgrids", cgrids, &count, &found); CHKERRQ(ierr);
136
137 // --- Assign the parsed values to the specific UserCtx struct passed in ---
138 user->IM = IMs[block_index];
139 user->JM = JMs[block_index];
140 user->KM = KMs[block_index];
141 user->Min_X = xMins[block_index];
142 user->Max_X = xMaxs[block_index];
143 user->rx = rxs[block_index];
144 user->Min_Y = yMins[block_index];
145 user->Max_Y = yMaxs[block_index];
146 user->ry = rys[block_index];
147 user->Min_Z = zMins[block_index];
148 user->Max_Z = zMaxs[block_index];
149 user->rz = rzs[block_index];
150 user->cgrid = cgrids[block_index];
151
152 LOG_ALLOW(LOCAL, LOG_DEBUG, "Rank %d: Block %d grid generation inputs set: IM=%d, JM=%d, KM=%d\n",
153 simCtx->rank, block_index, user->IM, user->JM, user->KM);
154 LOG_ALLOW(LOCAL, LOG_DEBUG, "Rank %d: Block %d bounds: X=[%.2f, %.2f], Y=[%.2f, %.2f], Z=[%.2f, %.2f]\n",
155 simCtx->rank, block_index, user->Min_X, user->Max_X, user->Min_Y, user->Max_Y, user->Min_Z, user->Max_Z);
156
157 // --- Clean up temporary storage ---
158 ierr = PetscFree4(IMs, JMs, KMs, cgrids); CHKERRQ(ierr);
159 ierr = PetscFree6(xMins, xMaxs, rxs, yMins, yMaxs, rys); CHKERRQ(ierr);
160 ierr = PetscFree3(zMins, zMaxs, rzs); CHKERRQ(ierr);
161
162 PROFILE_FUNCTION_END;
163 PetscFunctionReturn(0);
164}
LOCAL
#define LOCAL
Logging scope definitions for controlling message output.
Definition logging.h:44
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:199
PROFILE_FUNCTION_END
#define PROFILE_FUNCTION_END
Marks the end of a profiled code block.
Definition logging.h:827
LOG_DEBUG
@ LOG_DEBUG
Detailed debugging information.
Definition logging.h:31
PROFILE_FUNCTION_BEGIN
#define PROFILE_FUNCTION_BEGIN
Marks the beginning of a profiled code block (typically a function).
Definition logging.h:818
SimCtx::rank
PetscMPIInt rank
Definition variables.h:687
UserCtx::cgrid
PetscInt cgrid
Definition variables.h:891
SimCtx::block_number
PetscInt block_number
Definition variables.h:768
UserCtx::simCtx
SimCtx * simCtx
Back-pointer to the master simulation context.
Definition variables.h:879
UserCtx::Min_X
PetscReal Min_X
Definition variables.h:886
UserCtx::KM
PetscInt KM
Definition variables.h:885
UserCtx::_this
PetscInt _this
Definition variables.h:889
UserCtx::ry
PetscReal ry
Definition variables.h:890
UserCtx::Max_Y
PetscReal Max_Y
Definition variables.h:886
UserCtx::rz
PetscReal rz
Definition variables.h:890
UserCtx::JM
PetscInt JM
Definition variables.h:885
UserCtx::Min_Z
PetscReal Min_Z
Definition variables.h:886
UserCtx::Max_X
PetscReal Max_X
Definition variables.h:886
UserCtx::Min_Y
PetscReal Min_Y
Definition variables.h:886
UserCtx::IM
PetscInt IM
Definition variables.h:885
UserCtx::rx
PetscReal rx
Definition variables.h:890
UserCtx::Max_Z
PetscReal Max_Z
Definition variables.h:886
SimCtx
The master context for the entire simulation.
Definition variables.h:684
Here is the caller graph for this function:

◆ PopulateFinestUserGridResolutionFromOptions()

PetscErrorCode PopulateFinestUserGridResolutionFromOptions ( UserCtx *  finest_users,
PetscInt  nblk 
)

Internal helper implementation: PopulateFinestUserGridResolutionFromOptions().

Parses grid resolution arrays (-im, -jm, -km) once and applies them to all finest-grid blocks.

Local to this translation unit.

Definition at line 170 of file io.c.

171{
172 PetscErrorCode ierr;
173 PetscBool found;
174 PetscInt *IMs = NULL, *JMs = NULL, *KMs = NULL;
175 SimCtx *simCtx = NULL;
176
177 PetscFunctionBeginUser;
178
179 if (!finest_users) {
180 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "finest_users cannot be NULL.");
181 }
182 if (nblk <= 0) {
183 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "nblk must be positive. Got %d.", nblk);
184 }
185 simCtx = finest_users[0].simCtx;
186
187 ierr = PetscMalloc3(nblk, &IMs, nblk, &JMs, nblk, &KMs); CHKERRQ(ierr);
188 for (PetscInt i = 0; i < nblk; ++i) {
189 IMs[i] = 10; JMs[i] = 10; KMs[i] = 10;
190 }
191
192 PetscInt count;
193 count = nblk; ierr = PetscOptionsGetIntArray(NULL, NULL, "-im", IMs, &count, &found); CHKERRQ(ierr);
194 count = nblk; ierr = PetscOptionsGetIntArray(NULL, NULL, "-jm", JMs, &count, &found); CHKERRQ(ierr);
195 count = nblk; ierr = PetscOptionsGetIntArray(NULL, NULL, "-km", KMs, &count, &found); CHKERRQ(ierr);
196
197 for (PetscInt bi = 0; bi < nblk; ++bi) {
198 finest_users[bi].IM = IMs[bi];
199 finest_users[bi].JM = JMs[bi];
200 finest_users[bi].KM = KMs[bi];
201 if (simCtx) {
202 LOG_ALLOW(LOCAL, LOG_DEBUG,
203 "Rank %d: Preloaded analytical grid resolution for block %d: IM=%d, JM=%d, KM=%d\n",
204 simCtx->rank, bi, finest_users[bi].IM, finest_users[bi].JM, finest_users[bi].KM);
205 }
206 }
207
208 ierr = PetscFree3(IMs, JMs, KMs); CHKERRQ(ierr);
209 PetscFunctionReturn(0);
210}
Here is the caller graph for this function:

◆ ReadGridFile()

PetscErrorCode ReadGridFile ( UserCtx *  user)

Internal helper implementation: ReadGridFile().

Sets grid dimensions from a file for a SINGLE block using a one-time read cache.

Local to this translation unit.

Definition at line 219 of file io.c.

220{
221 PetscErrorCode ierr;
222 SimCtx *simCtx = user->simCtx;
223 PetscInt block_index = user->_this;
224
225 PetscFunctionBeginUser;
226 PROFILE_FUNCTION_BEGIN;
227
228 // --- One-Time Read and Broadcast Logic ---
229 if (!g_file_has_been_read) {
230 LOG_ALLOW_SYNC(GLOBAL, LOG_INFO, "First call to ReadGridFile. Reading and broadcasting grid file header from '%s'...\n", simCtx->grid_file);
231 PetscMPIInt rank = simCtx->rank;
232 PetscInt nblk = simCtx->block_number;
233
234 if (rank == 0) {
235 FILE *fd = fopen(simCtx->grid_file, "r");
236 if (!fd) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open file: %s", simCtx->grid_file);
237
238 // Read and validate the canonical PICGRID header.
239 char firstTok[32] = {0};
240 if (fscanf(fd, "%31s", firstTok) != 1)
241 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_READ, "Empty grid file: %s", simCtx->grid_file);
242 if (strcmp(firstTok, "PICGRID") != 0)
243 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_READ,
244 "Grid file %s must begin with the canonical PICGRID header.", simCtx->grid_file);
245 if (fscanf(fd, "%d", &g_nblk_from_file) != 1)
246 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_READ, "Expected number of blocks after \"PICGRID\" in %s", simCtx->grid_file);
247 if (g_nblk_from_file != nblk) {
248 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_UNEXPECTED, "Mismatch: -nblk is %d but grid file specifies %d blocks.", nblk, g_nblk_from_file);
249 }
250
251 ierr = PetscMalloc3(nblk, &g_IMs_from_file, nblk, &g_JMs_from_file, nblk, &g_KMs_from_file); CHKERRQ(ierr);
252 for (PetscInt i = 0; i < nblk; ++i) {
253 if (fscanf(fd, "%d %d %d\n", &g_IMs_from_file[i], &g_JMs_from_file[i], &g_KMs_from_file[i]) != 3) {
254 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_READ, "Expected 3 integers for block %d in %s", i, simCtx->grid_file);
255 }
256 }
257 fclose(fd);
258 }
259
260 // Broadcast nblk to verify (optional, good practice)
261 ierr = MPI_Bcast(&g_nblk_from_file, 1, MPI_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
262
263 // Allocate on other ranks before receiving the broadcast
264 if (rank != 0) {
265 ierr = PetscMalloc3(nblk, &g_IMs_from_file, nblk, &g_JMs_from_file, nblk, &g_KMs_from_file); CHKERRQ(ierr);
266 }
267
268 // Broadcast the data arrays
269 ierr = MPI_Bcast(g_IMs_from_file, nblk, MPI_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
270 ierr = MPI_Bcast(g_JMs_from_file, nblk, MPI_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
271 ierr = MPI_Bcast(g_KMs_from_file, nblk, MPI_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
272
273 g_file_has_been_read = PETSC_TRUE;
274 LOG_ALLOW(GLOBAL, LOG_INFO, "Grid file header read and broadcast complete.\n");
275 }
276
277 // --- Per-Block Assignment Logic (runs on every call) ---
278 user->IM = g_IMs_from_file[block_index];
279 user->JM = g_JMs_from_file[block_index];
280 user->KM = g_KMs_from_file[block_index];
281
282 LOG_ALLOW(LOCAL, LOG_DEBUG, "Rank %d: Set file inputs for Block %d: IM=%d, JM=%d, KM=%d\n",
283 simCtx->rank, block_index, user->IM, user->JM, user->KM);
284
285 PROFILE_FUNCTION_END;
286 PetscFunctionReturn(0);
287}
g_IMs_from_file
static PetscInt * g_IMs_from_file
Caches the IM dimensions for all blocks read from the grid file.
Definition io.c:19
g_nblk_from_file
static PetscInt g_nblk_from_file
Stores the number of blocks read from the grid file.
Definition io.c:17
g_file_has_been_read
static PetscBool g_file_has_been_read
A flag to ensure the grid file is read only once.
Definition io.c:25
g_KMs_from_file
static PetscInt * g_KMs_from_file
Caches the KM dimensions for all blocks read from the grid file.
Definition io.c:23
g_JMs_from_file
static PetscInt * g_JMs_from_file
Caches the JM dimensions for all blocks read from the grid file.
Definition io.c:21
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:252
GLOBAL
#define GLOBAL
Scope for global logging across all processes.
Definition logging.h:45
LOG_INFO
@ LOG_INFO
Informational messages about program execution.
Definition logging.h:30
SimCtx::grid_file
char grid_file[PETSC_MAX_PATH_LEN]
Definition variables.h:773
Here is the caller graph for this function:

◆ FreeBC_ParamList()

void FreeBC_ParamList ( BC_Param *  head)

Implementation of FreeBC_ParamList().

Frees an entire linked list of boundary-condition parameters.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
FreeBC_ParamList()

Definition at line 302 of file io.c.

302 {
303 BC_Param *current = head;
304 while (current != NULL) {
305 BC_Param *next = current->next;
306 PetscFree(current->key);
307 PetscFree(current->value);
308 PetscFree(current);
309 current = next;
310 }
311}
BC_Param_s::next
struct BC_Param_s * next
Definition variables.h:337
BC_Param_s::key
char * key
Definition variables.h:335
BC_Param_s::value
char * value
Definition variables.h:336
BC_Param_s
A node in a linked list for storing key-value parameters from the bcs.dat file.
Definition variables.h:334
Here is the caller graph for this function:

◆ StringToBCFace()

PetscErrorCode StringToBCFace ( const char *  str,
BCFace *  face_out 
)

Internal helper implementation: StringToBCFace().

Converts a face-token string (e.g., "-Xi", "+Eta") to the internal BCFace enum.

Local to this translation unit.

Definition at line 317 of file io.c.

317 {
318 if (strcasecmp(str, "-Xi") == 0) *face_out = BC_FACE_NEG_X;
319 else if (strcasecmp(str, "+Xi") == 0) *face_out = BC_FACE_POS_X;
320 else if (strcasecmp(str, "-Eta") == 0) *face_out = BC_FACE_NEG_Y;
321 else if (strcasecmp(str, "+Eta") == 0) *face_out = BC_FACE_POS_Y;
322 else if (strcasecmp(str, "-Zeta") == 0) *face_out = BC_FACE_NEG_Z;
323 else if (strcasecmp(str, "+Zeta") == 0) *face_out = BC_FACE_POS_Z;
324 else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown face specifier: %s", str);
325 return 0;
326}
BC_FACE_NEG_X
@ BC_FACE_NEG_X
Definition variables.h:260
BC_FACE_POS_Z
@ BC_FACE_POS_Z
Definition variables.h:262
BC_FACE_POS_Y
@ BC_FACE_POS_Y
Definition variables.h:261
BC_FACE_NEG_Z
@ BC_FACE_NEG_Z
Definition variables.h:262
BC_FACE_POS_X
@ BC_FACE_POS_X
Definition variables.h:260
BC_FACE_NEG_Y
@ BC_FACE_NEG_Y
Definition variables.h:261
Here is the caller graph for this function:

◆ StringToBCType()

PetscErrorCode StringToBCType ( const char *  str,
BCType *  type_out 
)

Internal helper implementation: StringToBCType().

Converts a mathematical BC type string (e.g., "PERIODIC", "WALL") to BCType.

Local to this translation unit.

Definition at line 332 of file io.c.

332 {
333 if (strcasecmp(str, "WALL") == 0) *type_out = WALL;
334 else if (strcasecmp(str, "SYMMETRY") == 0) *type_out = SYMMETRY;
335 else if (strcasecmp(str, "INLET") == 0) *type_out = INLET;
336 else if (strcasecmp(str, "OUTLET") == 0) *type_out = OUTLET;
337 else if (strcasecmp(str, "PERIODIC") == 0) *type_out = PERIODIC;
338 // ... add other BCTypes here ...
339 else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown BC Type string: %s", str);
340 return 0;
341}
INLET
@ INLET
Definition variables.h:288
SYMMETRY
@ SYMMETRY
Definition variables.h:286
OUTLET
@ OUTLET
Definition variables.h:287
PERIODIC
@ PERIODIC
Definition variables.h:290
WALL
@ WALL
Definition variables.h:284
Here is the caller graph for this function:

◆ StringToBCHandlerType()

PetscErrorCode StringToBCHandlerType ( const char *  str,
BCHandlerType *  handler_out 
)

Internal helper implementation: StringToBCHandlerType().

Converts a BC handler token (implementation strategy) to BCHandlerType.

Local to this translation unit.

Definition at line 347 of file io.c.

347 {
348 if (strcasecmp(str, "noslip") == 0) *handler_out = BC_HANDLER_WALL_NOSLIP;
349 else if (strcasecmp(str, "constant_velocity") == 0) *handler_out = BC_HANDLER_INLET_CONSTANT_VELOCITY;
350 else if (strcasecmp(str, "conservation") == 0) *handler_out = BC_HANDLER_OUTLET_CONSERVATION;
351 else if (strcasecmp(str, "parabolic") == 0) *handler_out = BC_HANDLER_INLET_PARABOLIC;
352 else if (strcasecmp(str, "prescribed_flow") == 0) *handler_out = BC_HANDLER_INLET_PROFILE_FROM_FILE;
353 else if (strcasecmp(str,"geometric") == 0) *handler_out = BC_HANDLER_PERIODIC_GEOMETRIC;
354 else if (strcasecmp(str,"constant_flux") == 0) *handler_out = BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX;
355 else if (strcasecmp(str,"initial_flux") == 0) *handler_out = BC_HANDLER_PERIODIC_DRIVEN_INITIAL_FLUX;
356 // ... add other BCHandlerTypes here ...
357 else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown BC Handler string: %s", str);
358 return 0;
359}
BC_HANDLER_PERIODIC_GEOMETRIC
@ BC_HANDLER_PERIODIC_GEOMETRIC
Definition variables.h:314
BC_HANDLER_INLET_PARABOLIC
@ BC_HANDLER_INLET_PARABOLIC
Definition variables.h:307
BC_HANDLER_INLET_CONSTANT_VELOCITY
@ BC_HANDLER_INLET_CONSTANT_VELOCITY
Definition variables.h:306
BC_HANDLER_PERIODIC_DRIVEN_INITIAL_FLUX
@ BC_HANDLER_PERIODIC_DRIVEN_INITIAL_FLUX
Definition variables.h:317
BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX
@ BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX
Definition variables.h:316
BC_HANDLER_INLET_PROFILE_FROM_FILE
@ BC_HANDLER_INLET_PROFILE_FROM_FILE
Definition variables.h:308
BC_HANDLER_WALL_NOSLIP
@ BC_HANDLER_WALL_NOSLIP
Definition variables.h:303
BC_HANDLER_OUTLET_CONSERVATION
@ BC_HANDLER_OUTLET_CONSERVATION
Definition variables.h:312
Here is the caller graph for this function:

◆ ValidateBCHandlerForBCType()

PetscErrorCode ValidateBCHandlerForBCType ( BCType  type,
BCHandlerType  handler 
)

Internal helper implementation: ValidateBCHandlerForBCType().

Validates that a selected handler is compatible with a mathematical BC type.

Local to this translation unit.

Definition at line 365 of file io.c.

365 {
366 switch (type) {
367 case OUTLET:
368 if(handler != BC_HANDLER_OUTLET_CONSERVATION) return PETSC_ERR_ARG_WRONG;
369 break;
370 case WALL:
371 if (handler != BC_HANDLER_WALL_NOSLIP && handler != BC_HANDLER_WALL_MOVING) return PETSC_ERR_ARG_WRONG;
372 break;
373 case INLET:
374 if (handler != BC_HANDLER_INLET_CONSTANT_VELOCITY &&
375 handler != BC_HANDLER_INLET_PARABOLIC &&
376 handler != BC_HANDLER_INLET_PROFILE_FROM_FILE) return PETSC_ERR_ARG_WRONG;
377 break;
378 case PERIODIC:
379 if(handler != BC_HANDLER_PERIODIC_GEOMETRIC && handler != BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX && handler != BC_HANDLER_PERIODIC_DRIVEN_INITIAL_FLUX) return PETSC_ERR_ARG_WRONG;
380 // ... add other validation cases here ...
381 default: break;
382 }
383 return 0; // Combination is valid
384}
BC_HANDLER_WALL_MOVING
@ BC_HANDLER_WALL_MOVING
Definition variables.h:304
Here is the caller graph for this function:

◆ GetBCParamReal()

PetscErrorCode GetBCParamReal ( BC_Param *  params,
const char *  key,
PetscReal *  value_out,
PetscBool *  found 
)

Internal helper implementation: GetBCParamReal().

Searches a BC_Param linked list for a key and returns its value as a double.

Local to this translation unit.

Definition at line 390 of file io.c.

390 {
391 *found = PETSC_FALSE;
392 *value_out = 0.0;
393 if (!key) return 0; // No key to search for
394
395 BC_Param *current = params;
396 while (current) {
397 if (strcasecmp(current->key, key) == 0) {
398 *value_out = atof(current->value);
399 *found = PETSC_TRUE;
400 return 0; // Found it, we're done
401 }
402 current = current->next;
403 }
404 return 0; // It's not an error to not find the key.
405}
Here is the caller graph for this function:

◆ GetBCParamBool()

PetscErrorCode GetBCParamBool ( BC_Param *  params,
const char *  key,
PetscBool *  value_out,
PetscBool *  found 
)

Internal helper implementation: GetBCParamBool().

Searches a BC_Param linked list for a key and returns its value as a bool.

Local to this translation unit.

Definition at line 411 of file io.c.

411 {
412 *found = PETSC_FALSE;
413 *value_out = PETSC_FALSE;
414 if (!key) return 0; // No key to search for
415
416 BC_Param *current = params;
417 while (current) {
418 if (strcasecmp(current->key, key) == 0) {
419 // Key was found.
420 *found = PETSC_TRUE;
421
422 // Check the value string. Default to FALSE if the value is NULL or doesn't match a "true" string.
423 if (current->value &&
424 (strcasecmp(current->value, "true") == 0 ||
425 strcmp(current->value, "1") == 0 ||
426 strcasecmp(current->value, "yes") == 0))
427 {
428 *value_out = PETSC_TRUE;
429 } else {
430 *value_out = PETSC_FALSE;
431 }
432 return 0; // Found it, we're done
433 }
434 current = current->next;
435 }
436 return 0; // It's not an error to not find the key.
437}
Here is the caller graph for this function:

◆ ParseAllBoundaryConditions()

PetscErrorCode ParseAllBoundaryConditions ( UserCtx *  user,
const char *  bcs_input_filename 
)

Internal helper implementation: ParseAllBoundaryConditions().

Parses the boundary conditions file to configure the type, handler, and any associated parameters for all 6 global faces of the domain.

Local to this translation unit.

Definition at line 450 of file io.c.

451{
452 PetscErrorCode ierr;
453 PetscMPIInt rank;
454
455 // Temporary storage for rank 0 to build the configuration before broadcasting.
456 BoundaryFaceConfig configs_rank0[6];
457
458 PetscFunctionBeginUser;
459 PROFILE_FUNCTION_BEGIN;
460 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
461
462 if (rank == 0) {
463 FILE *file;
464 char line_buffer[1024];
465
466 // Initialize the temporary config array with safe defaults on rank 0.
467 for (int i = 0; i < 6; i++) {
468 configs_rank0[i].face_id = (BCFace)i;
469 configs_rank0[i].mathematical_type = WALL;
470 configs_rank0[i].handler_type = BC_HANDLER_WALL_NOSLIP;
471 configs_rank0[i].params = NULL;
472 configs_rank0[i].handler = NULL; // Handler object is not created here.
473 }
474
475 LOG_ALLOW(GLOBAL, LOG_INFO, "Parsing BC configuration from '%s' on rank 0... \n", bcs_input_filename);
476 file = fopen(bcs_input_filename, "r");
477 if (!file) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Could not open BCs file '%s'.", bcs_input_filename);
478
479 while (fgets(line_buffer, sizeof(line_buffer), file)) {
480 char *current_pos = line_buffer;
481 while (isspace((unsigned char)*current_pos)) current_pos++; // Skip leading whitespace
482 if (*current_pos == '#' || *current_pos == '\0' || *current_pos == '\n' || *current_pos == '\r') continue;
483
484 char *face_str = strtok(current_pos, " \t\n\r");
485 char *type_str = strtok(NULL, " \t\n\r");
486 char *handler_str = strtok(NULL, " \t\n\r");
487
488 if (!face_str || !type_str || !handler_str) {
489 LOG_ALLOW(GLOBAL, LOG_WARNING, "Malformed line in bcs.dat, skipping: %s \n", line_buffer);
490 continue;
491 }
492
493 BCFace face_enum;
494 BCType type_enum;
495 BCHandlerType handler_enum;
496 const char* handler_name_for_log;
497
498 // --- Convert strings to enums and validate ---
499 ierr = StringToBCFace(face_str, &face_enum); CHKERRQ(ierr);
500 ierr = StringToBCType(type_str, &type_enum); CHKERRQ(ierr);
501 ierr = StringToBCHandlerType(handler_str, &handler_enum); CHKERRQ(ierr);
502 ierr = ValidateBCHandlerForBCType(type_enum, handler_enum);
503 if (ierr) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Validation failed: Handler '%s' is not valid for Type '%s' on Face '%s'.\n", handler_str, type_str, face_str);
504
505 // Store the core types for the corresponding face
506 configs_rank0[face_enum].mathematical_type = type_enum;
507 configs_rank0[face_enum].handler_type = handler_enum;
508 handler_name_for_log = BCHandlerTypeToString(handler_enum); // Assumes this utility exists
509 LOG_ALLOW(GLOBAL, LOG_DEBUG, " Parsed Face '%s': Type=%s, Handler=%s \n", face_str, type_str, handler_name_for_log);
510
511 // --- Parse optional key=value parameters for this face ---
512 FreeBC_ParamList(configs_rank0[face_enum].params); // Clear any previous (default) params
513 configs_rank0[face_enum].params = NULL;
514 BC_Param **param_next_ptr = &configs_rank0[face_enum].params; // Pointer to the 'next' pointer to build the list
515
516 char* token;
517 while ((token = strtok(NULL, " \t\n\r")) != NULL) {
518 char* equals_ptr = strchr(token, '=');
519 if (!equals_ptr) {
520 LOG_ALLOW(GLOBAL, LOG_WARNING, "Malformed parameter '%s' on face '%s', skipping. \n", token, face_str);
521 continue;
522 }
523
524 *equals_ptr = '\0'; // Temporarily split the string at '=' to separate key and value
525 char* key_str = token;
526 char* value_str = equals_ptr + 1;
527
528 BC_Param *new_param;
529 ierr = PetscMalloc1(1, &new_param); CHKERRQ(ierr);
530 ierr = PetscStrallocpy(key_str, &new_param->key); CHKERRQ(ierr);
531 ierr = PetscStrallocpy(value_str, &new_param->value); CHKERRQ(ierr);
532 new_param->next = NULL;
533
534 *param_next_ptr = new_param;
535 param_next_ptr = &new_param->next;
536 LOG_ALLOW(GLOBAL, LOG_TRACE, " - Found param: [%s] = [%s] \n", new_param->key, new_param->value);
537 }
538 }
539 fclose(file);
540 }
541
542 // =========================================================================
543 // BROADCASTING THE CONFIGURATION FROM RANK 0
544 // =========================================================================
545 // This is a critical step to ensure all processes have the same configuration.
546
547 LOG_ALLOW_SYNC(GLOBAL, LOG_DEBUG, "Rank %d broadcasting/receiving BC configuration.\n", rank);
548
549 for (int i = 0; i < 6; i++) {
550 // --- Broadcast simple enums ---
551 if (rank == 0) {
552 user->boundary_faces[i] = configs_rank0[i]; // Rank 0 populates its final struct
553 }
554 ierr = MPI_Bcast(&user->boundary_faces[i].mathematical_type, 1, MPI_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
555 ierr = MPI_Bcast(&user->boundary_faces[i].handler_type, 1, MPI_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
556
557 // --- Serialize and Broadcast the parameter linked list ---
558 PetscInt n_params = 0;
559 if (rank == 0) { // On rank 0, count the number of parameters to send
560 for (BC_Param *p = user->boundary_faces[i].params; p; p = p->next) n_params++;
561 }
562 ierr = MPI_Bcast(&n_params, 1, MPI_INT, 0, PETSC_COMM_WORLD);CHKERRQ(ierr);
563
564 if (rank != 0) { // Non-root ranks need to receive and build the list
565 FreeBC_ParamList(user->boundary_faces[i].params); // Ensure list is empty before building
566 user->boundary_faces[i].params = NULL;
567 }
568
569 BC_Param **param_next_ptr = &user->boundary_faces[i].params;
570
571 for (int j = 0; j < n_params; j++) {
572 char key_buf[256] = {0}, val_buf[256] = {0};
573 if (rank == 0) {
574 // On rank 0, navigate to the j-th param and copy its data to buffers
575 BC_Param *p = user->boundary_faces[i].params;
576 for (int k = 0; k < j; k++) p = p->next;
577 strncpy(key_buf, p->key, 255);
578 strncpy(val_buf, p->value, 255);
579 }
580
581 ierr = MPI_Bcast(key_buf, 256, MPI_CHAR, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
582 ierr = MPI_Bcast(val_buf, 256, MPI_CHAR, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
583
584 if (rank != 0) {
585 // On non-root ranks, deserialize: create a new node and append it
586 BC_Param *new_param;
587 ierr = PetscMalloc1(1, &new_param); CHKERRQ(ierr);
588 ierr = PetscStrallocpy(key_buf, &new_param->key); CHKERRQ(ierr);
589 ierr = PetscStrallocpy(val_buf, &new_param->value); CHKERRQ(ierr);
590 new_param->next = NULL;
591 *param_next_ptr = new_param;
592 param_next_ptr = &new_param->next;
593 } else {
594 // On rank 0, just advance the pointer for the next iteration
595 param_next_ptr = &((*param_next_ptr)->next);
596 }
597 }
598 user->boundary_faces[i].face_id = (BCFace)i; // Ensure face_id is set on all ranks
599 }
600
601 // --- Set particle inlet lookup fields used by the particle system ---
602 user->inletFaceDefined = PETSC_FALSE;
603 for (int i=0; i<6; i++) {
604
605 if (user->boundary_faces[i].mathematical_type == INLET && !user->inletFaceDefined) {
606 user->inletFaceDefined = PETSC_TRUE;
607 user->identifiedInletBCFace = (BCFace)i;
608 LOG_ALLOW(GLOBAL, LOG_INFO, "Inlet face for particle initialization identified as Face %d.\n", i);
609 break; // Found the first one, stop looking
610 }
611 }
612
613
614 if (rank == 0) {
615 // Rank 0 can now free the linked lists it created for the temporary storage.
616 // As written, user->boundary_faces was populated directly on rank 0, so no extra free is needed.
617 // for(int i=0; i<6; i++) FreeBC_ParamList(configs_rank0[i].params); // This would be needed if we used configs_rank0 exclusively
618 }
619
620 PROFILE_FUNCTION_END;
621 PetscFunctionReturn(0);
622}
StringToBCHandlerType
PetscErrorCode StringToBCHandlerType(const char *str, BCHandlerType *handler_out)
Internal helper implementation: StringToBCHandlerType().
Definition io.c:347
ValidateBCHandlerForBCType
PetscErrorCode ValidateBCHandlerForBCType(BCType type, BCHandlerType handler)
Internal helper implementation: ValidateBCHandlerForBCType().
Definition io.c:365
StringToBCFace
PetscErrorCode StringToBCFace(const char *str, BCFace *face_out)
Internal helper implementation: StringToBCFace().
Definition io.c:317
FreeBC_ParamList
void FreeBC_ParamList(BC_Param *head)
Implementation of FreeBC_ParamList().
Definition io.c:302
StringToBCType
PetscErrorCode StringToBCType(const char *str, BCType *type_out)
Internal helper implementation: StringToBCType().
Definition io.c:332
BCHandlerTypeToString
const char * BCHandlerTypeToString(BCHandlerType handler_type)
Converts a BCHandlerType enum to its string representation.
Definition logging.c:792
LOG_TRACE
@ LOG_TRACE
Very fine-grained tracing information for in-depth debugging.
Definition logging.h:32
LOG_WARNING
@ LOG_WARNING
Non-critical issues that warrant attention.
Definition logging.h:29
BCType
BCType
Defines the general mathematical/physical Category of a boundary.
Definition variables.h:281
UserCtx::inletFaceDefined
PetscBool inletFaceDefined
Definition variables.h:897
UserCtx::boundary_faces
BoundaryFaceConfig boundary_faces[6]
Definition variables.h:896
UserCtx::identifiedInletBCFace
BCFace identifiedInletBCFace
Definition variables.h:898
BCHandlerType
BCHandlerType
Defines the specific computational "strategy" for a boundary handler.
Definition variables.h:301
BoundaryFaceConfig::handler_type
BCHandlerType handler_type
Definition variables.h:367
BoundaryFaceConfig::params
BC_Param * params
Definition variables.h:368
BoundaryFaceConfig::mathematical_type
BCType mathematical_type
Definition variables.h:366
BCFace
BCFace
Identifies the six logical faces of a structured computational block.
Definition variables.h:259
BoundaryFaceConfig::handler
BoundaryCondition * handler
Definition variables.h:369
BoundaryFaceConfig
Holds the complete configuration for one of the six boundary faces.
Definition variables.h:364
Here is the call graph for this function:
Here is the caller graph for this function:

◆ DeterminePeriodicity()

PetscErrorCode DeterminePeriodicity ( SimCtx *  simCtx)

Internal helper implementation: DeterminePeriodicity().

Scans all block-specific boundary condition files to determine a globally consistent periodicity for each dimension, reusing the core type parser.

Local to this translation unit.

Definition at line 637 of file io.c.

638{
639 PetscErrorCode ierr;
640 PetscMPIInt rank;
641 PetscInt periodic_flags[3] = {0, 0, 0}; // Index 0:I, 1:J, 2:K
642
643 PetscFunctionBeginUser;
644 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
645
646 // --- Part 1: Collectively verify all BCS files exist before proceeding ---
647 for (PetscInt bi = 0; bi < simCtx->block_number; bi++) {
648 const char *bcs_filename = simCtx->bcs_files[bi];
649 if (!bcs_filename) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_NULL, "BCS filename for block %d is not set in SimCtx.", bi);
650 char desc_buf[256];
651 PetscBool file_exists;
652 snprintf(desc_buf, sizeof(desc_buf), "BCS file for block %d", bi);
653 ierr = VerifyPathExistence(bcs_filename, PETSC_FALSE, PETSC_FALSE, desc_buf, &file_exists); CHKERRQ(ierr);
654 }
655
656 // --- Part 2: Rank 0 does the parsing, since we know all files exist ---
657 if (rank == 0) {
658 PetscBool global_is_periodic[3] = {PETSC_FALSE, PETSC_FALSE, PETSC_FALSE};
659 PetscBool is_set = PETSC_FALSE;
660
661 for (PetscInt bi = 0; bi < simCtx->block_number; bi++) {
662 const char *bcs_filename = simCtx->bcs_files[bi];
663 FILE *file = fopen(bcs_filename, "r");
664
665 PetscBool face_is_periodic[6] = {PETSC_FALSE};
666 char line_buffer[1024];
667
668 while (fgets(line_buffer, sizeof(line_buffer), file)) {
669 char *current_pos = line_buffer;
670 while (isspace((unsigned char)*current_pos)) current_pos++;
671 if (*current_pos == '#' || *current_pos == '\0' || *current_pos == '\n') continue;
672
673 // --- Tokenize the line exactly like the main parser ---
674 char *face_str = strtok(current_pos, " \t\n\r");
675 char *type_str = strtok(NULL, " \t\n\r");
676
677 // If the line doesn't have at least two tokens, we can't determine the type.
678 if (!face_str || !type_str) continue;
679
680 // --- Perform a direct, non-erroring check on the mathematical type string ---
681 if (strcasecmp(type_str, "PERIODIC") == 0) {
682 BCFace face_enum;
683 // A malformed face string on a periodic line IS a fatal error.
684 ierr = StringToBCFace(face_str, &face_enum); CHKERRQ(ierr);
685 face_is_periodic[face_enum] = PETSC_TRUE;
686 }
687 // Any other type_str (e.g., "WALL", "INLET") is correctly and silently ignored.
688 }
689 fclose(file);
690
691 // --- Validate consistency within this file ---
692 if (face_is_periodic[BC_FACE_NEG_X] != face_is_periodic[BC_FACE_POS_X])
693 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Inconsistent X-periodicity in file '%s' for block %d. Both -Xi and +Xi must be periodic or neither.", bcs_filename, bi);
694 if (face_is_periodic[BC_FACE_NEG_Y] != face_is_periodic[BC_FACE_POS_Y])
695 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Inconsistent Y-periodicity in file '%s' for block %d. Both -Eta and +Eta must be periodic or neither.", bcs_filename, bi);
696 if (face_is_periodic[BC_FACE_NEG_Z] != face_is_periodic[BC_FACE_POS_Z])
697 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Inconsistent Z-periodicity in file '%s' for block %d. Both -Zeta and +Zeta must be periodic or neither.", bcs_filename, bi);
698
699 PetscBool local_is_periodic[3] = {face_is_periodic[BC_FACE_NEG_X], face_is_periodic[BC_FACE_NEG_Y], face_is_periodic[BC_FACE_NEG_Z]};
700
701 // --- Validate consistency across block files ---
702 if (!is_set) {
703 global_is_periodic[0] = local_is_periodic[0];
704 global_is_periodic[1] = local_is_periodic[1];
705 global_is_periodic[2] = local_is_periodic[2];
706 is_set = PETSC_TRUE;
707 } else {
708 if (global_is_periodic[0] != local_is_periodic[0] || global_is_periodic[1] != local_is_periodic[1] || global_is_periodic[2] != local_is_periodic[2]) {
709 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP,
710 "Periodicity mismatch between blocks. Block 0 requires (I:%d, J:%d, K:%d), but block %d (file '%s') has (I:%d, J:%d, K:%d).",
711 (int)global_is_periodic[0], (int)global_is_periodic[1], (int)global_is_periodic[2],
712 bi, bcs_filename,
713 (int)local_is_periodic[0], (int)local_is_periodic[1], (int)local_is_periodic[2]);
714 }
715 }
716 } // end loop over blocks
717
718 periodic_flags[0] = (global_is_periodic[0]) ? 1 : 0;
719 periodic_flags[1] = (global_is_periodic[1]) ? 1 : 0;
720 periodic_flags[2] = (global_is_periodic[2]) ? 1 : 0;
721
722 LOG_ALLOW(GLOBAL, LOG_INFO, "Global periodicity determined: I-periodic=%d, J-periodic=%d, K-periodic=%d\n",
723 periodic_flags[0], periodic_flags[1], periodic_flags[2]);
724 }
725
726 // --- Part 3: Broadcast the final flags from rank 0 to all other ranks ---
727 ierr = MPI_Bcast(periodic_flags, 3, MPIU_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
728
729 // --- All ranks now update their SimCtx ---
730 simCtx->i_periodic = periodic_flags[0];
731 simCtx->j_periodic = periodic_flags[1];
732 simCtx->k_periodic = periodic_flags[2];
733
734 PetscFunctionReturn(0);
735}
VerifyPathExistence
PetscErrorCode VerifyPathExistence(const char *path, PetscBool is_dir, PetscBool is_optional, const char *description, PetscBool *exists)
Internal helper implementation: VerifyPathExistence().
Definition io.c:741
SimCtx::bcs_files
char ** bcs_files
Definition variables.h:776
SimCtx::k_periodic
PetscInt k_periodic
Definition variables.h:769
SimCtx::i_periodic
PetscInt i_periodic
Definition variables.h:769
SimCtx::j_periodic
PetscInt j_periodic
Definition variables.h:769
Here is the call graph for this function:
Here is the caller graph for this function:

◆ VerifyPathExistence()

PetscErrorCode VerifyPathExistence ( const char *  path,
PetscBool  is_dir,
PetscBool  is_optional,
const char *  description,
PetscBool *  exists 
)

Internal helper implementation: VerifyPathExistence().

A parallel-safe helper to verify the existence of a generic file or directory path.

Local to this translation unit.

Definition at line 741 of file io.c.

742{
743 PetscErrorCode ierr;
744 PetscMPIInt rank;
745 MPI_Comm comm = PETSC_COMM_WORLD;
746
747 PetscFunctionBeginUser;
748 ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr);
749
750 if (rank == 0) {
751 if (is_dir) {
752 ierr = PetscTestDirectory(path, 'r', exists); CHKERRQ(ierr);
753 } else {
754 ierr = PetscTestFile(path, 'r', exists); CHKERRQ(ierr);
755 }
756
757 if (!(*exists)) {
758 if (is_optional) {
759 LOG_ALLOW(GLOBAL, LOG_WARNING, "Optional %s not found at: %s (using defaults/ignoring).\n", description, path);
760 } else {
761 LOG_ALLOW(GLOBAL, LOG_ERROR, "Mandatory %s not found at: %s\n", description, path);
762 }
763 } else {
764 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Found %s: %s\n", description, path);
765 }
766 }
767
768 // Broadcast the result from Rank 0
769 PetscMPIInt exists_int = (rank == 0) ? (PetscMPIInt)(*exists) : 0;
770 ierr = MPI_Bcast(&exists_int, 1, MPI_INT, 0, comm); CHKERRMPI(ierr);
771 *exists = (PetscBool)exists_int;
772
773 // Collective error for mandatory files
774 if (!(*exists) && !is_optional) {
775 SETERRQ(comm, PETSC_ERR_FILE_OPEN, "Mandatory %s not found. Rank 0 expected it at '%s'. Check path and permissions.", description, path);
776 }
777
778 PetscFunctionReturn(0);
779}
LOG_ERROR
@ LOG_ERROR
Critical errors that may halt the program.
Definition logging.h:28
Here is the caller graph for this function:

◆ CheckDataFile()

static PetscErrorCode CheckDataFile ( UserCtx *  user,
PetscInt  ti,
const char *  fieldName,
const char *  ext,
PetscBool *  fileExists 
)
static

Internal helper implementation: CheckDataFile().

Local to this translation unit.

Definition at line 785 of file io.c.

786{
787 PetscErrorCode ierr;
788 PetscMPIInt rank;
789 MPI_Comm comm = PETSC_COMM_WORLD;
790 PetscInt placeholder_int = 0;
791 SimCtx *simCtx = user->simCtx;
792
793 PetscFunctionBeginUser;
794
795 if(!simCtx->current_io_directory) {
796 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_NULL, "I/O context directory is NULL. Ensure it is set before calling CheckDataFile().");
797 }
798
799 ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr);
800
801 if (rank == 0) {
802 char filename[PETSC_MAX_PATH_LEN];
803 // Use the same standardized, rank-independent filename format
804
805 ierr = PetscSNPrintf(filename, sizeof(filename), "%s/%s%05"PetscInt_FMT"_%d.%s",simCtx->current_io_directory,fieldName, ti, placeholder_int, ext);
806 ierr = PetscTestFile(filename, 'r', fileExists); CHKERRQ(ierr);
807 if (!(*fileExists)) {
808 LOG_ALLOW(GLOBAL, LOG_WARNING, "(Rank 0) - Optional data file '%s' not found.\n", filename);
809 }
810 }
811
812 // Broadcast the result from Rank 0 to all other ranks.
813 // We cast the PetscBool to a PetscMPIInt for MPI_Bcast.
814 PetscMPIInt fileExists_int = (rank == 0) ? (PetscMPIInt)(*fileExists) : 0;
815 ierr = MPI_Bcast(&fileExists_int, 1, MPI_INT, 0, comm); CHKERRMPI(ierr);
816 *fileExists = (PetscBool)fileExists_int;
817
818 PetscFunctionReturn(0);
819}
SimCtx::current_io_directory
char * current_io_directory
Definition variables.h:711
Here is the caller graph for this function:

◆ ReadOptionalField()

static PetscErrorCode ReadOptionalField ( UserCtx *  user,
const char *  field_name,
const char *  field_label,
Vec  field_vec,
PetscInt  ti,
const char *  ext 
)
static

Internal helper implementation: ReadOptionalField().

Local to this translation unit.

Definition at line 825 of file io.c.

826{
827 PetscErrorCode ierr;
828 PetscBool fileExists;
829
830 PetscFunctionBeginUser;
831
832 /* Check if the data file for this optional field exists. */
833 ierr = CheckDataFile(user,ti, field_name, ext, &fileExists); CHKERRQ(ierr);
834
835 if (fileExists) {
836 /* File exists, so we MUST be able to read it. */
837 LOG_ALLOW(GLOBAL, LOG_DEBUG, "File for %s found, attempting to read...\n", field_label);
838 ierr = ReadFieldData(user, field_name, field_vec, ti, ext);
839
840 if (ierr) {
841 /* Any error here is fatal. A PETSC_ERR_FILE_OPEN would mean a race condition or
842 permissions issue. Other errors could be size mismatch or corruption. */
843 SETERRQ(PETSC_COMM_WORLD, ierr, "Failed to read data for %s from existing file for step %d. The file may be corrupt or have an incorrect size.", field_label, ti);
844 } else {
845 LOG_ALLOW(GLOBAL, LOG_INFO, "Successfully read %s field for step %d.\n", field_label, ti);
846 }
847 } else {
848 /* File does not exist, which is acceptable for an optional field. */
849 LOG_ALLOW(GLOBAL, LOG_WARNING, "Optional %s file for step %d not found. Skipping.\n", field_label, ti);
850 }
851
852 PetscFunctionReturn(0);
853}
ReadFieldData
PetscErrorCode ReadFieldData(UserCtx *user, const char *field_name, Vec field_vec, PetscInt ti, const char *ext)
Internal helper implementation: ReadFieldData().
Definition io.c:901
CheckDataFile
static PetscErrorCode CheckDataFile(UserCtx *user, PetscInt ti, const char *fieldName, const char *ext, PetscBool *fileExists)
Internal helper implementation: CheckDataFile().
Definition io.c:785
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadOptionalSwarmField()

static PetscErrorCode ReadOptionalSwarmField ( UserCtx *  user,
const char *  field_name,
const char *  field_label,
PetscInt  ti,
const char *  ext 
)
static

Internal helper implementation: ReadOptionalSwarmField().

Local to this translation unit.

Definition at line 859 of file io.c.

860{
861 PetscErrorCode ierr;
862 PetscBool fileExists;
863
864 PetscFunctionBeginUser;
865
866 /* Check if the data file for this optional field exists. */
867 ierr = CheckDataFile(user,ti, field_name, ext, &fileExists); CHKERRQ(ierr);
868
869 if (fileExists) {
870 /* File exists, so we MUST be able to read it. */
871 LOG_ALLOW(GLOBAL, LOG_DEBUG, "File for %s found, attempting to read...\n", field_label);
872 if(strcasecmp(field_name,"DMSwarm_CellID") == 0 || strcasecmp(field_name,"DMSwarm_pid")== 0 || strcasecmp(field_name,"DMSwarm_location_status")== 0 ) {
873 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Reading integer swarm field '%s'.\n", field_name);
874 ierr = ReadSwarmIntField(user,field_name,ti,ext);
875 }
876 else{
877 ierr = ReadSwarmField(user, field_name, ti, ext);
878 }
879
880 if (ierr) {
881 /* Any error here is fatal. A PETSC_ERR_FILE_OPEN would mean a race condition or
882 permissions issue. Other errors could be size mismatch or corruption. */
883 SETERRQ(PETSC_COMM_WORLD, ierr, "Failed to read data for %s from existing file for step %d. The file may be corrupt or have an incorrect size.", field_label, ti);
884 } else {
885 LOG_ALLOW(GLOBAL, LOG_INFO, "Successfully read %s field for step %d.\n", field_label, ti);
886 }
887 } else {
888 /* File does not exist, which is acceptable for an optional field. */
889 LOG_ALLOW(GLOBAL, LOG_WARNING, "Optional %s file for step %d not found. Skipping.\n", field_label, ti);
890 }
891
892 PetscFunctionReturn(0);
893}
ReadSwarmIntField
PetscErrorCode ReadSwarmIntField(UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
Internal helper implementation: ReadSwarmIntField().
Definition io.c:1316
ReadSwarmField
PetscErrorCode ReadSwarmField(UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
Internal helper implementation: ReadSwarmField().
Definition io.c:1286
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadFieldData()

PetscErrorCode ReadFieldData ( UserCtx *  user,
const char *  field_name,
Vec  field_vec,
PetscInt  ti,
const char *  ext 
)

Internal helper implementation: ReadFieldData().

Reads data for a specific field from a file into the provided vector.

Local to this translation unit.

Definition at line 901 of file io.c.

906{
907 PetscErrorCode ierr;
908 char filename[PETSC_MAX_PATH_LEN];
909 MPI_Comm comm;
910 PetscMPIInt rank,size;
911 SimCtx *simCtx = user->simCtx;
912
913
914 PetscFunctionBeginUser;
915 PROFILE_FUNCTION_BEGIN;
916
917 if(!simCtx->current_io_directory){
918 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "I/O context directory was not set before calling ReadFieldData().");
919 }
920
921
922 ierr = PetscObjectGetComm((PetscObject)field_vec,&comm);CHKERRQ(ierr);
923 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr);
924 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr);
925
926 const char *source_path = NULL;
927 source_path = simCtx->current_io_directory;
928
929 if(!source_path){
930 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "source_path was not set for the current execution mode.");
931 }
932 /* ---------------------------------------------------------------------
933 * Compose <path-to-file>/<field_name><step with 5 digits>_0.<ext>
934 * (all restart files are written by rank-0 with that naming scheme).
935 * ------------------------------------------------------------------ */
936 ierr = PetscSNPrintf(filename,sizeof(filename),
937 "%s/%s%05" PetscInt_FMT "_0.%s",
938 source_path,field_name,ti,ext);CHKERRQ(ierr);
939
940 LOG_ALLOW(GLOBAL,LOG_DEBUG,
941 "Attempting to read <%s> on rank %d/%d\n",
942 filename,(int)rank,(int)size);
943
944 /* ======================================================================
945 * 1. SERIAL JOB – just hand the Vec to VecLoad()
946 * ==================================================================== */
947 if(size==1)
948 {
949 PetscViewer viewer;
950 PetscBool found;
951 Vec temp_vec;
952 PetscInt expectedSize,loadedSize;
953
954 ierr = PetscTestFile(filename,'r',&found);CHKERRQ(ierr);
955 if(!found) SETERRQ(comm,PETSC_ERR_FILE_OPEN,
956 "Restart/Source file not found: %s",filename);
957
958 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,filename,FILE_MODE_READ,&viewer);CHKERRQ(ierr);
959// ---- START MODIFICATION ----
960 // DO NOT load directly into field_vec, as this can resize it, which is
961 // illegal for DMSwarm "view" vectors. Instead, load into a temporary vector.
962 ierr = VecCreate(PETSC_COMM_SELF, &temp_vec); CHKERRQ(ierr);
963 ierr = VecLoad(temp_vec,viewer);CHKERRQ(ierr);
964 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
965
966 // Sanity check: ensure the file size matches the expected vector size.
967 ierr = VecGetSize(field_vec, &expectedSize);CHKERRQ(ierr);
968 ierr = VecGetSize(temp_vec, &loadedSize);CHKERRQ(ierr);
969 if (loadedSize != expectedSize) {
970 SETERRQ(comm,PETSC_ERR_FILE_UNEXPECTED,
971 "File %s holds %d entries – expected %d for field '%s'",
972 filename, loadedSize, expectedSize, field_name);
973 }
974
975 // Now, safely copy the data from the temporary vector to the final destination.
976 ierr = VecCopy(temp_vec, field_vec);CHKERRQ(ierr);
977
978 // Clean up the temporary vector.
979 ierr = VecDestroy(&temp_vec);CHKERRQ(ierr);
980
981 // ---- END MODIFICATION ----
982
983 /* create EMPTY sequential Vec – VecLoad() will size it correctly */
984 /*
985 ierr = VecCreate(PETSC_COMM_SELF,&seq_vec);CHKERRQ(ierr);
986 ierr = VecSetType(seq_vec,VECSEQ);CHKERRQ(ierr);
987
988 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,filename,
989 FILE_MODE_READ,&viewer);CHKERRQ(ierr);
990
991 ierr = VecLoad(field_vec,viewer);CHKERRQ(ierr);
992 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
993 */
994 LOG_ALLOW(GLOBAL,LOG_INFO,
995 "Loaded <%s> (serial path)\n",filename);
996
997 PROFILE_FUNCTION_END;
998 PetscFunctionReturn(0);
999 }
1000
1001 /* ======================================================================
1002 * 2. PARALLEL JOB
1003 * ==================================================================== */
1004 PetscInt globalSize;
1005 ierr = VecGetSize(field_vec,&globalSize);CHKERRQ(ierr);
1006
1007 DM dm = NULL;
1008 const char *dmtype = NULL;
1009 Vec nat = NULL; /* Natural-ordered vector for DMDA */
1010
1011 /* -------------------- rank-0 : read the sequential file -------------- */
1012 Vec seq_vec = NULL; /* only valid on rank-0 */
1013 const PetscScalar *seqArray = NULL; /* borrowed pointer on rank-0 only */
1014
1015 if(rank==0)
1016 {
1017 PetscViewer viewer;
1018 PetscBool found;
1019
1020 ierr = PetscTestFile(filename,'r',&found);CHKERRQ(ierr);
1021 if(!found) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,
1022 "Restart file not found: %s",filename);
1023
1024 /* create EMPTY sequential Vec – VecLoad() will size it correctly */
1025 ierr = VecCreate(PETSC_COMM_SELF,&seq_vec);CHKERRQ(ierr);
1026 ierr = VecSetType(seq_vec,VECSEQ);CHKERRQ(ierr);
1027
1028 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,filename,
1029 FILE_MODE_READ,&viewer);CHKERRQ(ierr);
1030 ierr = VecLoad(seq_vec,viewer);CHKERRQ(ierr);
1031 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
1032
1033 /* size sanity-check */
1034 PetscInt loaded;
1035 ierr = VecGetSize(seq_vec,&loaded);CHKERRQ(ierr);
1036 if(loaded != globalSize)
1037 SETERRQ(comm,PETSC_ERR_FILE_UNEXPECTED,
1038 "File %s holds %d entries – expected %d",
1039 filename,loaded,globalSize);
1040
1041 /* borrow array for later Bcast */
1042 ierr = VecGetArrayRead(seq_vec,&seqArray);CHKERRQ(ierr);
1043
1044 LOG_ALLOW(GLOBAL,LOG_TRACE,
1045 "Rank 0 successfully loaded <%s>\n",filename);
1046 }
1047
1048 /* -------------------- Check if this is a DMDA vector ----------------- */
1049 ierr = VecGetDM(field_vec, &dm); CHKERRQ(ierr);
1050 if (dm) { ierr = DMGetType(dm, &dmtype); CHKERRQ(ierr); }
1051
1052 if (dmtype && !strcmp(dmtype, DMDA)) {
1053 /* ==================================================================
1054 * DMDA PATH: File is in natural ordering, need to convert to global
1055 * ================================================================== */
1056
1057 /* Create natural vector */
1058 ierr = DMDACreateNaturalVector(dm, &nat); CHKERRQ(ierr);
1059
1060 /* Scatter from rank 0's seq_vec to all ranks' natural vector */
1061 VecScatter scatter;
1062 Vec nat_seq = NULL; /* Sequential natural vector on rank 0 */
1063
1064 ierr = VecScatterCreateToZero(nat, &scatter, &nat_seq); CHKERRQ(ierr);
1065
1066 /* Reverse scatter: from rank 0 to all ranks */
1067 ierr = VecScatterBegin(scatter, (rank == 0 ? seq_vec : nat_seq), nat,
1068 INSERT_VALUES, SCATTER_REVERSE); CHKERRQ(ierr);
1069 ierr = VecScatterEnd(scatter, (rank == 0 ? seq_vec : nat_seq), nat,
1070 INSERT_VALUES, SCATTER_REVERSE); CHKERRQ(ierr);
1071
1072 /* Convert natural → global ordering */
1073 ierr = DMDANaturalToGlobalBegin(dm, nat, INSERT_VALUES, field_vec); CHKERRQ(ierr);
1074 ierr = DMDANaturalToGlobalEnd(dm, nat, INSERT_VALUES, field_vec); CHKERRQ(ierr);
1075
1076 /* Cleanup */
1077 ierr = VecScatterDestroy(&scatter); CHKERRQ(ierr);
1078 ierr = VecDestroy(&nat_seq); CHKERRQ(ierr);
1079 ierr = VecDestroy(&nat); CHKERRQ(ierr);
1080
1081 } else {
1082 /* ==================================================================
1083 * NON-DMDA PATH: Use broadcast and direct copy (assumes global ordering)
1084 * ================================================================== */
1085
1086 PetscScalar *buffer = NULL;
1087 if (rank == 0) {
1088 buffer = (PetscScalar *)seqArray;
1089 } else {
1090 ierr = PetscMalloc1(globalSize, &buffer); CHKERRQ(ierr);
1091 }
1092
1093 ierr = MPI_Bcast(buffer, (int)globalSize, MPIU_SCALAR, 0, comm); CHKERRQ(ierr);
1094
1095 /* Copy slice based on ownership range */
1096 PetscInt rstart, rend, loc;
1097 PetscScalar *locArray;
1098
1099 ierr = VecGetOwnershipRange(field_vec, &rstart, &rend); CHKERRQ(ierr);
1100 loc = rend - rstart;
1101
1102 ierr = VecGetArray(field_vec, &locArray); CHKERRQ(ierr);
1103 ierr = PetscMemcpy(locArray, buffer + rstart, loc * sizeof(PetscScalar)); CHKERRQ(ierr);
1104 ierr = VecRestoreArray(field_vec, &locArray); CHKERRQ(ierr);
1105
1106 if (rank != 0) {
1107 ierr = PetscFree(buffer); CHKERRQ(ierr);
1108 }
1109 }
1110
1111 /* -------------------- tidy up ---------------------------------------- */
1112 if (rank == 0) {
1113 ierr = VecRestoreArrayRead(seq_vec, &seqArray); CHKERRQ(ierr);
1114 ierr = VecDestroy(&seq_vec); CHKERRQ(ierr);
1115 }
1116
1117 LOG_ALLOW(GLOBAL,LOG_INFO,
1118 "Loaded <%s> (parallel path)\n",filename);
1119
1120 PROFILE_FUNCTION_END;
1121 PetscFunctionReturn(0);
1122}
Here is the caller graph for this function:

◆ ReadSimulationFields()

PetscErrorCode ReadSimulationFields ( UserCtx *  user,
PetscInt  ti 
)

Internal helper implementation: ReadSimulationFields().

Reads binary field data for velocity, pressure, and other required vectors.

Local to this translation unit.

Definition at line 1129 of file io.c.

1130{
1131 PetscErrorCode ierr;
1132
1133 SimCtx *simCtx = user->simCtx;
1134 const char *source_path = NULL;
1135 const char *eulerian_ext = "dat";
1136
1137 if(simCtx->exec_mode == EXEC_MODE_POSTPROCESSOR){
1138 source_path = simCtx->pps->source_dir;
1139 if (simCtx->pps->eulerianExt[0] != '\0') {
1140 eulerian_ext = simCtx->pps->eulerianExt;
1141 }
1142 } else if(simCtx->exec_mode == EXEC_MODE_SOLVER){
1143 source_path = simCtx->restart_dir;
1144 } else{
1145 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Invalid execution mode for reading simulation fields.");
1146 }
1147
1148 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to read simulation fields.\n");
1149
1150 // Set the current I/O directory context
1151 ierr = PetscSNPrintf(simCtx->_io_context_buffer, sizeof(simCtx->_io_context_buffer), "%s/%s",source_path,simCtx->euler_subdir); CHKERRQ(ierr);
1152 simCtx->current_io_directory = simCtx->_io_context_buffer;
1153
1154 // Read Cartesian velocity field
1155 ierr = ReadFieldData(user, "ufield", user->Ucat, ti, eulerian_ext); CHKERRQ(ierr);
1156
1157 // Read contravariant velocity field
1158 ierr = ReadFieldData(user, "vfield", user->Ucont, ti, eulerian_ext); CHKERRQ(ierr);
1159
1160 // Read pressure field
1161 ierr = ReadFieldData(user, "pfield", user->P, ti, eulerian_ext); CHKERRQ(ierr);
1162
1163 // Read node state field (nvert)
1164 ierr = ReadFieldData(user, "nvfield", user->Nvert, ti, eulerian_ext); CHKERRQ(ierr);
1165
1166 LOG_ALLOW(GLOBAL,LOG_INFO,"Successfully read all mandatory fields. \n");
1167
1168 if(simCtx->np>0){
1169 // Read Particle Count field
1170 if(!user->ParticleCount){
1171 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "ParticleCount Vec is NULL but np>0");
1172 }
1173 ierr = ReadOptionalField(user, "ParticleCount", "Particle Count", user->ParticleCount, ti, eulerian_ext); CHKERRQ(ierr);
1174
1175 ierr = ReadOptionalField(user, "psifield", "Scalar Psi Field", user->Psi, ti, eulerian_ext); CHKERRQ(ierr);
1176 }
1177 else{
1178 LOG_ALLOW(GLOBAL, LOG_INFO, "No particles in simulation, skipping Particle fields reading.\n");
1179 }
1180 // Process LES fields if enabled
1181 if (simCtx->les) {
1182 ierr = ReadLESFields(user,ti); CHKERRQ(ierr);
1183 }
1184
1185 // Process RANS fields if enabled
1186 if (simCtx->rans) {
1187 ierr = ReadRANSFields(user,ti); CHKERRQ(ierr);
1188 }
1189
1190 // Process statistical fields if averaging is enabled
1191 if (simCtx->averaging) {
1192 ierr = ReadStatisticalFields(user,ti); CHKERRQ(ierr);
1193 }
1194
1195 simCtx->current_io_directory = NULL; // Clear the I/O context after reading
1196
1197 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished reading simulation fields.\n");
1198
1199 return 0;
1200}
ReadStatisticalFields
PetscErrorCode ReadStatisticalFields(UserCtx *user, PetscInt ti)
Implementation of ReadStatisticalFields().
Definition io.c:1208
ReadRANSFields
PetscErrorCode ReadRANSFields(UserCtx *user, PetscInt ti)
Internal helper implementation: ReadRANSFields().
Definition io.c:1257
ReadLESFields
PetscErrorCode ReadLESFields(UserCtx *user, PetscInt ti)
Internal helper implementation: ReadLESFields().
Definition io.c:1232
ReadOptionalField
static PetscErrorCode ReadOptionalField(UserCtx *user, const char *field_name, const char *field_label, Vec field_vec, PetscInt ti, const char *ext)
Internal helper implementation: ReadOptionalField().
Definition io.c:825
SimCtx::euler_subdir
char euler_subdir[PETSC_MAX_PATH_LEN]
Definition variables.h:707
SimCtx::rans
PetscInt rans
Definition variables.h:789
SimCtx::np
PetscInt np
Definition variables.h:796
SimCtx::averaging
PetscBool averaging
Definition variables.h:793
UserCtx::Ucont
Vec Ucont
Definition variables.h:904
PostProcessParams::eulerianExt
char eulerianExt[8]
Definition variables.h:619
PostProcessParams::source_dir
char source_dir[PETSC_MAX_PATH_LEN]
Definition variables.h:596
UserCtx::Ucat
Vec Ucat
Definition variables.h:904
UserCtx::ParticleCount
Vec ParticleCount
Definition variables.h:952
SimCtx::pps
PostProcessParams * pps
Definition variables.h:860
EXEC_MODE_SOLVER
@ EXEC_MODE_SOLVER
Definition variables.h:657
EXEC_MODE_POSTPROCESSOR
@ EXEC_MODE_POSTPROCESSOR
Definition variables.h:658
SimCtx::_io_context_buffer
char _io_context_buffer[PETSC_MAX_PATH_LEN]
Definition variables.h:710
SimCtx::les
PetscInt les
Definition variables.h:789
UserCtx::Nvert
Vec Nvert
Definition variables.h:904
SimCtx::exec_mode
ExecutionMode exec_mode
Definition variables.h:703
SimCtx::restart_dir
char restart_dir[PETSC_MAX_PATH_LEN]
Definition variables.h:705
UserCtx::Psi
Vec Psi
Definition variables.h:953
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadStatisticalFields()

PetscErrorCode ReadStatisticalFields ( UserCtx *  user,
PetscInt  ti 
)

Implementation of ReadStatisticalFields().

Reads statistical fields used for time-averaged simulations.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
ReadStatisticalFields()

Definition at line 1208 of file io.c.

1209{
1210 PetscErrorCode ierr;
1211 const char *eulerian_ext = "dat";
1212 if (user->simCtx->exec_mode == EXEC_MODE_POSTPROCESSOR && user->simCtx->pps->eulerianExt[0] != '\0') {
1213 eulerian_ext = user->simCtx->pps->eulerianExt;
1214 }
1215
1216 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to read statistical fields.\n");
1217
1218 ierr = ReadOptionalField(user, "su0", "Velocity Sum", user->Ucat_sum, ti, eulerian_ext); CHKERRQ(ierr);
1219 ierr = ReadOptionalField(user, "su1", "Velocity Cross Sum", user->Ucat_cross_sum, ti, eulerian_ext); CHKERRQ(ierr);
1220 ierr = ReadOptionalField(user, "su2", "Velocity Square Sum",user->Ucat_square_sum, ti, eulerian_ext); CHKERRQ(ierr);
1221 ierr = ReadOptionalField(user, "sp", "Pressure Sum", user->P_sum, ti, eulerian_ext); CHKERRQ(ierr);
1222
1223 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished reading statistical fields.\n");
1224
1225 return 0;
1226}
UserCtx::Ucat_square_sum
Vec Ucat_square_sum
Definition variables.h:938
UserCtx::P_sum
Vec P_sum
Definition variables.h:938
UserCtx::Ucat_sum
Vec Ucat_sum
Definition variables.h:938
UserCtx::Ucat_cross_sum
Vec Ucat_cross_sum
Definition variables.h:938
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadLESFields()

PetscErrorCode ReadLESFields ( UserCtx *  user,
PetscInt  ti 
)

Internal helper implementation: ReadLESFields().

Reads LES-related fields used in turbulence modeling.

Local to this translation unit.

Definition at line 1232 of file io.c.

1233{
1234 PetscErrorCode ierr;
1235 const char *eulerian_ext = "dat";
1236 if (user->simCtx->exec_mode == EXEC_MODE_POSTPROCESSOR && user->simCtx->pps->eulerianExt[0] != '\0') {
1237 eulerian_ext = user->simCtx->pps->eulerianExt;
1238 }
1239
1240 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to read LES fields.\n");
1241
1242 ierr = ReadOptionalField(user, "Nu_t", "Turbulent Viscosity", user->Nu_t, ti, eulerian_ext); CHKERRQ(ierr);
1243 ierr = ReadOptionalField(user, "cs", "Smagorinsky Constant (Cs)", user->CS, ti, eulerian_ext); CHKERRQ(ierr);
1244
1245 ierr = UpdateLocalGhosts(user, "CS"); CHKERRQ(ierr);
1246 ierr = UpdateLocalGhosts(user, "Nu_t"); CHKERRQ(ierr);
1247
1248 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished reading LES fields.\n");
1249
1250 return 0;
1251}
UpdateLocalGhosts
PetscErrorCode UpdateLocalGhosts(UserCtx *user, const char *fieldName)
Updates the local vector (including ghost points) from its corresponding global vector.
Definition setup.c:1755
UserCtx::Nu_t
Vec Nu_t
Definition variables.h:935
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadRANSFields()

PetscErrorCode ReadRANSFields ( UserCtx *  user,
PetscInt  ti 
)

Internal helper implementation: ReadRANSFields().

Reads RANS-related fields for turbulence modeling.

Local to this translation unit.

Definition at line 1257 of file io.c.

1258{
1259 PetscErrorCode ierr;
1260 const char *eulerian_ext = "dat";
1261 if (user->simCtx->exec_mode == EXEC_MODE_POSTPROCESSOR && user->simCtx->pps->eulerianExt[0] != '\0') {
1262 eulerian_ext = user->simCtx->pps->eulerianExt;
1263 }
1264
1265 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to read RANS fields.\n");
1266
1267 ierr = ReadOptionalField(user, "kfield", "K-Omega RANS", user->K_Omega, ti, eulerian_ext); CHKERRQ(ierr);
1268 ierr = ReadOptionalField(user, "Nu_t", "Turbulent Viscosity", user->Nu_t, ti, eulerian_ext); CHKERRQ(ierr);
1269
1270 VecCopy(user->K_Omega, user->K_Omega_o);
1271
1272 ierr = UpdateLocalGhosts(user, "K_Omega"); CHKERRQ(ierr);
1273 ierr = UpdateLocalGhosts(user, "Nu_t"); CHKERRQ(ierr);
1274 ierr = UpdateLocalGhosts(user, "K_Omega_o"); CHKERRQ(ierr);
1275
1276 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished reading RANS fields.\n");
1277
1278 return 0;
1279}
UserCtx::K_Omega_o
Vec K_Omega_o
Definition variables.h:935
UserCtx::K_Omega
Vec K_Omega
Definition variables.h:935
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadSwarmField()

PetscErrorCode ReadSwarmField ( UserCtx *  user,
const char *  field_name,
PetscInt  ti,
const char *  ext 
)

Internal helper implementation: ReadSwarmField().

Reads data from a file into a specified field of a PETSc DMSwarm.

Local to this translation unit.

Definition at line 1286 of file io.c.

1287{
1288 PetscErrorCode ierr;
1289 DM swarm;
1290 Vec fieldVec;
1291
1292 PetscFunctionBegin;
1293
1294 swarm = user->swarm;
1295
1296 LOG_ALLOW(GLOBAL,LOG_DEBUG," ReadSwarmField Begins \n");
1297
1298 /* 2) Create a global vector that references the specified Swarm field. */
1299 ierr = DMSwarmCreateGlobalVectorFromField(swarm, field_name, &fieldVec);CHKERRQ(ierr);
1300
1301 LOG_ALLOW(GLOBAL,LOG_DEBUG," Vector created from Field \n");
1302
1303 /* 3) Use the ReadFieldData() function to read data into fieldVec. */
1304 ierr = ReadFieldData(user, field_name, fieldVec, ti, ext);CHKERRQ(ierr);
1305
1306 /* 4) Destroy the global vector reference. */
1307 ierr = DMSwarmDestroyGlobalVectorFromField(swarm, field_name, &fieldVec);CHKERRQ(ierr);
1308
1309 PetscFunctionReturn(0);
1310}
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadSwarmIntField()

PetscErrorCode ReadSwarmIntField ( UserCtx *  user,
const char *  field_name,
PetscInt  ti,
const char *  ext 
)

Internal helper implementation: ReadSwarmIntField().

Reads integer swarm data by using ReadFieldData and casting the result.

Local to this translation unit.

Definition at line 1316 of file io.c.

1317{
1318 PetscErrorCode ierr;
1319 DM swarm = user->swarm;
1320 Vec temp_vec;
1321 PetscInt nlocal, nglobal, bs, i;
1322 const PetscScalar *scalar_array; // Read-only pointer from the temp Vec
1323 void *field_array_void;
1324
1325
1326 PetscFunctionBeginUser;
1327
1328 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Reading '%s' via temporary Vec.\n", field_name);
1329
1330 // Get the properties of the swarm field to determine the expected layout
1331 ierr = DMSwarmGetLocalSize(swarm, &nlocal); CHKERRQ(ierr);
1332 ierr = DMSwarmGetSize(swarm, &nglobal); CHKERRQ(ierr);
1333 // We get the block size but not the data pointer yet
1334 ierr = DMSwarmGetField(swarm, field_name, &bs, NULL, NULL); CHKERRQ(ierr);
1335 ierr = DMSwarmRestoreField(swarm, field_name, &bs, NULL, NULL); CHKERRQ(ierr);
1336
1337 // Create a temporary Vec with the CORRECT layout to receive the data
1338 ierr = VecCreate(PETSC_COMM_WORLD, &temp_vec); CHKERRQ(ierr);
1339 ierr = VecSetType(temp_vec, VECMPI); CHKERRQ(ierr);
1340 ierr = VecSetSizes(temp_vec, nlocal * bs, nglobal * bs); CHKERRQ(ierr);
1341 ierr = VecSetBlockSize(temp_vec, bs); CHKERRQ(ierr);
1342 ierr = VecSetUp(temp_vec); CHKERRQ(ierr);
1343
1344 // Call your existing reader to populate the temporary Vec
1345 ierr = ReadFieldData(user, field_name, temp_vec, ti, ext); CHKERRQ(ierr);
1346
1347 // Get local pointers
1348 ierr = VecGetArrayRead(temp_vec, &scalar_array); CHKERRQ(ierr);
1349 ierr = DMSwarmGetField(swarm, field_name, NULL, NULL, &field_array_void); CHKERRQ(ierr);
1350
1351 // Perform the cast back, using the correct loop size (nlocal * bs)
1352 if (strcmp(field_name, "DMSwarm_pid") == 0) {
1353 PetscInt64 *int64_array = (PetscInt64 *)field_array_void;
1354 for (i = 0; i < nlocal * bs; i++) {
1355 int64_array[i] = (PetscInt64)scalar_array[i];
1356 }
1357 } else {
1358 PetscInt *int_array = (PetscInt *)field_array_void;
1359 for (i = 0; i < nlocal * bs; i++) {
1360 int_array[i] = (PetscInt)scalar_array[i];
1361 }
1362 }
1363
1364 // Restore access
1365 ierr = DMSwarmRestoreField(swarm, field_name, NULL, NULL, &field_array_void); CHKERRQ(ierr);
1366 ierr = VecRestoreArrayRead(temp_vec, &scalar_array); CHKERRQ(ierr);
1367
1368 // 6. Clean up
1369 ierr = VecDestroy(&temp_vec); CHKERRQ(ierr);
1370
1371 PetscFunctionReturn(0);
1372}
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReadAllSwarmFields()

PetscErrorCode ReadAllSwarmFields ( UserCtx *  user,
PetscInt  ti 
)

Internal helper implementation: ReadAllSwarmFields().

Reads multiple fields (positions, velocity, CellID, and weight) into a DMSwarm.

Local to this translation unit.

Definition at line 1378 of file io.c.

1379{
1380 PetscErrorCode ierr;
1381 PetscInt nGlobal;
1382 SimCtx *simCtx = user->simCtx;
1383 const char *source_path = NULL;
1384 const char *particle_ext = "dat";
1385
1386 PetscFunctionBeginUser;
1387 ierr = DMSwarmGetSize(user->swarm, &nGlobal); CHKERRQ(ierr);
1388 LOG_ALLOW(GLOBAL, LOG_INFO, "Reading DMSwarm fields for timestep %d (swarm size is %d).\n", ti, nGlobal);
1389
1390 if (nGlobal == 0) {
1391 LOG_ALLOW(GLOBAL, LOG_INFO, "Swarm is empty for timestep %d. Nothing to read.\n", ti);
1392 PetscFunctionReturn(0);
1393 }
1394
1395 // First, determine the top-level source directory based on the execution mode.
1396 if (simCtx->exec_mode == EXEC_MODE_SOLVER) {
1397 source_path = simCtx->restart_dir;
1398 } else if (simCtx->exec_mode == EXEC_MODE_POSTPROCESSOR) {
1399 source_path = simCtx->pps->source_dir;
1400 if (simCtx->pps->particleExt[0] != '\0') {
1401 particle_ext = simCtx->pps->particleExt;
1402 }
1403 } else {
1404 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "Invalid execution mode for reading simulation fields.");
1405 }
1406
1407 // Set the current I/O directory context
1408 ierr = PetscSNPrintf(simCtx->_io_context_buffer, sizeof(simCtx->_io_context_buffer),
1409 "%s/%s", source_path, simCtx->particle_subdir); CHKERRQ(ierr);
1410
1411 simCtx->current_io_directory = simCtx->_io_context_buffer;
1412
1413 /* 1) Read positions (REQUIRED) */
1414 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Reading mandatory position field...\n");
1415 ierr = ReadSwarmField(user, "position", ti, particle_ext);
1416 if (ierr) {
1417 SETERRQ(PETSC_COMM_WORLD, ierr, "Failed to read MANDATORY 'position' field for step %d. Cannot continue.", ti);
1418 }
1419 LOG_ALLOW(GLOBAL, LOG_INFO, "Successfully read mandatory position field for step %d.\n", ti);
1420
1421 /* 2) Read all OPTIONAL fields using the helper function. */
1422 /* The helper will print a warning and continue if a file is not found. */
1423 ierr = ReadOptionalSwarmField(user, "velocity", "Velocity", ti, particle_ext); CHKERRQ(ierr);
1424 ierr = ReadOptionalSwarmField(user, "DMSwarm_pid", "Particle ID", ti, particle_ext); CHKERRQ(ierr);
1425 ierr = ReadOptionalSwarmField(user, "DMSwarm_CellID", "Cell ID", ti, particle_ext); CHKERRQ(ierr);
1426 ierr = ReadOptionalSwarmField(user, "weight", "Particle Weight", ti, particle_ext); CHKERRQ(ierr);
1427 ierr = ReadOptionalSwarmField(user, "Psi", "Scalar Psi", ti, particle_ext); CHKERRQ(ierr);
1428 ierr = ReadOptionalSwarmField(user, "DMSwarm_location_status", "Migration Status", ti, particle_ext); CHKERRQ(ierr);
1429
1430 simCtx->current_io_directory = NULL; // Clear the I/O context after reading
1431
1432 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished reading DMSwarm fields for timestep %d.\n", ti);
1433 PetscFunctionReturn(0);
1434}
ReadOptionalSwarmField
static PetscErrorCode ReadOptionalSwarmField(UserCtx *user, const char *field_name, const char *field_label, PetscInt ti, const char *ext)
Internal helper implementation: ReadOptionalSwarmField().
Definition io.c:859
SimCtx::particle_subdir
char particle_subdir[PETSC_MAX_PATH_LEN]
Definition variables.h:708
PostProcessParams::particleExt
char particleExt[8]
Definition variables.h:620
Here is the call graph for this function:
Here is the caller graph for this function:

◆ WriteFieldData()

PetscErrorCode WriteFieldData ( UserCtx *  user,
const char *  field_name,
Vec  field_vec,
PetscInt  ti,
const char *  ext 
)

Internal helper implementation: WriteFieldData().

Writes data from a specific PETSc vector to a file.

Local to this translation unit.

Definition at line 1443 of file io.c.

1448{
1449 PetscErrorCode ierr;
1450 MPI_Comm comm;
1451 PetscMPIInt rank, size;
1452
1453 const PetscInt placeholder_int = 0; /* keep legacy name */
1454 char filename[PETSC_MAX_PATH_LEN];
1455 SimCtx *simCtx=user->simCtx;
1456
1457 PetscFunctionBeginUser;
1458 PROFILE_FUNCTION_BEGIN;
1459
1460 if(!simCtx->current_io_directory){
1461 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "I/O context directory was not set before calling WriteFieldData().");
1462 }
1463
1464 /* ------------------------------------------------------------ */
1465 /* Basic communicator information */
1466 /* ------------------------------------------------------------ */
1467 ierr = PetscObjectGetComm((PetscObject)field_vec,&comm);CHKERRQ(ierr);
1468 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr);
1469 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr);
1470
1471 ierr = PetscSNPrintf(filename,sizeof(filename),
1472 "%s/%s%05" PetscInt_FMT "_%d.%s",
1473 simCtx->current_io_directory,field_name,ti,placeholder_int,ext);CHKERRQ(ierr);
1474
1475 LOG_ALLOW(GLOBAL,LOG_DEBUG,
1476 " Preparing to write <%s> on rank %d/%d\n",
1477 filename,rank,size);
1478
1479 /* ------------------------------------------------------------ */
1480 /* 1. Serial path */
1481 /* ------------------------------------------------------------ */
1482 if (size == 1) {
1483 PetscViewer viewer;
1484
1485 ierr = PetscViewerBinaryOpen(comm,filename,
1486 FILE_MODE_WRITE,&viewer);CHKERRQ(ierr);
1487 ierr = VecView(field_vec,viewer);CHKERRQ(ierr);
1488 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
1489
1490 LOG_ALLOW(GLOBAL,LOG_INFO,
1491 " Wrote <%s> (serial path)\n",filename);
1492 PetscFunctionReturn(0);
1493 }
1494
1495 /* ------------------------------------------------------------ */
1496 /* 2. Parallel path */
1497 /* ------------------------------------------------------------ */
1498 VecScatter scatter;
1499 Vec seq_vec=NULL; /* created by PETSc, lives only on rank 0 */
1500 DM dm = NULL;
1501 const char *dmtype = NULL;
1502 Vec nat = NULL; /* Natural-ordered vector for DMDA */
1503
1504 /* 2.1 Check if vector has DMDA and convert to natural ordering if needed */
1505 ierr = VecGetDM(field_vec, &dm); CHKERRQ(ierr);
1506 if (dm) { ierr = DMGetType(dm, &dmtype); CHKERRQ(ierr); }
1507
1508 if (dmtype && !strcmp(dmtype, DMDA)) {
1509 /* DMDA path: convert to natural ordering first */
1510 ierr = DMDACreateNaturalVector(dm, &nat); CHKERRQ(ierr);
1511 ierr = DMDAGlobalToNaturalBegin(dm, field_vec, INSERT_VALUES, nat); CHKERRQ(ierr);
1512 ierr = DMDAGlobalToNaturalEnd(dm, field_vec, INSERT_VALUES, nat); CHKERRQ(ierr);
1513
1514 /* Gather the natural-ordered vector */
1515 ierr = VecScatterCreateToZero(nat, &scatter, &seq_vec); CHKERRQ(ierr);
1516 } else {
1517 /* Non-DMDA path: direct gather in global ordering */
1518 ierr = VecScatterCreateToZero(field_vec, &scatter, &seq_vec); CHKERRQ(ierr);
1519 }
1520
1521 /* 2.2 Gather distributed → sequential (on rank 0) */
1522 ierr = VecScatterBegin(scatter, (nat ? nat : field_vec), seq_vec,
1523 INSERT_VALUES, SCATTER_FORWARD); CHKERRQ(ierr);
1524 ierr = VecScatterEnd(scatter, (nat ? nat : field_vec), seq_vec,
1525 INSERT_VALUES, SCATTER_FORWARD); CHKERRQ(ierr);
1526
1527 /* 2.3 Rank 0 writes the file */
1528 if (rank == 0) {
1529 PetscViewer viewer;
1530
1531 /* (optional) value diagnostics */
1532 PetscReal vmin,vmax;
1533 ierr = VecMin(seq_vec,NULL,&vmin);CHKERRQ(ierr);
1534 ierr = VecMax(seq_vec,NULL,&vmax);CHKERRQ(ierr);
1535 LOG_ALLOW(GLOBAL,LOG_DEBUG,
1536 " <%s> range = [%.4e … %.4e]\n",
1537 field_name,(double)vmin,(double)vmax);
1538
1539 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,filename,
1540 FILE_MODE_WRITE,&viewer);CHKERRQ(ierr);
1541 ierr = VecView(seq_vec,viewer);CHKERRQ(ierr);
1542 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
1543
1544 LOG_ALLOW(GLOBAL,LOG_INFO,
1545 " Wrote <%s> (parallel path)\n",filename);
1546 }
1547
1548 /* 2.4 Cleanup */
1549 ierr = VecScatterDestroy(&scatter);CHKERRQ(ierr);
1550 ierr = VecDestroy(&seq_vec);CHKERRQ(ierr);
1551 if (nat) { ierr = VecDestroy(&nat); CHKERRQ(ierr); }
1552
1553 PROFILE_FUNCTION_END;
1554 PetscFunctionReturn(0);
1555}
Here is the caller graph for this function:

◆ WriteSimulationFields()

PetscErrorCode WriteSimulationFields ( UserCtx *  user)

Implementation of WriteSimulationFields().

Writes simulation fields to files.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
WriteSimulationFields()

Definition at line 1563 of file io.c.

1564{
1565 PetscErrorCode ierr;
1566
1567 SimCtx *simCtx = user->simCtx;
1568
1569 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to write simulation fields.\n");
1570
1571 // Set the current IO directory
1572 ierr = PetscSNPrintf(simCtx->_io_context_buffer, sizeof(simCtx->_io_context_buffer),
1573 "%s/%s", simCtx->output_dir, simCtx->euler_subdir);CHKERRQ(ierr);
1574 simCtx->current_io_directory = simCtx->_io_context_buffer;
1575
1576 // Write contravariant velocity field
1577 ierr = WriteFieldData(user, "vfield", user->Ucont, simCtx->step, "dat"); CHKERRQ(ierr);
1578
1579 // Write Cartesian velocity field
1580 ierr = WriteFieldData(user, "ufield", user->Ucat, simCtx->step, "dat"); CHKERRQ(ierr);
1581
1582 // Write pressure field
1583 ierr = WriteFieldData(user, "pfield", user->P, simCtx->step, "dat"); CHKERRQ(ierr);
1584
1585 // Write node state field (nvert)
1586 ierr = WriteFieldData(user, "nvfield", user->Nvert, simCtx->step, "dat"); CHKERRQ(ierr);
1587
1588 // Write ParticleCountPerCell if enabled.
1589 if(simCtx->np>0){
1590 ierr = WriteFieldData(user, "ParticleCount",user->ParticleCount,simCtx->step,"dat"); CHKERRQ(ierr);
1591 ierr = WriteFieldData(user, "psifield", user->Psi, simCtx->step, "dat"); CHKERRQ(ierr);
1592 }
1593
1594 // Write LES fields if enabled
1595 if (simCtx->les) {
1596 ierr = WriteLESFields(user); CHKERRQ(ierr);
1597 }
1598
1599 // Write RANS fields if enabled
1600 if (simCtx->rans) {
1601 ierr = WriteRANSFields(user); CHKERRQ(ierr);
1602 }
1603
1604 // Write statistical fields if averaging is enabled
1605 if (simCtx->averaging) {
1606 ierr = WriteStatisticalFields(user); CHKERRQ(ierr);
1607 }
1608
1609 simCtx->current_io_directory = NULL;
1610
1611 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished writing simulation fields.\n");
1612
1613 return 0;
1614}
WriteStatisticalFields
PetscErrorCode WriteStatisticalFields(UserCtx *user)
Implementation of WriteStatisticalFields().
Definition io.c:1622
WriteFieldData
PetscErrorCode WriteFieldData(UserCtx *user, const char *field_name, Vec field_vec, PetscInt ti, const char *ext)
Internal helper implementation: WriteFieldData().
Definition io.c:1443
WriteRANSFields
PetscErrorCode WriteRANSFields(UserCtx *user)
Implementation of WriteRANSFields().
Definition io.c:1674
WriteLESFields
PetscErrorCode WriteLESFields(UserCtx *user)
Internal helper implementation: WriteLESFields().
Definition io.c:1644
SimCtx::output_dir
char output_dir[PETSC_MAX_PATH_LEN]
Definition variables.h:706
SimCtx::step
PetscInt step
Definition variables.h:692
Here is the call graph for this function:
Here is the caller graph for this function:

◆ WriteStatisticalFields()

PetscErrorCode WriteStatisticalFields ( UserCtx *  user)

Implementation of WriteStatisticalFields().

Writes statistical fields for averaging purposes.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
WriteStatisticalFields()

Definition at line 1622 of file io.c.

1623{
1624 PetscErrorCode ierr;
1625
1626 SimCtx *simCtx = user->simCtx;
1627
1628 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to write statistical fields.\n");
1629
1630 ierr = WriteFieldData(user, "su0", user->Ucat_sum, simCtx->step, "dat"); CHKERRQ(ierr);
1631 ierr = WriteFieldData(user, "su1", user->Ucat_cross_sum, simCtx->step, "dat"); CHKERRQ(ierr);
1632 ierr = WriteFieldData(user, "su2", user->Ucat_square_sum, simCtx->step, "dat"); CHKERRQ(ierr);
1633 ierr = WriteFieldData(user, "sp", user->P_sum, simCtx->step, "dat"); CHKERRQ(ierr);
1634
1635 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished writing statistical fields.\n");
1636
1637 return 0;
1638}
Here is the call graph for this function:
Here is the caller graph for this function:

◆ WriteLESFields()

PetscErrorCode WriteLESFields ( UserCtx *  user)

Internal helper implementation: WriteLESFields().

Writes LES-related fields.

Local to this translation unit.

Definition at line 1644 of file io.c.

1645{
1646 PetscErrorCode ierr;
1647
1648 SimCtx *simCtx = user->simCtx;
1649
1650 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to write LES fields.\n");
1651
1652
1653 DMLocalToGlobalBegin(user->da, user->lCs, INSERT_VALUES, user->CS);
1654 DMLocalToGlobalEnd(user->da, user->lCs, INSERT_VALUES, user->CS);
1655
1656 DMLocalToGlobalBegin(user->da, user->lNu_t, INSERT_VALUES, user->Nu_t);
1657 DMLocalToGlobalEnd(user->da, user->lNu_t, INSERT_VALUES, user->Nu_t);
1658
1659 ierr = WriteFieldData(user, "Nu_t", user->Nu_t, simCtx->step, "dat"); CHKERRQ(ierr);
1660 ierr = WriteFieldData(user, "cs", user->CS, simCtx->step, "dat"); CHKERRQ(ierr);
1661
1662
1663 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished writing LES fields.\n");
1664
1665 return 0;
1666}
UserCtx::lCs
Vec lCs
Definition variables.h:935
UserCtx::lNu_t
Vec lNu_t
Definition variables.h:935
Here is the call graph for this function:
Here is the caller graph for this function:

◆ WriteRANSFields()

PetscErrorCode WriteRANSFields ( UserCtx *  user)

Implementation of WriteRANSFields().

Writes RANS-related fields.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
WriteRANSFields()

Definition at line 1674 of file io.c.

1675{
1676 PetscErrorCode ierr;
1677
1678 SimCtx *simCtx = user->simCtx;
1679
1680 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to write RANS fields.\n");
1681
1682 ierr = WriteFieldData(user, "kfield", user->K_Omega, simCtx->step, "dat"); CHKERRQ(ierr);
1683
1684 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished writing RANS fields.\n");
1685
1686 return 0;
1687}
Here is the call graph for this function:
Here is the caller graph for this function:

◆ WriteSwarmField()

PetscErrorCode WriteSwarmField ( UserCtx *  user,
const char *  field_name,
PetscInt  ti,
const char *  ext 
)

Implementation of WriteSwarmField().

Writes data from a specific field in a PETSc Swarm to a file.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
WriteSwarmField()

Definition at line 1695 of file io.c.

1696{
1697 PetscErrorCode ierr;
1698 Vec fieldVec;
1699 DM swarm;
1700
1701 PetscFunctionBeginUser; /* PETSc macro indicating start of function */
1702
1703 /*
1704 * 1) Retrieve the PetscSwarm from the user context.
1705 * Ensure user->swarm is initialized and not NULL.
1706 */
1707 swarm = user->swarm;
1708
1709 /*
1710 * 2) Create a global vector from the specified swarm field.
1711 * This function is available in PETSc 3.14.4.
1712 * It provides a read/write "view" of the swarm field as a global Vec.
1713 */
1714 LOG_ALLOW(GLOBAL, LOG_DEBUG,
1715 "Attempting to create global vector from field: %s\n",
1716 field_name);
1717 ierr = DMSwarmCreateGlobalVectorFromField(swarm, field_name, &fieldVec);CHKERRQ(ierr);
1718
1719 /*
1720 * 3) Use your existing WriteFieldData() to write the global vector to a file.
1721 * The field name, time index, and extension are passed along for naming.
1722 */
1723 LOG_ALLOW(GLOBAL, LOG_DEBUG,
1724 "Calling WriteFieldData for field: %s\n",
1725 field_name);
1726 ierr = WriteFieldData(user, field_name, fieldVec, ti, ext);CHKERRQ(ierr);
1727
1728 /*
1729 * 4) Destroy the global vector once the data is successfully written.
1730 * This step is crucial for avoiding memory leaks.
1731 * DMSwarmDestroyGlobalVectorFromField() is also available in PETSc 3.14.4.
1732 */
1733 LOG_ALLOW(GLOBAL, LOG_DEBUG,
1734 "Destroying the global vector for field: %s\n",
1735 field_name);
1736 ierr = DMSwarmDestroyGlobalVectorFromField(swarm, field_name, &fieldVec);CHKERRQ(ierr);
1737
1738 /* Log and return success. */
1739 LOG_ALLOW(GLOBAL, LOG_INFO,
1740 "Successfully wrote swarm data for field: %s\n",
1741 field_name);
1742
1743 PetscFunctionReturn(0); /* PETSc macro indicating end of function */
1744}
Here is the call graph for this function:
Here is the caller graph for this function:

◆ WriteSwarmIntField()

PetscErrorCode WriteSwarmIntField ( UserCtx *  user,
const char *  field_name,
PetscInt  ti,
const char *  ext 
)

Internal helper implementation: WriteSwarmIntField().

Writes integer data from a specific PETSc Swarm field to a file.

Local to this translation unit.

Definition at line 1750 of file io.c.

1751{
1752 PetscErrorCode ierr;
1753 DM swarm = user->swarm;
1754 Vec temp_vec; // Temporary Vec to hold casted data
1755 PetscInt nlocal, nglobal,bs,i;
1756 void *field_array_void;
1757 PetscScalar *scalar_array; // Pointer to the temporary Vec's scalar data
1758
1759 PetscFunctionBeginUser;
1760
1761 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Casting '%s' to Vec for writing.\n", field_name);
1762
1763 // Get the swarm field properties
1764 ierr = DMSwarmGetLocalSize(swarm, &nlocal); CHKERRQ(ierr);
1765 ierr = DMSwarmGetSize(swarm, &nglobal); CHKERRQ(ierr);
1766 ierr = DMSwarmGetField(swarm, field_name, &bs, NULL, &field_array_void); CHKERRQ(ierr);
1767
1768 // Create Temporary parallel Vec wit the CORRECT layout
1769 ierr = VecCreate(PETSC_COMM_WORLD, &temp_vec); CHKERRQ(ierr);
1770 ierr = VecSetType(temp_vec, VECMPI); CHKERRQ(ierr);
1771 ierr = VecSetSizes(temp_vec, nlocal*bs, nglobal*bs); CHKERRQ(ierr);
1772 ierr = VecSetUp(temp_vec); CHKERRQ(ierr);
1773
1774 // Defining Vector field to mandatory field 'position'
1775 DMSwarmVectorDefineField(swarm,"position");
1776
1777 ierr = VecGetArray(temp_vec, &scalar_array); CHKERRQ(ierr);
1778
1779 if(strcasecmp(field_name,"DMSwarm_pid") == 0){
1780 PetscInt64 *int64_array = (PetscInt64 *)field_array_void;
1781 // Perform the cast from PetscInt64 to PetscScalar
1782 for (i = 0; i < nlocal*bs; i++) {
1783 scalar_array[i] = (PetscScalar)int64_array[i];
1784 }
1785 }else{
1786 PetscInt *int_array = (PetscInt *)field_array_void;
1787 //Perform the cast from PetscInt to PetscScalar
1788 for (i = 0; i < nlocal*bs; i++) {
1789 scalar_array[i] = (PetscScalar)int_array[i];
1790 }
1791 }
1792
1793 // Restore access to both arrays
1794 ierr = VecRestoreArray(temp_vec, &scalar_array); CHKERRQ(ierr);
1795 ierr = DMSwarmRestoreField(swarm, field_name, &bs, NULL, &field_array_void); CHKERRQ(ierr);
1796
1797 // Call your existing writer with the temporary, populated Vec
1798 ierr = WriteFieldData(user, field_name, temp_vec, ti, ext); CHKERRQ(ierr);
1799
1800 // Clean up
1801 ierr = VecDestroy(&temp_vec); CHKERRQ(ierr);
1802
1803 PetscFunctionReturn(0);
1804}
Here is the call graph for this function:
Here is the caller graph for this function:

◆ WriteAllSwarmFields()

PetscErrorCode WriteAllSwarmFields ( UserCtx *  user)

Internal helper implementation: WriteAllSwarmFields().

Writes a predefined set of PETSc Swarm fields to files.

Local to this translation unit.

Definition at line 1810 of file io.c.

1811{
1812 PetscErrorCode ierr;
1813 SimCtx *simCtx = user->simCtx;
1814
1815 PetscFunctionBeginUser;
1816
1817 // If no swarm is configured or there are no particles, do nothing and return.
1818 if (!user->swarm || simCtx->np <= 0) {
1819 PetscFunctionReturn(0);
1820 }
1821
1822 LOG_ALLOW(GLOBAL, LOG_INFO, "Starting to write swarm fields.\n");
1823
1824 // Ensure the current IO directory is set
1825 ierr = PetscSNPrintf(simCtx->_io_context_buffer, sizeof(simCtx->_io_context_buffer),
1826 "%s/%s", simCtx->output_dir, simCtx->particle_subdir);CHKERRQ(ierr);
1827 simCtx->current_io_directory = simCtx->_io_context_buffer;
1828
1829 // Write particle position field
1830 ierr = WriteSwarmField(user, "position", simCtx->step, "dat"); CHKERRQ(ierr);
1831
1832 // Write particle velocity field
1833 ierr = WriteSwarmField(user, "velocity", simCtx->step, "dat"); CHKERRQ(ierr);
1834
1835 // Write particle weight field
1836 ierr = WriteSwarmField(user, "weight", simCtx->step, "dat"); CHKERRQ(ierr);
1837
1838 // Write custom particle field "Psi"
1839 ierr = WriteSwarmField(user, "Psi", simCtx->step, "dat"); CHKERRQ(ierr);
1840
1841 // Integer fields require special handling
1842
1843 // Write the background mesh cell ID for each particle
1844 ierr = WriteSwarmIntField(user, "DMSwarm_CellID", simCtx->step, "dat"); CHKERRQ(ierr);
1845
1846 // Write the particle location status (e.g., inside or outside the domain)
1847 ierr = WriteSwarmIntField(user, "DMSwarm_location_status", simCtx->step, "dat"); CHKERRQ(ierr);
1848
1849 // Write the unique particle ID
1850 ierr = WriteSwarmIntField(user, "DMSwarm_pid", simCtx->step, "dat"); CHKERRQ(ierr);
1851
1852 simCtx->current_io_directory = NULL;
1853
1854 LOG_ALLOW(GLOBAL, LOG_INFO, "Finished writing swarm fields.\n");
1855
1856 PetscFunctionReturn(0);
1857}
WriteSwarmField
PetscErrorCode WriteSwarmField(UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
Implementation of WriteSwarmField().
Definition io.c:1695
WriteSwarmIntField
PetscErrorCode WriteSwarmIntField(UserCtx *user, const char *field_name, PetscInt ti, const char *ext)
Internal helper implementation: WriteSwarmIntField().
Definition io.c:1750
Here is the call graph for this function:
Here is the caller graph for this function:

◆ VecToArrayOnRank0()

PetscErrorCode VecToArrayOnRank0 ( Vec  inVec,
PetscInt *  N,
double **  arrayOut 
)

Internal helper implementation: VecToArrayOnRank0().

Gathers the contents of a distributed PETSc Vec into a single array on rank 0.

Local to this translation unit.

Definition at line 1863 of file io.c.

1864{
1865 PetscErrorCode ierr;
1866 MPI_Comm comm;
1867 PetscMPIInt rank;
1868 PetscInt globalSize;
1869 DM dm = NULL;
1870 const char *dmtype = NULL;
1871
1872 /* For DMDA path */
1873 Vec nat = NULL, seqNat = NULL;
1874 VecScatter scatNat = NULL;
1875 const PetscScalar *nar = NULL;
1876 PetscScalar *buf = NULL;
1877
1878 /* For generic (no DM) path */
1879 Vec seq = NULL;
1880 VecScatter scat = NULL;
1881 const PetscScalar *sar = NULL;
1882
1883 PetscFunctionBeginUser;
1884
1885 ierr = PetscObjectGetComm((PetscObject)inVec, &comm); CHKERRQ(ierr);
1886 ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr);
1887 ierr = VecGetSize(inVec, &globalSize); CHKERRQ(ierr);
1888 *N = globalSize;
1889 *arrayOut = NULL;
1890
1891 ierr = VecGetDM(inVec, &dm); CHKERRQ(ierr);
1892 if (dm) { ierr = DMGetType(dm, &dmtype); CHKERRQ(ierr); }
1893
1894 if (dmtype && !strcmp(dmtype, DMDA)) {
1895 /* --- DMDA path: go to NATURAL ordering, then gather to rank 0 --- */
1896 ierr = DMDACreateNaturalVector(dm, &nat); CHKERRQ(ierr);
1897 ierr = DMDAGlobalToNaturalBegin(dm, inVec, INSERT_VALUES, nat); CHKERRQ(ierr);
1898 ierr = DMDAGlobalToNaturalEnd (dm, inVec, INSERT_VALUES, nat); CHKERRQ(ierr);
1899
1900 ierr = VecScatterCreateToZero(nat, &scatNat, &seqNat); CHKERRQ(ierr);
1901 ierr = VecScatterBegin(scatNat, nat, seqNat, INSERT_VALUES, SCATTER_FORWARD); CHKERRQ(ierr);
1902 ierr = VecScatterEnd (scatNat, nat, seqNat, INSERT_VALUES, SCATTER_FORWARD); CHKERRQ(ierr);
1903
1904 if (rank == 0) {
1905 PetscInt nseq;
1906 ierr = VecGetLocalSize(seqNat, &nseq); CHKERRQ(ierr);
1907 ierr = VecGetArrayRead(seqNat, &nar); CHKERRQ(ierr);
1908
1909 ierr = PetscMalloc1(nseq, &buf); CHKERRQ(ierr);
1910 ierr = PetscMemcpy(buf, nar, (size_t)nseq * sizeof(PetscScalar)); CHKERRQ(ierr);
1911
1912 ierr = VecRestoreArrayRead(seqNat, &nar); CHKERRQ(ierr);
1913 *arrayOut = (double*)buf; /* hand back as double* to match the helper signature */
1914 }
1915
1916 ierr = VecScatterDestroy(&scatNat); CHKERRQ(ierr);
1917 ierr = VecDestroy(&seqNat); CHKERRQ(ierr);
1918 ierr = VecDestroy(&nat); CHKERRQ(ierr);
1919 } else {
1920 /* --- No DM attached: plain gather in Vec’s global (parallel) layout order --- */
1921 ierr = VecScatterCreateToZero(inVec, &scat, &seq); CHKERRQ(ierr);
1922 ierr = VecScatterBegin(scat, inVec, seq, INSERT_VALUES, SCATTER_FORWARD); CHKERRQ(ierr);
1923 ierr = VecScatterEnd (scat, inVec, seq, INSERT_VALUES, SCATTER_FORWARD); CHKERRQ(ierr);
1924
1925 if (rank == 0) {
1926 PetscInt nseq;
1927 ierr = VecGetLocalSize(seq, &nseq); CHKERRQ(ierr);
1928 ierr = VecGetArrayRead(seq, &sar); CHKERRQ(ierr);
1929
1930 ierr = PetscMalloc1(nseq, &buf); CHKERRQ(ierr);
1931 ierr = PetscMemcpy(buf, sar, (size_t)nseq * sizeof(PetscScalar)); CHKERRQ(ierr);
1932
1933 ierr = VecRestoreArrayRead(seq, &sar); CHKERRQ(ierr);
1934 *arrayOut = (double*)buf;
1935 }
1936
1937 ierr = VecScatterDestroy(&scat); CHKERRQ(ierr);
1938 ierr = VecDestroy(&seq); CHKERRQ(ierr);
1939 }
1940
1941 PetscFunctionReturn(0);
1942}
Here is the caller graph for this function:

◆ SwarmFieldToArrayOnRank0()

PetscErrorCode SwarmFieldToArrayOnRank0 ( DM  swarm,
const char *  field_name,
PetscInt *  n_total_particles,
PetscInt *  n_components,
void **  gathered_array 
)

Internal helper implementation: SwarmFieldToArrayOnRank0().

Gathers any DMSwarm field from all ranks to a single, contiguous array on rank 0.

Local to this translation unit.

Definition at line 1948 of file io.c.

1949{
1950 PetscErrorCode ierr;
1951 PetscMPIInt rank, size;
1952 PetscInt nlocal, nglobal, bs;
1953 void *local_array_void;
1954 size_t element_size = 0;
1955
1956 PetscFunctionBeginUser;
1957
1958 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
1959 ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size); CHKERRQ(ierr);
1960
1961 // All ranks get swarm properties to determine send/receive counts
1962 ierr = DMSwarmGetLocalSize(swarm, &nlocal); CHKERRQ(ierr);
1963 ierr = DMSwarmGetSize(swarm, &nglobal); CHKERRQ(ierr);
1964 ierr = DMSwarmGetField(swarm, field_name, &bs, NULL, &local_array_void); CHKERRQ(ierr);
1965
1966 // Determine the size of one element of the field's data type
1967 if (strcasecmp(field_name, "DMSwarm_pid") == 0) {
1968 element_size = sizeof(PetscInt64);
1969 } else if (strcasecmp(field_name, "DMSwarm_CellID") == 0 || strcasecmp(field_name, "DMSwarm_location_status") == 0) {
1970 element_size = sizeof(PetscInt);
1971 } else {
1972 element_size = sizeof(PetscScalar);
1973 }
1974
1975 if (rank == 0) {
1976 *n_total_particles = nglobal;
1977 *n_components = bs;
1978 *gathered_array = NULL; // Initialize output
1979 }
1980
1981 if (size == 1) { // Serial case is a simple copy
1982 if (rank == 0) {
1983 ierr = PetscMalloc(nglobal * bs * element_size, gathered_array); CHKERRQ(ierr);
1984 ierr = PetscMemcpy(*gathered_array, local_array_void, nglobal * bs * element_size); CHKERRQ(ierr);
1985 }
1986 } else { // Parallel case: use MPI_Gatherv
1987 PetscInt *recvcounts = NULL, *displs = NULL;
1988 if (rank == 0) {
1989 ierr = PetscMalloc1(size, &recvcounts); CHKERRQ(ierr);
1990 ierr = PetscMalloc1(size, &displs); CHKERRQ(ierr);
1991 }
1992 PetscInt sendcount = nlocal * bs;
1993
1994 // Gather the number of elements (not bytes) from each rank
1995 ierr = MPI_Gather(&sendcount, 1, MPIU_INT, recvcounts, 1, MPIU_INT, 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
1996
1997 if (rank == 0) {
1998 displs[0] = 0;
1999 // Convert counts and calculate displacements in terms of BYTES
2000 for (PetscMPIInt i = 0; i < size; i++) recvcounts[i] *= element_size;
2001 for (PetscMPIInt i = 1; i < size; i++) displs[i] = displs[i-1] + recvcounts[i-1];
2002
2003 ierr = PetscMalloc(nglobal * bs * element_size, gathered_array); CHKERRQ(ierr);
2004 }
2005
2006 // Use Gatherv with MPI_BYTE to handle any data type generically
2007 ierr = MPI_Gatherv(local_array_void, nlocal * bs * element_size, MPI_BYTE,
2008 *gathered_array, recvcounts, displs, MPI_BYTE,
2009 0, PETSC_COMM_WORLD); CHKERRQ(ierr);
2010
2011 if (rank == 0) {
2012 ierr = PetscFree(recvcounts); CHKERRQ(ierr);
2013 ierr = PetscFree(displs); CHKERRQ(ierr);
2014 }
2015 }
2016
2017 ierr = DMSwarmRestoreField(swarm, field_name, &bs, NULL, &local_array_void); CHKERRQ(ierr);
2018
2019 PetscFunctionReturn(0);
2020}
Here is the caller graph for this function:

◆ DisplayBanner()

PetscErrorCode DisplayBanner ( SimCtx *  simCtx)

Internal helper implementation: DisplayBanner().

Displays a structured banner summarizing the simulation configuration.

Local to this translation unit.

Definition at line 2026 of file io.c.

2027{
2028 PetscErrorCode ierr;
2029 PetscMPIInt rank;
2030 Cmpnts global_min_coords, global_max_coords;
2031 PetscReal StartTime;
2032 PetscInt StartStep,StepsToRun,total_num_particles;
2033 PetscMPIInt num_mpi_procs;
2034
2035 // SimCtx *simCtx = user->simCtx;
2036 UserCtx *user = simCtx->usermg.mgctx[simCtx->usermg.mglevels - 1].user;
2037 num_mpi_procs = simCtx->size;
2038 StartTime = simCtx->StartTime;
2039 StartStep = simCtx->StartStep;
2040 StepsToRun = simCtx->StepsToRun;
2041 total_num_particles = simCtx->np;
2042 BoundingBox *bboxlist_on_rank0 = simCtx->bboxlist;
2043
2044
2045 PetscFunctionBeginUser;
2046
2047 if (!user) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "DisplayBanner - UserCtx pointer is NULL.");
2048 global_min_coords = user->bbox.min_coords;
2049 global_max_coords = user->bbox.max_coords;
2050 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
2051
2052 if (rank == 0) {
2053 // If global_domain_bbox is not pre-populated in UserCtx, compute it here from bboxlist_on_rank0
2054 // This assumes bboxlist_on_rank0 is valid and contains all local bounding boxes on rank 0.
2055 if (bboxlist_on_rank0 && num_mpi_procs > 0) {
2056 global_min_coords = bboxlist_on_rank0[0].min_coords;
2057 global_max_coords = bboxlist_on_rank0[0].max_coords;
2058 for (PetscMPIInt p = 1; p < num_mpi_procs; ++p) {
2059 global_min_coords.x = PetscMin(global_min_coords.x, bboxlist_on_rank0[p].min_coords.x);
2060 global_min_coords.y = PetscMin(global_min_coords.y, bboxlist_on_rank0[p].min_coords.y);
2061 global_min_coords.z = PetscMin(global_min_coords.z, bboxlist_on_rank0[p].min_coords.z);
2062 global_max_coords.x = PetscMax(global_max_coords.x, bboxlist_on_rank0[p].max_coords.x);
2063 global_max_coords.y = PetscMax(global_max_coords.y, bboxlist_on_rank0[p].max_coords.y);
2064 global_max_coords.z = PetscMax(global_max_coords.z, bboxlist_on_rank0[p].max_coords.z);
2065 }
2066 // Optionally store this in user->global_domain_bbox if it's useful elsewhere
2067 // user->global_domain_bbox.min_coords = global_min_coords;
2068 // user->global_domain_bbox.max_coords = global_max_coords;
2069 } else {
2070 // Fallback or warning if bboxlist is not available for global calculation
2071 LOG_ALLOW(LOCAL, LOG_WARNING, "(Rank 0) - bboxlist not provided or num_mpi_procs <=0; using user->bbox for domain bounds.\n");
2072 // global_min_coords = user->bbox.min_coords; // Use local bbox of rank 0 as fallback
2073 // global_max_coords = user->bbox.max_coords;
2074 }
2075
2076 ierr = PetscPrintf(PETSC_COMM_SELF, "\n"); CHKERRQ(ierr);
2077 ierr = PetscPrintf(PETSC_COMM_SELF, "=============================================================\n"); CHKERRQ(ierr);
2078 ierr = PetscPrintf(PETSC_COMM_SELF, " CASE SUMMARY \n"); CHKERRQ(ierr);
2079 ierr = PetscPrintf(PETSC_COMM_SELF, "=============================================================\n"); CHKERRQ(ierr);
2080 ierr = PetscPrintf(PETSC_COMM_SELF, " Grid Points : %d X %d X %d\n", user->IM, user->JM, user->KM); CHKERRQ(ierr);
2081 ierr = PetscPrintf(PETSC_COMM_SELF, " Cells : %d X %d X %d\n", user->IM - 1, user->JM - 1, user->KM - 1); CHKERRQ(ierr);
2082 ierr = PetscPrintf(PETSC_COMM_SELF, " Global Domain Bounds (X) : %.6f to %.6f\n", (double)global_min_coords.x, (double)global_max_coords.x); CHKERRQ(ierr);
2083 ierr = PetscPrintf(PETSC_COMM_SELF, " Global Domain Bounds (Y) : %.6f to %.6f\n", (double)global_min_coords.y, (double)global_max_coords.y); CHKERRQ(ierr);
2084 ierr = PetscPrintf(PETSC_COMM_SELF, " Global Domain Bounds (Z) : %.6f to %.6f\n", (double)global_min_coords.z, (double)global_max_coords.z); CHKERRQ(ierr);
2085 ierr = PetscPrintf(PETSC_COMM_SELF, " Periodic Axes (BC-derived) : I=%s, J=%s, K=%s\n",
2086 simCtx->i_periodic ? "YES" : "NO",
2087 simCtx->j_periodic ? "YES" : "NO",
2088 simCtx->k_periodic ? "YES" : "NO"); CHKERRQ(ierr);
2089 for (PetscInt axis = 0; axis < 3; axis++) {
2090 if (!user->periodic_translation_valid[axis]) continue;
2091 ierr = PetscPrintf(PETSC_COMM_SELF,
2092 " Periodic %c Translation : (%.6e, %.6e, %.6e)\n",
2093 "IJK"[axis],
2094 (double)user->periodic_translation[axis].x,
2095 (double)user->periodic_translation[axis].y,
2096 (double)user->periodic_translation[axis].z); CHKERRQ(ierr);
2097 }
2098 if (total_num_particles > 0 &&
2099 (simCtx->i_periodic || simCtx->j_periodic || simCtx->k_periodic)) {
2100 ierr = PetscPrintf(PETSC_COMM_SELF,
2101 " Particle Periodic Wrapping : UNSUPPORTED (Eulerian periodicity only)\n"); CHKERRQ(ierr);
2102 }
2103 if(strcmp(simCtx->eulerianSource,"load")==0 || strcmp(simCtx->eulerianSource,"solve")==0){
2104 ierr = PetscPrintf(PETSC_COMM_SELF, "-------------------- Boundary Conditions --------------------\n"); CHKERRQ(ierr);
2105 const int face_name_width = 17; // Adjusted for longer names (Zeta,Eta,Xi)
2106 for (PetscInt i_face = 0; i_face < 6; ++i_face) {
2107 BCFace current_face = (BCFace)i_face;
2108 // The BCFaceToString will now return the Xi, Eta, Zeta versions
2109 const char* face_str = BCFaceToString(current_face);
2110 const char* bc_type_str = BCTypeToString(user->boundary_faces[current_face].mathematical_type);
2111 const char* bc_handler_type_str = BCHandlerTypeToString(user->boundary_faces[current_face].handler_type);
2112 if(user->boundary_faces[current_face].mathematical_type == INLET){
2113 if(user->boundary_faces[current_face].handler_type == BC_HANDLER_INLET_CONSTANT_VELOCITY){
2114 Cmpnts inlet_velocity = {0.0,0.0,0.0};
2115 PetscBool found;
2116 ierr = GetBCParamReal(user->boundary_faces[current_face].params,"vx",&inlet_velocity.x,&found); CHKERRQ(ierr);
2117 ierr = GetBCParamReal(user->boundary_faces[current_face].params,"vy",&inlet_velocity.y,&found); CHKERRQ(ierr);
2118 ierr = GetBCParamReal(user->boundary_faces[current_face].params,"vz",&inlet_velocity.z,&found); CHKERRQ(ierr);
2119 ierr = PetscPrintf(PETSC_COMM_SELF, " Face %-*s : %s - %s - [%.4f,%.4f,%.4f]\n",
2120 face_name_width, face_str, bc_type_str, bc_handler_type_str,inlet_velocity.x,inlet_velocity.y,inlet_velocity.z); CHKERRQ(ierr);
2121 } else if(user->boundary_faces[current_face].handler_type == BC_HANDLER_INLET_PARABOLIC){
2122 PetscReal v_max = 0.0;
2123 PetscBool found;
2124 ierr = GetBCParamReal(user->boundary_faces[current_face].params,"v_max",&v_max,&found); CHKERRQ(ierr);
2125 ierr = PetscPrintf(PETSC_COMM_SELF, " Face %-*s : %s - %s - v_max=%.4f\n",
2126 face_name_width, face_str, bc_type_str, bc_handler_type_str, v_max); CHKERRQ(ierr);
2127 } else if(user->boundary_faces[current_face].handler_type == BC_HANDLER_INLET_PROFILE_FROM_FILE){
2128 const char *source_file = "(missing)";
2129 for (BC_Param *param = user->boundary_faces[current_face].params; param; param = param->next) {
2130 if (strcasecmp(param->key, "source_file") == 0 && param->value) {
2131 source_file = param->value;
2132 break;
2133 }
2134 }
2135 ierr = PetscPrintf(PETSC_COMM_SELF, " Face %-*s : %s - %s - source_file=%s\n",
2136 face_name_width, face_str, bc_type_str, bc_handler_type_str, source_file); CHKERRQ(ierr);
2137 }
2138 } else if(user->boundary_faces[current_face].handler_type == BC_HANDLER_PERIODIC_DRIVEN_CONSTANT_FLUX){
2139 PetscReal flux;
2140 PetscBool trimflag,foundflux,foundtrimflag;
2141 ierr = GetBCParamReal(user->boundary_faces[current_face].params,"target_flux",&flux,&foundflux); CHKERRQ(ierr);
2142 ierr = GetBCParamBool(user->boundary_faces[current_face].params,"apply_trim",&trimflag,&foundtrimflag); CHKERRQ(ierr);
2143 ierr = PetscPrintf(PETSC_COMM_SELF, " Face %-*s : %s - %s - [%.4f] - %s\n",
2144 face_name_width, face_str, bc_type_str, bc_handler_type_str,flux,trimflag?"Trim Ucont":"No Trim"); CHKERRQ(ierr);
2145 } else{
2146 ierr = PetscPrintf(PETSC_COMM_SELF, " Face %-*s : %s - %s\n",
2147 face_name_width, face_str, bc_type_str,bc_handler_type_str); CHKERRQ(ierr);
2148 }
2149 }
2150 }
2151 ierr = PetscPrintf(PETSC_COMM_SELF, "-------------------------------------------------------------\n"); CHKERRQ(ierr);
2152 ierr = PetscPrintf(PETSC_COMM_SELF, " Run Mode : %s\n", simCtx->OnlySetup ? "SETUP ONLY" : "Full Simulation"); CHKERRQ(ierr);
2153 ierr = PetscPrintf(PETSC_COMM_SELF, " Start Time : %.4f\n", (double)StartTime); CHKERRQ(ierr);
2154 ierr = PetscPrintf(PETSC_COMM_SELF, " Timestep Size : %.4f\n", (double)simCtx->dt); CHKERRQ(ierr);
2155 ierr = PetscPrintf(PETSC_COMM_SELF, " Starting Step : %d\n", StartStep); CHKERRQ(ierr);
2156 ierr = PetscPrintf(PETSC_COMM_SELF, " Total Steps to Run : %d\n", StepsToRun); CHKERRQ(ierr);
2157 ierr = PetscPrintf(PETSC_COMM_SELF, " Ending Step : %d\n", StartStep + StepsToRun); CHKERRQ(ierr);
2158 if (simCtx->tiout > 0) {
2159 ierr = PetscPrintf(PETSC_COMM_SELF, " Field/Restart Cadence : every %d step(s)\n", simCtx->tiout); CHKERRQ(ierr);
2160 } else {
2161 ierr = PetscPrintf(PETSC_COMM_SELF, " Field/Restart Cadence : DISABLED\n"); CHKERRQ(ierr);
2162 }
2163 ierr = PetscPrintf(PETSC_COMM_SELF, " Immersed Boundary : %s\n", simCtx->immersed ? "ENABLED" : "DISABLED"); CHKERRQ(ierr);
2164 if (simCtx->walltimeGuardEnabled) {
2165 ierr = PetscPrintf(
2166 PETSC_COMM_SELF,
2167 " Runtime Walltime Guard : %s (warmup=%d, multiplier=%.2f, min=%.1f s, alpha=%.2f)\n",
2168 simCtx->walltimeGuardActive ? "ENABLED" : "CONFIGURED BUT INACTIVE",
2169 simCtx->walltimeGuardWarmupSteps,
2170 (double)simCtx->walltimeGuardMultiplier,
2171 (double)simCtx->walltimeGuardMinSeconds,
2172 (double)simCtx->walltimeGuardEstimatorAlpha
2173 ); CHKERRQ(ierr);
2174 } else {
2175 ierr = PetscPrintf(PETSC_COMM_SELF, " Runtime Walltime Guard : DISABLED\n"); CHKERRQ(ierr);
2176 }
2177 ierr = PetscPrintf(PETSC_COMM_SELF, " Number of MPI Processes : %d\n", num_mpi_procs); CHKERRQ(ierr);
2178 ierr = PetscPrintf(PETSC_COMM_WORLD," Number of Particles : %d\n", total_num_particles); CHKERRQ(ierr);
2179 if (simCtx->np > 0) {
2180 const char *particle_init_str = ParticleInitializationToString(simCtx->ParticleInitialization);
2181
2182 if (simCtx->particleConsoleOutputFreq > 0) {
2183 ierr = PetscPrintf(PETSC_COMM_SELF, " Particle Console Cadence : every %d step(s)\n", simCtx->particleConsoleOutputFreq); CHKERRQ(ierr);
2184 } else {
2185 ierr = PetscPrintf(PETSC_COMM_SELF, " Particle Console Cadence : DISABLED\n"); CHKERRQ(ierr);
2186 }
2187 ierr = PetscPrintf(PETSC_COMM_SELF, " Particle Log Row Sampling : every %d particle(s)\n", simCtx->LoggingFrequency); CHKERRQ(ierr);
2188 if (simCtx->StartStep > 0) {
2189 ierr = PetscPrintf(PETSC_COMM_SELF, " Particle Restart Mode : %s\n", simCtx->particleRestartMode); CHKERRQ(ierr);
2190 }
2191 ierr = PetscPrintf(PETSC_COMM_SELF, " Particle Initialization Mode: %s\n", particle_init_str); CHKERRQ(ierr);
2192 ierr = PetscPrintf(PETSC_COMM_SELF, " Interpolation Method : %s\n",
2193 simCtx->interpolationMethod == INTERP_TRILINEAR ? "Trilinear (direct cell-center)" : "CornerAveraged (legacy)"); CHKERRQ(ierr);
2194 if (simCtx->ParticleInitialization == PARTICLE_INIT_SURFACE_RANDOM ||
2195 simCtx->ParticleInitialization == PARTICLE_INIT_SURFACE_EDGES) {
2196 if (user->inletFaceDefined) {
2197 ierr = PetscPrintf(PETSC_COMM_SELF, " Particles Initialized At : %s (Enum Val: %d)\n", BCFaceToString(user->identifiedInletBCFace), user->identifiedInletBCFace); CHKERRQ(ierr);
2198 } else {
2199 ierr = PetscPrintf(PETSC_COMM_SELF, " Particles Initialized At : --- (No INLET face identified)\n"); CHKERRQ(ierr);
2200 }
2201 }
2202 }
2203 if(strcmp(simCtx->eulerianSource,"solve")==0 || strcmp(simCtx->eulerianSource,"load")==0){
2204 ierr = PetscPrintf(PETSC_COMM_WORLD," Reynolds Number : %le\n", simCtx->ren); CHKERRQ(ierr);
2205 //ierr = PetscPrintf(PETSC_COMM_WORLD," Von-Neumann Number : %le\n", simCtx->vnn); CHKERRQ(ierr);
2206 if(strcmp(simCtx->eulerianSource,"solve")==0){
2207 //ierr = PetscPrintf(PETSC_COMM_WORLD," Stanton Number : %le\n", simCtx->st); CHKERRQ(ierr);
2208 ierr = PetscPrintf(PETSC_COMM_WORLD," Momentum Equation Solver : %s\n", MomentumSolverTypeToString(simCtx->mom_solver_type)); CHKERRQ(ierr);
2209 ierr = PetscPrintf(PETSC_COMM_WORLD," Initial Pseudo-CFL (Courant): %le\n", simCtx->pseudo_cfl); CHKERRQ(ierr);
2210 ierr = PetscPrintf(PETSC_COMM_WORLD," Large Eddy Simulation Model : %s\n", LESModelToString(simCtx->les)); CHKERRQ(ierr);
2211 }
2212 if (strcmp(simCtx->eulerianSource, "load") == 0) {
2213 ierr = PetscPrintf(PETSC_COMM_SELF, " Eulerian State Source : load (%s)\n",
2214 simCtx->restart_dir); CHKERRQ(ierr);
2215 } else if (simCtx->StartStep > 0) {
2216 ierr = PetscPrintf(PETSC_COMM_SELF, " Eulerian State Source : restart step %d (%s)\n",
2217 simCtx->StartStep, simCtx->restart_dir); CHKERRQ(ierr);
2218 } else {
2219 const char* field_init_str = InitialConditionModeToString(simCtx->initialConditionMode);
2220 ierr = PetscPrintf(PETSC_COMM_SELF, " Eulerian State Source : initial condition (%s)\n",
2221 field_init_str); CHKERRQ(ierr);
2222 }
2223 if (strcmp(simCtx->eulerianSource, "solve") == 0 && simCtx->StartStep == 0 &&
2224 simCtx->initialConditionMode == IC_MODE_CONSTANT_CARTESIAN) {
2225 ierr = PetscPrintf(PETSC_COMM_SELF,
2226 " Constant Velocity (Cart.) : x=%.4f y=%.4f z=%.4f\n",
2227 (double)simCtx->InitialConstantContra.x,
2228 (double)simCtx->InitialConstantContra.y,
2229 (double)simCtx->InitialConstantContra.z); CHKERRQ(ierr);
2230 } else if (strcmp(simCtx->eulerianSource, "solve") == 0 && simCtx->StartStep == 0 &&
2231 simCtx->initialConditionMode == IC_MODE_CONSTANT_STREAMWISE) {
2232 ierr = PetscPrintf(PETSC_COMM_SELF,
2233 " Constant Velocity (Curv.) : speed=%.4f direction=%s\n",
2234 (double)simCtx->icVelocityPhysical,
2235 FlowDirectionToString(simCtx->flowDirection)); CHKERRQ(ierr);
2236 } else if (strcmp(simCtx->eulerianSource, "solve") == 0 && simCtx->StartStep == 0 &&
2237 simCtx->initialConditionMode == IC_MODE_POISEUILLE) {
2238 ierr = PetscPrintf(PETSC_COMM_SELF,
2239 " Poiseuille Peak Velocity : speed=%.4f direction=%s\n",
2240 (double)simCtx->icVelocityPhysical,
2241 FlowDirectionToString(simCtx->flowDirection)); CHKERRQ(ierr);
2242 } else if (strcmp(simCtx->eulerianSource, "solve") == 0 && simCtx->StartStep == 0 &&
2243 simCtx->initialConditionMode == IC_MODE_FILE) {
2244 ierr = PetscPrintf(PETSC_COMM_SELF,
2245 " Initial Velocity File : field=%s directory=%s\n",
2246 simCtx->initialConditionField == IC_FIELD_UCAT ? "Ucat" : "Ucont",
2247 simCtx->initialConditionDirectory); CHKERRQ(ierr);
2248 }
2249 } else if(strcmp(simCtx->eulerianSource,"analytical")==0){
2250 ierr = PetscPrintf(PETSC_COMM_WORLD," Analytical Solution Type : %s\n", simCtx->AnalyticalSolutionType); CHKERRQ(ierr);
2251 }
2252 ierr = PetscPrintf(PETSC_COMM_SELF, "=============================================================\n"); CHKERRQ(ierr);
2253 ierr = PetscPrintf(PETSC_COMM_SELF, "\n"); CHKERRQ(ierr);
2254 }
2255 PetscFunctionReturn(0);
2256}
GetBCParamReal
PetscErrorCode GetBCParamReal(BC_Param *params, const char *key, PetscReal *value_out, PetscBool *found)
Internal helper implementation: GetBCParamReal().
Definition io.c:390
GetBCParamBool
PetscErrorCode GetBCParamBool(BC_Param *params, const char *key, PetscBool *value_out, PetscBool *found)
Internal helper implementation: GetBCParamBool().
Definition io.c:411
BCFaceToString
const char * BCFaceToString(BCFace face)
Helper function to convert BCFace enum to a string representation.
Definition logging.c:669
BCTypeToString
const char * BCTypeToString(BCType type)
Helper function to convert BCType enum to a string representation.
Definition logging.c:772
FlowDirectionToString
const char * FlowDirectionToString(FlowDirection fd)
Convert a FlowDirection enum value to its YAML token string.
Definition logging.c:704
InitialConditionModeToString
const char * InitialConditionModeToString(InitialConditionMode mode)
Convert an initial-condition mode to a string representation.
Definition logging.c:687
LESModelToString
const char * LESModelToString(LESModelType LESFlag)
Helper function to convert LES Flag to a string representation.
Definition logging.c:740
MomentumSolverTypeToString
const char * MomentumSolverTypeToString(MomentumSolverType SolverFlag)
Helper function to convert Momentum Solver flag to a string representation.
Definition logging.c:756
ParticleInitializationToString
const char * ParticleInitializationToString(ParticleInitializationType ParticleInitialization)
Helper function to convert ParticleInitialization to a string representation.
Definition logging.c:723
SimCtx::icVelocityPhysical
PetscReal icVelocityPhysical
Definition variables.h:747
MGCtx::user
UserCtx * user
Definition variables.h:569
SimCtx::walltimeGuardActive
PetscBool walltimeGuardActive
Definition variables.h:839
SimCtx::initialConditionMode
InitialConditionMode initialConditionMode
Definition variables.h:742
PARTICLE_INIT_SURFACE_RANDOM
@ PARTICLE_INIT_SURFACE_RANDOM
Random placement on the inlet face.
Definition variables.h:550
PARTICLE_INIT_SURFACE_EDGES
@ PARTICLE_INIT_SURFACE_EDGES
Deterministic placement at inlet face edges.
Definition variables.h:553
SimCtx::StartTime
PetscReal StartTime
Definition variables.h:698
SimCtx::flowDirection
FlowDirection flowDirection
Definition variables.h:746
SimCtx::usermg
UserMG usermg
Definition variables.h:821
SimCtx::ren
PetscReal ren
Definition variables.h:732
BoundingBox::max_coords
Cmpnts max_coords
Maximum x, y, z coordinates of the bounding box.
Definition variables.h:171
SimCtx::dt
PetscReal dt
Definition variables.h:699
SimCtx::StepsToRun
PetscInt StepsToRun
Definition variables.h:695
SimCtx::StartStep
PetscInt StartStep
Definition variables.h:694
BoundingBox::min_coords
Cmpnts min_coords
Minimum x, y, z coordinates of the bounding box.
Definition variables.h:170
SimCtx::OnlySetup
PetscBool OnlySetup
Definition variables.h:700
Cmpnts::x
PetscScalar x
Definition variables.h:101
SimCtx::interpolationMethod
InterpolationMethod interpolationMethod
Definition variables.h:801
SimCtx::particleRestartMode
char particleRestartMode[16]
Definition variables.h:802
SimCtx::bboxlist
BoundingBox * bboxlist
Definition variables.h:799
SimCtx::eulerianSource
char eulerianSource[PETSC_MAX_PATH_LEN]
Definition variables.h:704
SimCtx::walltimeGuardEnabled
PetscBool walltimeGuardEnabled
Definition variables.h:838
SimCtx::ParticleInitialization
ParticleInitializationType ParticleInitialization
Definition variables.h:800
Cmpnts::z
PetscScalar z
Definition variables.h:101
INTERP_TRILINEAR
@ INTERP_TRILINEAR
Definition variables.h:563
UserMG::mglevels
PetscInt mglevels
Definition variables.h:576
SimCtx::AnalyticalSolutionType
char AnalyticalSolutionType[PETSC_MAX_PATH_LEN]
Definition variables.h:717
IC_MODE_CONSTANT_CARTESIAN
@ IC_MODE_CONSTANT_CARTESIAN
Definition variables.h:151
IC_MODE_POISEUILLE
@ IC_MODE_POISEUILLE
Definition variables.h:152
IC_MODE_CONSTANT_STREAMWISE
@ IC_MODE_CONSTANT_STREAMWISE
Definition variables.h:153
IC_MODE_FILE
@ IC_MODE_FILE
Definition variables.h:154
SimCtx::particleConsoleOutputFreq
PetscInt particleConsoleOutputFreq
Definition variables.h:697
SimCtx::InitialConstantContra
Cmpnts InitialConstantContra
Definition variables.h:745
Cmpnts::y
PetscScalar y
Definition variables.h:101
IC_FIELD_UCAT
@ IC_FIELD_UCAT
Definition variables.h:159
SimCtx::size
PetscMPIInt size
Definition variables.h:688
UserCtx::periodic_translation
Cmpnts periodic_translation[3]
Definition variables.h:892
UserMG::mgctx
MGCtx * mgctx
Definition variables.h:579
UserCtx::periodic_translation_valid
PetscBool periodic_translation_valid[3]
Definition variables.h:893
SimCtx::initialConditionField
InitialConditionField initialConditionField
Definition variables.h:743
UserCtx::bbox
BoundingBox bbox
Definition variables.h:887
SimCtx::mom_solver_type
MomentumSolverType mom_solver_type
Definition variables.h:724
SimCtx::immersed
PetscInt immersed
Definition variables.h:714
SimCtx::pseudo_cfl
PetscReal pseudo_cfl
Definition variables.h:732
SimCtx::LoggingFrequency
PetscInt LoggingFrequency
Definition variables.h:826
BoundingBox
Defines a 3D axis-aligned bounding box.
Definition variables.h:169
Cmpnts
A 3D point or vector with PetscScalar components.
Definition variables.h:100
UserCtx
User-defined context containing data specific to a single computational grid level.
Definition variables.h:876
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ParsePostProcessingSettings()

PetscErrorCode ParsePostProcessingSettings ( SimCtx *  simCtx)

Internal helper implementation: ParsePostProcessingSettings().

Initializes post-processing settings from a config file and command-line overrides.

Local to this translation unit.

Definition at line 2264 of file io.c.

2265{
2266 FILE *file;
2267 char line[1024];
2268 PetscBool startTimeSet, endTimeSet, timeStepSet;
2269
2270 PetscFunctionBeginUser;
2271 PROFILE_FUNCTION_BEGIN;
2272
2273 if (!simCtx || !simCtx->pps) {
2274 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_NULL, "SimCtx or its pps member is NULL in ParsePostProcessingSettings.");
2275 }
2276
2277 char *configFile = simCtx->PostprocessingControlFile;
2278 PostProcessParams *pps = simCtx->pps;
2279
2280
2281 // --- 1. Set Sane Defaults First ---
2282 pps->startTime = 0;
2283 pps->endTime = 0;
2284 pps->timeStep = 1;
2285 pps->outputParticles = PETSC_FALSE;
2286 pps->particle_output_freq = simCtx->LoggingFrequency; // Default to logging frequency;
2287 strcpy(pps->process_pipeline, "");
2288 strcpy(pps->output_fields_instantaneous, "Ucat,P");
2289 strcpy(pps->output_fields_averaged, "");
2290 strcpy(pps->output_prefix, "Field");
2291 strcpy(pps->particle_output_prefix,"Particle");
2292 strcpy(pps->particle_fields,"velocity,CellID,weight,pid");
2293 strcpy(pps->particle_pipeline,"");
2294 strncpy(pps->statistics_pipeline, "", MAX_PIPELINE_LENGTH - 1);
2295 strncpy(pps->statistics_output_prefix, "Stats", MAX_FILENAME_LENGTH - 1);
2296 strcpy(pps->particleExt,"dat"); // The input file format for particles.
2297 strcpy(pps->eulerianExt,"dat"); // The input file format for Eulerian fields.
2298 pps->reference[0] = pps->reference[1] = pps->reference[2] = 1;
2299 strncpy(pps->source_dir, simCtx->output_dir, sizeof(pps->source_dir) - 1);
2300 pps->source_dir[sizeof(pps->source_dir) - 1] = '\0'; // Ensure null-termination
2301
2302 // --- 2. Parse the Configuration File (overrides defaults) ---
2303 file = fopen(configFile, "r");
2304 if (file) {
2305 LOG_ALLOW(GLOBAL, LOG_INFO, "Parsing post-processing config file: %s\n", configFile);
2306 while (fgets(line, sizeof(line), file)) {
2307 char *key, *value, *comment;
2308 comment = strchr(line, '#'); if (comment) *comment = '\0';
2309 TrimWhitespace(line); if (strlen(line) == 0) continue;
2310 key = strtok(line, "="); value = strtok(NULL, "=");
2311 if (key && value) {
2312 TrimWhitespace(key); TrimWhitespace(value);
2313 if (strcmp(key, "startTime") == 0) pps->startTime = atoi(value);
2314 else if (strcmp(key, "endTime") == 0) pps->endTime = atoi(value);
2315 else if (strcmp(key, "timeStep") == 0) pps->timeStep = atoi(value);
2316 else if (strcmp(key, "output_particles") == 0) {
2317 if (strcasecmp(value, "true") == 0) pps->outputParticles = PETSC_TRUE;
2318 }
2319 else if (strcasecmp(key, "process_pipeline") == 0) {
2320 strncpy(pps->process_pipeline, value, MAX_PIPELINE_LENGTH - 1);
2321 pps->process_pipeline[MAX_PIPELINE_LENGTH - 1] = '\0'; // Ensure null-termination
2322 } else if (strcasecmp(key, "output_fields_instantaneous") == 0) {
2323 strncpy(pps->output_fields_instantaneous, value, MAX_FIELD_LIST_LENGTH - 1);
2324 pps->output_fields_instantaneous[MAX_FIELD_LIST_LENGTH - 1] = '\0';
2325 } else if (strcasecmp(key, "output_fields_averaged") == 0) {
2326 strncpy(pps->output_fields_averaged, value, MAX_FIELD_LIST_LENGTH - 1);
2327 pps->output_fields_averaged[MAX_FIELD_LIST_LENGTH - 1] = '\0';
2328 } else if (strcasecmp(key, "output_prefix") == 0) {
2329 strncpy(pps->output_prefix, value, MAX_FILENAME_LENGTH - 1);
2330 pps->output_prefix[MAX_FILENAME_LENGTH - 1] = '\0';
2331 } else if (strcasecmp(key, "particle_output_prefix") == 0) {
2332 strncpy(pps->particle_output_prefix, value, MAX_FILENAME_LENGTH - 1);
2333 pps->particle_output_prefix[MAX_FILENAME_LENGTH - 1] = '\0';
2334 } else if (strcasecmp(key, "particle_fields_instantaneous") == 0) {
2335 strncpy(pps->particle_fields, value, MAX_FIELD_LIST_LENGTH - 1);
2336 pps->particle_fields[MAX_FIELD_LIST_LENGTH - 1] = '\0';
2337 } else if (strcasecmp(key, "particle_pipeline") == 0) {
2338 strncpy(pps->particle_pipeline, value, MAX_PIPELINE_LENGTH - 1);
2339 pps->particle_pipeline[MAX_PIPELINE_LENGTH - 1] = '\0';
2340 } else if (strcasecmp(key, "particle_output_freq") == 0) {
2341 pps->particle_output_freq = atoi(value);
2342 } else if (strcasecmp(key, "statistics_pipeline") == 0) {
2343 strncpy(pps->statistics_pipeline, value, MAX_PIPELINE_LENGTH - 1);
2344 pps->statistics_pipeline[MAX_PIPELINE_LENGTH - 1] = '\0';
2345 } else if (strcasecmp(key, "statistics_output_prefix") == 0) {
2346 strncpy(pps->statistics_output_prefix, value, MAX_FILENAME_LENGTH - 1);
2347 pps->statistics_output_prefix[MAX_FILENAME_LENGTH - 1] = '\0';
2348 } else if (strcasecmp(key, "particleExt") == 0) {
2349 strncpy(pps->particleExt, value, sizeof(pps->particleExt) - 1);
2350 pps->particleExt[sizeof(pps->particleExt) - 1] = '\0';
2351 } else if (strcasecmp(key, "eulerianExt") == 0) {
2352 strncpy(pps->eulerianExt, value, sizeof(pps->eulerianExt) - 1);
2353 pps->eulerianExt[sizeof(pps->eulerianExt) - 1] = '\0';
2354 } else if (strcmp(key, "reference_ip") == 0) {pps->reference[0] = atoi(value);
2355 } else if (strcmp(key, "reference_jp") == 0) {pps->reference[1] = atoi(value);
2356 } else if (strcmp(key, "reference_kp") == 0) {pps->reference[2] = atoi(value);
2357 } else if (strcasecmp(key, "source_directory") == 0) {
2358 strncpy(pps->source_dir, value, sizeof(pps->source_dir) - 1);
2359 pps->source_dir[sizeof(pps->source_dir) - 1] = '\0';
2360 } else {
2361 LOG_ALLOW(GLOBAL, LOG_WARNING, "Unknown key '%s' in post-processing config file. Ignoring.\n", key);
2362 }
2363 // Add parsing for pipeline, fields, etc. in later phases
2364 }
2365 }
2366 fclose(file);
2367 } else {
2368 LOG_ALLOW(GLOBAL, LOG_WARNING, "Could not open post-processing config file '%s'. Using defaults and command-line overrides.\n", configFile);
2369 }
2370
2371 // --- 3. Parse Command-Line Options (overrides file settings and defaults) ---
2372 PetscOptionsGetInt(NULL, NULL, "-startTime", &pps->startTime, &startTimeSet);
2373 PetscOptionsGetInt(NULL, NULL, "-endTime", &pps->endTime, &endTimeSet);
2374 PetscOptionsGetInt(NULL, NULL, "-timeStep", &pps->timeStep, &timeStepSet);
2375 PetscOptionsGetBool(NULL, NULL, "-output_particles", &pps->outputParticles, NULL);
2376
2377 if(pps->endTime==-1){
2378 pps->endTime = simCtx->StartStep + simCtx->StepsToRun; // Total steps if endTime is set to -1.
2379 }
2380
2381 // If only startTime is given on command line, run for a single step
2382 if (startTimeSet && !endTimeSet) {
2383 pps->endTime = pps->startTime;
2384 }
2385
2386 LOG_ALLOW(GLOBAL, LOG_INFO, "Post-processing configured to run from t=%d to t=%d with step %d. Particle output: %s.\n",
2387 pps->startTime, pps->endTime, pps->timeStep, pps->outputParticles ? "TRUE" : "FALSE");
2388
2389 LOG_ALLOW(GLOBAL, LOG_INFO, "Process Pipeline: %s\n", pps->process_pipeline);
2390 LOG_ALLOW(GLOBAL, LOG_INFO, "Instantaneous Output Fields: %s\n", pps->output_fields_instantaneous);
2391 LOG_ALLOW(GLOBAL, LOG_INFO, "Output Prefix: %s\n", pps->output_prefix);
2392 LOG_ALLOW(GLOBAL, LOG_INFO, "Particle Output Prefix: %s\n", pps->particle_output_prefix);
2393 LOG_ALLOW(GLOBAL, LOG_INFO, "Particle Fields: %s\n", pps->particle_fields);
2394 LOG_ALLOW(GLOBAL, LOG_INFO, "Particle Pipeline: %s\n", pps->particle_pipeline);
2395 LOG_ALLOW(GLOBAL, LOG_INFO, "Particle Output Frequency: %d\n", pps->particle_output_freq);
2396 LOG_ALLOW(GLOBAL, LOG_INFO, "Averaged Output Fields (reserved): %s\n", pps->output_fields_averaged);
2397 LOG_ALLOW(GLOBAL, LOG_INFO, "Post input extensions: Eulerian='.%s', Particle='.%s'\n", pps->eulerianExt, pps->particleExt);
2398
2399 PROFILE_FUNCTION_END;
2400 PetscFunctionReturn(0);
2401}
TrimWhitespace
void TrimWhitespace(char *str)
Implementation of TrimWhitespace().
Definition io.c:38
PostProcessParams::statistics_output_prefix
char statistics_output_prefix[256]
basename for CSV output, e.g.
Definition variables.h:616
PostProcessParams::particle_output_prefix
char particle_output_prefix[256]
Definition variables.h:611
PostProcessParams::output_prefix
char output_prefix[256]
Definition variables.h:608
MAX_FIELD_LIST_LENGTH
#define MAX_FIELD_LIST_LENGTH
Definition variables.h:585
PostProcessParams::output_fields_averaged
char output_fields_averaged[1024]
Definition variables.h:607
MAX_PIPELINE_LENGTH
#define MAX_PIPELINE_LENGTH
Definition variables.h:584
PostProcessParams::reference
PetscInt reference[3]
Definition variables.h:623
PostProcessParams::timeStep
PetscInt timeStep
Definition variables.h:601
MAX_FILENAME_LENGTH
#define MAX_FILENAME_LENGTH
Definition variables.h:586
PostProcessParams::statistics_pipeline
char statistics_pipeline[1024]
e.g.
Definition variables.h:615
PostProcessParams::output_fields_instantaneous
char output_fields_instantaneous[1024]
Definition variables.h:606
PostProcessParams::particle_pipeline
char particle_pipeline[1024]
Definition variables.h:609
PostProcessParams::process_pipeline
char process_pipeline[1024]
Definition variables.h:605
PostProcessParams::particle_output_freq
PetscInt particle_output_freq
Definition variables.h:612
PostProcessParams::particle_fields
char particle_fields[1024]
Definition variables.h:610
PostProcessParams::outputParticles
PetscBool outputParticles
Definition variables.h:602
PostProcessParams::startTime
PetscInt startTime
Definition variables.h:599
SimCtx::PostprocessingControlFile
char PostprocessingControlFile[PETSC_MAX_PATH_LEN]
Definition variables.h:859
PostProcessParams
Holds all configuration parameters for a post-processing run.
Definition variables.h:594
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ParseScalingInformation()

PetscErrorCode ParseScalingInformation ( SimCtx *  simCtx)

Implementation of ParseScalingInformation().

Parses physical scaling parameters from command-line options.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
ParseScalingInformation()

Definition at line 2412 of file io.c.

2413{
2414 PetscErrorCode ierr;
2415 PetscBool flg;
2416
2417 PetscFunctionBeginUser;
2418 PROFILE_FUNCTION_BEGIN;
2419
2420 if (!simCtx) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "SimCtx is NULL in ParseScalingInformation");
2421
2422 // --- 1. Set default values to 1.0 ---
2423 // This represents a purely non-dimensional run if no scaling is provided.
2424 simCtx->scaling.L_ref = 1.0;
2425 simCtx->scaling.U_ref = 1.0;
2426 simCtx->scaling.rho_ref = 1.0;
2427
2428 // --- 2. Read overrides from the command line / control file ---
2429 ierr = PetscOptionsGetReal(NULL, NULL, "-scaling_L_ref", &simCtx->scaling.L_ref, &flg); CHKERRQ(ierr);
2430 ierr = PetscOptionsGetReal(NULL, NULL, "-scaling_U_ref", &simCtx->scaling.U_ref, &flg); CHKERRQ(ierr);
2431 ierr = PetscOptionsGetReal(NULL, NULL, "-scaling_rho_ref", &simCtx->scaling.rho_ref, &flg); CHKERRQ(ierr);
2432
2433 // --- 3. Calculate derived scaling factors ---
2434 // Check for division by zero to be safe, though U_ref should be positive.
2435 if (simCtx->scaling.U_ref <= 0.0) {
2436 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Reference velocity U_ref must be positive. Got %g", (double)simCtx->scaling.U_ref);
2437 }
2438 simCtx->scaling.P_ref = simCtx->scaling.rho_ref * simCtx->scaling.U_ref * simCtx->scaling.U_ref;
2439
2440 // --- 4. Log the final, effective scales for verification ---
2441 LOG(GLOBAL, LOG_INFO, "---------------- Physical Scales Initialized -----------------\n");
2442 LOG(GLOBAL, LOG_INFO, " L_ref: %.4f, U_ref: %.4f, rho_ref: %.4f, P_ref: %.4f\n",
2443 simCtx->scaling.L_ref, simCtx->scaling.U_ref, simCtx->scaling.rho_ref, simCtx->scaling.P_ref);
2444 LOG(GLOBAL, LOG_INFO, "--------------------------------------------------------------\n");
2445
2446 PROFILE_FUNCTION_END;
2447 PetscFunctionReturn(0);
2448}
LOG
#define LOG(scope, level, fmt,...)
Logging macro for PETSc-based applications with scope control.
Definition logging.h:83
ScalingCtx::L_ref
PetscReal L_ref
Definition variables.h:666
SimCtx::scaling
ScalingCtx scaling
Definition variables.h:763
ScalingCtx::P_ref
PetscReal P_ref
Definition variables.h:669
ScalingCtx::rho_ref
PetscReal rho_ref
Definition variables.h:668
ScalingCtx::U_ref
PetscReal U_ref
Definition variables.h:667
Here is the caller graph for this function:

◆ ReadDataFileToArray()

PetscInt ReadDataFileToArray ( const char *  filename,
double **  data_out,
PetscInt *  Nout,
MPI_Comm  comm 
)

Implementation of ReadDataFileToArray().

Reads a simple ASCII data file containing one numeric value per line.

Full API contract (arguments, ownership, side effects) is documented with the header declaration in include/io.h.

See also
ReadDataFileToArray()

Definition at line 2456 of file io.c.

2460{
2461 /* STEP 0: Prepare local variables & log function entry */
2462 PetscMPIInt rank, size;
2463 PetscErrorCode ierr;
2464 FILE *fp = NULL;
2465 PetscInt N = 0; /* number of lines/values read on rank 0 */
2466 double *array = NULL; /* pointer to local array on each rank */
2467 PetscInt fileExistsFlag = 0; /* 0 = doesn't exist, 1 = does exist */
2468
2469 LOG_ALLOW(GLOBAL, LOG_DEBUG,
2470 "Start reading from file: %s\n",
2471 filename);
2472
2473 /* Basic error checking: data_out, Nout must be non-null. */
2474 if (!filename || !data_out || !Nout) {
2475 LOG_ALLOW(GLOBAL, LOG_WARNING,
2476 "Null pointer argument provided.\n");
2477 return 1;
2478 }
2479
2480 /* Determine rank/size for coordinating I/O. */
2481 MPI_Comm_rank(comm, &rank);
2482 MPI_Comm_size(comm, &size);
2483
2484 /* STEP 1: On rank 0, check if file can be opened. */
2485 if (!rank) {
2486 fp = fopen(filename, "r");
2487 if (fp) {
2488 fileExistsFlag = 1;
2489 fclose(fp);
2490 }
2491 }
2492
2493 /* STEP 2: Broadcast file existence to all ranks. */
2494 // In ReadDataFileToArray:
2495 ierr = MPI_Bcast(&fileExistsFlag, 1, MPI_INT, 0, comm); CHKERRQ(ierr);
2496
2497 if (!fileExistsFlag) {
2498 /* If file does not exist, log & return. */
2499 if (!rank) {
2500 LOG_ALLOW(GLOBAL, LOG_WARNING,
2501 "File '%s' not found.\n",
2502 filename);
2503 }
2504 return 2;
2505 }
2506
2507 /* STEP 3: Rank 0 re-opens and reads the file, counting lines, etc. */
2508 if (!rank) {
2509 fp = fopen(filename, "r");
2510 if (!fp) {
2511 LOG_ALLOW(GLOBAL, LOG_WARNING,
2512 "File '%s' could not be opened for reading.\n",
2513 filename);
2514 return 3;
2515 }
2516
2517 /* (3a) Count lines first. */
2518 {
2519 char line[256];
2520 while (fgets(line, sizeof(line), fp)) {
2521 N++;
2522 }
2523 }
2524
2525 LOG_ALLOW(GLOBAL, LOG_DEBUG,
2526 "File '%s' has %d lines.\n",
2527 filename, N);
2528
2529 /* (3b) Allocate array on rank 0. */
2530 array = (double*)malloc(N * sizeof(double));
2531 if (!array) {
2532 fclose(fp);
2533 LOG_ALLOW(GLOBAL, LOG_WARNING,
2534 "malloc failed for array.\n");
2535 return 4;
2536 }
2537
2538 /* (3c) Rewind & read values into array. */
2539 rewind(fp);
2540 {
2541 PetscInt i = 0;
2542 char line[256];
2543 while (fgets(line, sizeof(line), fp)) {
2544 double val;
2545 if (sscanf(line, "%lf", &val) == 1) {
2546 array[i++] = val;
2547 }
2548 }
2549 }
2550 fclose(fp);
2551
2552 LOG_ALLOW(GLOBAL, LOG_INFO,
2553 "Successfully read %d values from '%s'.\n",
2554 N, filename);
2555 }
2556
2557 /* STEP 4: Broadcast the integer N to all ranks. */
2558 ierr = MPI_Bcast(&N, 1, MPI_INT, 0, comm); CHKERRQ(ierr);
2559
2560 /* STEP 5: Each rank allocates an array to receive the broadcast if rank>0. */
2561 if (rank) {
2562 array = (double*)malloc(N * sizeof(double));
2563 if (!array) {
2564 LOG_ALLOW(GLOBAL, LOG_WARNING,
2565 "malloc failed on rank %d.\n",
2566 rank);
2567 return 5;
2568 }
2569 }
2570
2571 /* STEP 6: Broadcast the actual data from rank 0 to all. */
2572 ierr = MPI_Bcast(array, N, MPI_DOUBLE, 0, comm); CHKERRQ(ierr);
2573
2574 /* STEP 7: Assign outputs on all ranks. */
2575 *data_out = array;
2576 *Nout = N;
2577
2578 LOG_ALLOW(GLOBAL, LOG_DEBUG,
2579 "Done. Provided array of length=%d to all ranks.\n",
2580 N);
2581 return 0; /* success */
2582}

◆ ReadPositionsFromFile()

PetscErrorCode ReadPositionsFromFile ( PetscInt  timeIndex,
UserCtx *  user,
double **  coordsArray,
PetscInt *  Ncoords 
)

Internal helper implementation: ReadPositionsFromFile().

Reads coordinate data (for particles) from file into a PETSc Vec, then gathers it to rank 0.

Local to this translation unit.

Definition at line 2588 of file io.c.

2592{
2593 PetscFunctionBeginUser;
2594
2595 PetscErrorCode ierr;
2596 Vec coordsVec;
2597
2598 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Creating coords Vec.\n");
2599 ierr = VecCreate(PETSC_COMM_WORLD, &coordsVec);CHKERRQ(ierr);
2600 ierr = VecSetFromOptions(coordsVec);CHKERRQ(ierr);
2601
2602 // For example: "position" is the name of the coordinate data
2603 ierr = ReadFieldData(user, "position", coordsVec, timeIndex, "dat");
2604 if (ierr) {
2605 LOG_ALLOW(GLOBAL, LOG_ERROR,
2606 "Error reading position data (ti=%d).\n",
2607 timeIndex);
2608 PetscFunctionReturn(ierr);
2609 }
2610
2611 LOG_ALLOW(GLOBAL, LOG_DEBUG, "ReadPositions - Gathering coords Vec to rank 0.\n");
2612 ierr = VecToArrayOnRank0(coordsVec, Ncoords, coordsArray);CHKERRQ(ierr);
2613
2614 ierr = VecDestroy(&coordsVec);CHKERRQ(ierr);
2615
2616 LOG_ALLOW(GLOBAL, LOG_DEBUG,
2617 "Successfully gathered coordinates. Ncoords=%d.\n", *Ncoords);
2618 PetscFunctionReturn(0);
2619}
VecToArrayOnRank0
PetscErrorCode VecToArrayOnRank0(Vec inVec, PetscInt *N, double **arrayOut)
Internal helper implementation: VecToArrayOnRank0().
Definition io.c:1863
Here is the call graph for this function:

◆ ReadFieldDataToRank0()

PetscErrorCode ReadFieldDataToRank0 ( PetscInt  timeIndex,
const char *  fieldName,
UserCtx *  user,
double **  scalarArray,
PetscInt *  Nscalars 
)

Internal helper implementation: ReadFieldDataToRank0().

Reads a named field from file into a PETSc Vec, then gathers it to rank 0.

Local to this translation unit.

Definition at line 2626 of file io.c.

2631{
2632 PetscFunctionBeginUser;
2633
2634 PetscErrorCode ierr;
2635 Vec fieldVec;
2636
2637 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Creating field Vec.\n");
2638 ierr = VecCreate(PETSC_COMM_WORLD, &fieldVec);CHKERRQ(ierr);
2639 ierr = VecSetFromOptions(fieldVec);CHKERRQ(ierr);
2640
2641 ierr = ReadFieldData(user, fieldName, fieldVec, timeIndex, "dat");
2642 if (ierr) {
2643 LOG_ALLOW(GLOBAL, LOG_ERROR,
2644 "Error reading field '%s' (ti=%d).\n",
2645 fieldName, timeIndex);
2646 PetscFunctionReturn(ierr);
2647 }
2648
2649 LOG_ALLOW(GLOBAL, LOG_DEBUG, "Gathering field Vec to rank 0.\n");
2650 ierr = VecToArrayOnRank0(fieldVec, Nscalars, scalarArray);CHKERRQ(ierr);
2651
2652 ierr = VecDestroy(&fieldVec);CHKERRQ(ierr);
2653
2654 LOG_ALLOW(GLOBAL, LOG_DEBUG,
2655 "Successfully gathered field '%s'. Nscalars=%d.\n",
2656 fieldName, *Nscalars);
2657 PetscFunctionReturn(0);
2658}
Here is the call graph for this function:

Variable Documentation

◆ g_nblk_from_file

PetscInt g_nblk_from_file = 0
static

Stores the number of blocks read from the grid file.

Definition at line 17 of file io.c.

◆ g_IMs_from_file

PetscInt* g_IMs_from_file = NULL
static

Caches the IM dimensions for all blocks read from the grid file.

Definition at line 19 of file io.c.

◆ g_JMs_from_file

PetscInt* g_JMs_from_file = NULL
static

Caches the JM dimensions for all blocks read from the grid file.

Definition at line 21 of file io.c.

◆ g_KMs_from_file

PetscInt* g_KMs_from_file = NULL
static

Caches the KM dimensions for all blocks read from the grid file.

Definition at line 23 of file io.c.

◆ g_file_has_been_read

PetscBool g_file_has_been_read = PETSC_FALSE
static

A flag to ensure the grid file is read only once.

Definition at line 25 of file io.c.

Generated by doxygen 1.9.8