PICurv 0.1.0
A Parallel Particle-In-Cell Solver for Curvilinear LES
Loading...
Searching...
No Matches
walkingsearch.h File Reference

Header file for particle location functions using the walking search algorithm. More...

#include <petsc.h>
#include <stdbool.h>
#include <math.h>
#include "variables.h"
#include "logging.h"
#include "setup.h"
Include dependency graph for walkingsearch.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Functions

PetscErrorCode GetCellCharacteristicSize (const Cell *cell, PetscReal *cellSize)
 Estimates a characteristic length of the cell for threshold scaling.
 
PetscErrorCode ComputeSignedDistanceToPlane (const Cmpnts v1, const Cmpnts v2, const Cmpnts v3, const Cmpnts v4, const Cmpnts cell_centroid, const Cmpnts p_target, PetscReal *d_signed, const PetscReal threshold)
 Computes the signed distance from a point to the plane approximating a quadrilateral face.
 
PetscErrorCode CalculateDistancesToCellFaces (const Cmpnts p, const Cell *cell, PetscReal *d, const PetscReal threshold)
 Computes the signed distances from a point to each face of a cubic cell.
 
PetscErrorCode DeterminePointPosition (PetscReal *d, PetscInt *result)
 Classifies a point based on precomputed face distances.
 
PetscErrorCode GetCellVerticesFromGrid (Cmpnts ***coor, PetscInt idx, PetscInt idy, PetscInt idz, Cell *cell)
 Retrieves the coordinates of the eight vertices of a cell based on grid indices.
 
PetscErrorCode InitializeTraversalParameters (UserCtx *user, Particle *particle, PetscInt *idx, PetscInt *idy, PetscInt *idz, PetscInt *traversal_steps)
 Initializes traversal parameters for locating a particle.
 
PetscErrorCode CheckCellWithinLocalGrid (UserCtx *user, PetscInt idx, PetscInt idy, PetscInt idz, PetscBool *is_within)
 Checks if the current cell indices are within the local grid boundaries.
 
PetscErrorCode RetrieveCurrentCell (UserCtx *user, PetscInt idx, PetscInt idy, PetscInt idz, Cell *cell)
 Retrieves the coordinates of the eight vertices of the current cell.
 
PetscErrorCode EvaluateParticlePosition (const Cell *cell, PetscReal *d, const Cmpnts p, PetscInt *position, const PetscReal threshold)
 Determines the spatial relationship of a particle relative to a cubic cell.
 
PetscErrorCode LocateParticleInGrid (UserCtx *user, Particle *particle)
 Locates the grid cell containing a particle using a walking search.
 
PetscErrorCode FinalizeTraversal (UserCtx *user, Particle *particle, PetscInt traversal_steps, PetscBool cell_found, PetscInt idx, PetscInt idy, PetscInt idz)
 Finalizes the traversal by reporting the results.
 
PetscErrorCode FindOwnerOfCell (UserCtx *user, PetscInt i, PetscInt j, PetscInt k, PetscMPIInt *owner_rank)
 Finds the MPI rank that owns a given global cell index.
 
PetscErrorCode LocateParticleOrFindMigrationTarget (UserCtx *user, Particle *particle, ParticleLocationStatus *status_out)
 Locates a particle's host cell or identifies its migration target using a robust walk search.
 
PetscErrorCode ReportSearchOutcome (const Particle *particle, ParticleLocationStatus status, PetscInt traversal_steps)
 Logs the final outcome of the particle location search.
 
PetscErrorCode UpdateCellIndicesBasedOnDistances (PetscReal d[NUM_FACES], PetscInt *idx, PetscInt *idy, PetscInt *idz)
 Updates the cell indices based on the signed distances to each face.
 

Detailed Description

Header file for particle location functions using the walking search algorithm.

This file contains declarations of functions that implement the walking search algorithm to locate particles within a computational grid. It includes functions for computing distances to cell faces, determining particle positions relative to cells, and updating traversal parameters.

Definition in file walkingsearch.h.

Function Documentation

◆ GetCellCharacteristicSize()

PetscErrorCode GetCellCharacteristicSize ( const Cell cell,
PetscReal *  cellSize 
)

Estimates a characteristic length of the cell for threshold scaling.

For a hexahedral cell with vertices cell->vertices[0..7], we approximate the cell size by some measure—e.g. average edge length or diagonal.

Parameters
[in]cellA pointer to the Cell structure
[out]cellSizeA pointer to a PetscReal where the characteristic size is stored.
Returns
PetscErrorCode 0 on success, non-zero on failure.

Definition at line 26 of file walkingsearch.c.

27{
28 PetscErrorCode ierr = 0;
29 PetscFunctionBeginUser;
31 if (!cell || !cellSize) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
32 "GetCellCharacteristicSize: Null pointer(s).");
33
34 // A simple approach: compute the average of the distances between
35 // all pairs of adjacent vertices that share an edge. That gives
36 // a typical measure of cell dimension. For a uniform grid, you
37 // could do something simpler.
38
39 PetscInt edges[12][2] = {
40 {0,1},{1,2},{2,3},{3,0}, // bottom face edges
41 {4,5},{5,6},{6,7},{7,4}, // top face edges
42 {0,7},{1,6},{2,5},{3,4} // vertical edges
43 };
44
45 PetscReal totalEdgeLen = 0.0;
46 for (PetscInt i=0; i<12; i++) {
47 PetscInt vA = edges[i][0];
48 PetscInt vB = edges[i][1];
49 Cmpnts A = cell->vertices[vA];
50 Cmpnts B = cell->vertices[vB];
51 PetscReal dx = B.x - A.x;
52 PetscReal dy = B.y - A.y;
53 PetscReal dz = B.z - A.z;
54 PetscReal edgeLen = sqrt(dx*dx + dy*dy + dz*dz);
55 totalEdgeLen += edgeLen;
56 }
57
58 // Average edge length
59 *cellSize = totalEdgeLen / 12.0;
60
62 PetscFunctionReturn(ierr);
63}
#define PROFILE_FUNCTION_END
Marks the end of a profiled code block.
Definition logging.h:746
#define PROFILE_FUNCTION_BEGIN
Marks the beginning of a profiled code block (typically a function).
Definition logging.h:737
PetscScalar x
Definition variables.h:101
PetscScalar z
Definition variables.h:101
PetscScalar y
Definition variables.h:101
Cmpnts vertices[8]
Coordinates of the eight vertices of the cell.
Definition variables.h:161
A 3D point or vector with PetscScalar components.
Definition variables.h:100
Here is the caller graph for this function:

◆ ComputeSignedDistanceToPlane()

PetscErrorCode ComputeSignedDistanceToPlane ( const Cmpnts  v1,
const Cmpnts  v2,
const Cmpnts  v3,
const Cmpnts  v4,
const Cmpnts  cell_centroid,
const Cmpnts  p_target,
PetscReal *  d_signed,
const PetscReal  threshold 
)

Computes the signed distance from a point to the plane approximating a quadrilateral face.

This function calculates the signed distance from a given point p_target to the plane approximating a quadrilateral face defined by points v1, v2, v3, and v4. The plane's initial normal is determined using two edge vectors emanating from v1 (i.e., v2-v1 and v4-v1).

Normal Orientation: The initial normal's orientation is checked against the cell's interior. A vector from the face centroid to the cell centroid (vec_face_to_cell_centroid) is computed. If the dot product of the initial normal and vec_face_to_cell_centroid is positive, it means the initial normal is pointing towards the cell's interior (it's an inward normal relative to the cell). In this case, the normal is flipped to ensure it points outward.

Signed Distance Convention: The signed distance is the projection of the vector from p_target to the face's centroid onto the (now guaranteed) outward-pointing plane's normal vector.

  • d_signed < 0: p_target is "outside" and "beyond" the face, in the direction of the outward normal. This is the primary case indicating the particle has crossed this face.
  • d_signed > 0: p_target is "outside" but "behind" the face (on the same side as the cell interior relative to this face plane).
  • d_signed = 0: p_target is on the face plane (within threshold).
Parameters
[in]v1First vertex defining the face (used as origin for initial normal calculation).
[in]v2Second vertex defining the face (used for first edge vector: v2-v1).
[in]v3Third vertex defining the face (used for face centroid calculation).
[in]v4Fourth vertex defining the face (used for second edge vector: v4-v1).
[in]cell_centroidThe centroid of the 8-vertex cell to which this face belongs.
[in]p_targetThe point from which the distance to the plane is calculated.
[out]d_signedPointer to store the computed signed distance.
[in]thresholdThe threshold below which the absolute distance is considered zero.
Note
  • The order of vertices v1, v2, v3, v4 is used for calculating the face centroid. The order of v1, v2, v4 is used for the initial normal calculation.
  • The cell_centroid is crucial for robustly determining the outward direction of the normal.
Returns
PetscErrorCode Returns 0 on success, PETSC_ERR_ARG_NULL if d_signed is NULL, or PETSC_ERR_USER if a degenerate plane (near-zero normal) is detected.

Definition at line 109 of file walkingsearch.c.

112{
113 PetscErrorCode ierr = 0;
114 PetscMPIInt rank;
115
116 PetscFunctionBeginUser;
118
119 if (d_signed == NULL) {
120 ierr = MPI_Comm_rank(MPI_COMM_WORLD, &rank); CHKERRQ(ierr);
121 LOG_ALLOW(LOCAL, LOG_ERROR, "Output pointer 'd_signed' is NULL on rank %d.\n", rank);
122 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Output pointer 'd_signed' must not be NULL.");
123 }
124
125 LOG_ALLOW(LOCAL, LOG_VERBOSE, "Target point: (%.6e, %.6e, %.6e)\n", p_target.x, p_target.y, p_target.z);
126 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Cell Centroid: (%.6e, %.6e, %.6e)\n", cell_centroid.x, cell_centroid.y, cell_centroid.z);
127 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Face Vertices:\n");
128 LOG_ALLOW(LOCAL, LOG_VERBOSE, " v1: (%.6e, %.6e, %.6e)\n", v1.x, v1.y, v1.z);
129 LOG_ALLOW(LOCAL, LOG_VERBOSE, " v2: (%.6e, %.6e, %.6e)\n", v2.x, v2.y, v2.z);
130 LOG_ALLOW(LOCAL, LOG_VERBOSE, " v3: (%.6e, %.6e, %.6e)\n", v3.x, v3.y, v3.z);
131 LOG_ALLOW(LOCAL, LOG_VERBOSE, " v4: (%.6e, %.6e, %.6e)\n", v4.x, v4.y, v4.z);
132
133 // --- Calculate Edge Vectors for Initial Normal Computation ---
134 PetscReal edge1_x = v2.x - v1.x;
135 PetscReal edge1_y = v2.y - v1.y;
136 PetscReal edge1_z = v2.z - v1.z;
137
138 PetscReal edge2_x = v4.x - v1.x;
139 PetscReal edge2_y = v4.y - v1.y;
140 PetscReal edge2_z = v4.z - v1.z;
141 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Edge1 (v2-v1): (%.6e, %.6e, %.6e)\n", edge1_x, edge1_y, edge1_z);
142 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Edge2 (v4-v1): (%.6e, %.6e, %.6e)\n", edge2_x, edge2_y, edge2_z);
143
144 // --- Compute Initial Normal Vector (Cross Product: edge1 x edge2) ---
145 PetscReal normal_x_initial = edge1_y * edge2_z - edge1_z * edge2_y;
146 PetscReal normal_y_initial = edge1_z * edge2_x - edge1_x * edge2_z;
147 PetscReal normal_z_initial = edge1_x * edge2_y - edge1_y * edge2_x;
148 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Initial Raw Normal (edge1 x edge2): (%.6e, %.6e, %.6e)\n", normal_x_initial, normal_y_initial, normal_z_initial);
149
150 PetscReal normal_magnitude = sqrt(normal_x_initial * normal_x_initial +
151 normal_y_initial * normal_y_initial +
152 normal_z_initial * normal_z_initial);
153 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Initial Normal Magnitude: %.6e\n", normal_magnitude);
154
155 if (normal_magnitude < 1.0e-12) {
156 ierr = MPI_Comm_rank(MPI_COMM_WORLD, &rank); CHKERRQ(ierr);
158 "Degenerate plane detected on rank %d. Normal magnitude (%.3e) is too small.\n",
159 rank, normal_magnitude);
160 LOG_ALLOW(LOCAL, LOG_ERROR, " Offending vertices for normal: v1(%.3e,%.3e,%.3e), v2(%.3e,%.3e,%.3e), v4(%.3e,%.3e,%.3e)\n",
161 v1.x,v1.y,v1.z, v2.x,v2.y,v2.z, v4.x,v4.y,v4.z);
162 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_USER, "Degenerate plane detected (normal vector is near zero).");
163 }
164
165 // --- Calculate the Centroid of the Four Face Vertices (v1, v2, v3, v4) ---
166 PetscReal face_centroid_x = 0.25 * (v1.x + v2.x + v3.x + v4.x);
167 PetscReal face_centroid_y = 0.25 * (v1.y + v2.y + v3.y + v4.y);
168 PetscReal face_centroid_z = 0.25 * (v1.z + v2.z + v3.z + v4.z);
169 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Face Centroid: (%.6e, %.6e, %.6e)\n", face_centroid_x, face_centroid_y, face_centroid_z);
170
171 // --- Orient the Normal to Point Outward from the Cell ---
172 // Vector from face centroid to cell centroid
173 PetscReal vec_fc_to_cc_x = cell_centroid.x - face_centroid_x;
174 PetscReal vec_fc_to_cc_y = cell_centroid.y - face_centroid_y;
175 PetscReal vec_fc_to_cc_z = cell_centroid.z - face_centroid_z;
176 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Vec (FaceCentroid -> CellCentroid): (%.6e, %.6e, %.6e)\n", vec_fc_to_cc_x, vec_fc_to_cc_y, vec_fc_to_cc_z);
177
178 // Dot product of initial normal with vector from face centroid to cell centroid
179 PetscReal dot_prod_orientation = normal_x_initial * vec_fc_to_cc_x +
180 normal_y_initial * vec_fc_to_cc_y +
181 normal_z_initial * vec_fc_to_cc_z;
182 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Dot Product for Orientation (N_initial . Vec_FC_to_CC): %.6e\n", dot_prod_orientation);
183
184 PetscReal normal_x = normal_x_initial;
185 PetscReal normal_y = normal_y_initial;
186 PetscReal normal_z = normal_z_initial;
187
188 // If dot_prod_orientation > 0, initial normal points towards cell interior, so flip it.
189 if (dot_prod_orientation > 1.0e-9) { // Use a small epsilon to avoid issues if dot product is extremely close to zero
190 normal_x = -normal_x_initial;
191 normal_y = -normal_y_initial;
192 normal_z = -normal_z_initial;
193 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Initial normal was inward (dot_prod > 0). Flipped normal.\n");
194 } else if (dot_prod_orientation == 0.0 && normal_magnitude > 1e-12) {
195 // This case is ambiguous or face plane contains cell centroid.
196 // This might happen for highly symmetric cells or if face_centroid IS cell_centroid (e.g. 2D cell).
197 // For now, we keep the original normal direction based on v1,v2,v4 ordering.
198 // A more robust solution for this edge case might be needed if it occurs often.
199 ierr = MPI_Comm_rank(MPI_COMM_WORLD, &rank); CHKERRQ(ierr);
200 LOG_ALLOW(LOCAL, LOG_WARNING, "Rank %d: Dot product for normal orientation is zero. Normal direction might be ambiguous. Keeping initial normal direction from (v2-v1)x(v4-v1).\n", rank);
201 }
202 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Oriented Raw Normal: (%.6e, %.6e, %.6e)\n", normal_x, normal_y, normal_z);
203
204
205 // --- Normalize the (Now Outward-Pointing) Normal Vector ---
206 // Note: normal_magnitude was calculated from initial normal.
207 // If we flipped, magnitude is the same.
208 normal_x /= normal_magnitude;
209 normal_y /= normal_magnitude;
210 normal_z /= normal_magnitude;
211 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Normalized Outward Normal: (%.6f, %.6f, %.6f)\n", normal_x, normal_y, normal_z);
212
213
214 // --- Compute Vector from Target Point to Face Centroid ---
215 PetscReal vec_p_to_fc_x = face_centroid_x - p_target.x;
216 PetscReal vec_p_to_fc_y = face_centroid_y - p_target.y;
217 PetscReal vec_p_to_fc_z = face_centroid_z - p_target.z;
218
219 // --- Compute Signed Distance ---
220 *d_signed = vec_p_to_fc_x * normal_x +
221 vec_p_to_fc_y * normal_y +
222 vec_p_to_fc_z * normal_z;
223
224 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Raw Signed Distance (using outward normal): %.15e\n", *d_signed);
225
226 // --- Apply Threshold ---
227 if (fabs(*d_signed) < threshold) {
228 LOG_ALLOW(LOCAL, LOG_VERBOSE, " Distance %.15e is less than threshold %.1e. Setting to 0.0.\n", *d_signed, threshold);
229 *d_signed = 0.0;
230 }
231
232 LOG_ALLOW(LOCAL, LOG_VERBOSE, "Final Signed Distance: %.15e\n", *d_signed);
233
235 PetscFunctionReturn(0);
236}
#define LOCAL
Logging scope definitions for controlling message output.
Definition logging.h:45
#define LOG_ALLOW(scope, level, fmt,...)
Logging macro that checks both the log level and whether the calling function is in the allowed-funct...
Definition logging.h:200
@ LOG_ERROR
Critical errors that may halt the program.
Definition logging.h:28
@ LOG_WARNING
Non-critical issues that warrant attention.
Definition logging.h:29
@ LOG_VERBOSE
Extremely detailed logs, typically for development use only.
Definition logging.h:34
Here is the caller graph for this function:

◆ CalculateDistancesToCellFaces()

PetscErrorCode CalculateDistancesToCellFaces ( const Cmpnts  p,
const Cell cell,
PetscReal *  d,
const PetscReal  threshold 
)

Computes the signed distances from a point to each face of a cubic cell.

This function calculates the signed distances from a specified point p to each of the six faces of a cubic cell. The cell is defined by its eight vertices, and the distances are stored in the memory location pointed to by d.

Parameters
[in]pThe target point in 3D space.
[in]cellPointer to a Cell structure that defines the cubic cell via its vertices.
[out]dA pointer to an array of six PetscReal values where the computed signed distances will be stored.
[in]thresholdA PetscReal value that specifies the minimum distance below which a computed distance is considered to be zero.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

This function calculates the signed distances from a specified point p to each of the six faces of a cubic cell. The cell is defined by its eight vertices, and the distances are stored in the memory location pointed to by d. Each distance corresponds to a specific face of the cell, as enumerated by the Face enumeration. A threshold value is used to determine when a distance should be considered effectively zero, accounting for floating-point precision limitations.

Parameters
[in]pThe target point in 3D space for which distances to the cell's faces are to be calculated. Represented by the Cmpnts structure.
[in]cellA pointer to a Cell structure that defines the cubic cell via its eight vertices. The vertices must be ordered consistently to ensure correct normal directions for each face.
[out]dA pointer to an array of six PetscReal values where the computed signed distances will be stored. Each index in the array corresponds to a specific face as defined by the Face enumeration:
  • d[LEFT]: Distance to the Left Face of the cell.
  • d[RIGHT]: Distance to the Right Face of the cell.
  • d[BOTTOM]: Distance to the Bottom Face of the cell.
  • d[TOP]: Distance to the Top Face of the cell.
  • d[FRONT]: Distance to the Front Face of the cell.
  • d[BACK]: Distance to the Back Face of the cell.
[in]thresholdA PetscReal value that specifies the minimum distance below which a computed distance is considered to be zero. This helps in mitigating issues arising from floating-point arithmetic inaccuracies.
Returns
PetscErrorCode Returns 0 if the function executes successfully. If an error occurs, a non-zero error code is returned, indicating the type of failure.
Note
  • The vertices of the cell must be ordered in a counter-clockwise manner when viewed from outside the cell. This ordering ensures that the normal vectors computed for each face point outward, which is essential for correctly determining the sign of the distances.
  • All four vertices defining a face must lie on a non-degenerate (i.e., non-flat) plane. Degenerate planes can lead to undefined normal vectors and incorrect distance calculations.
  • The threshold parameter provides flexibility in handling cases where the point p is extremely close to a face, allowing the user to define what constitutes "close enough" to be treated as zero distance based on the specific requirements of their application.

Definition at line 282 of file walkingsearch.c.

283{
284 PetscErrorCode ierr;
285 PetscFunctionBeginUser;
287
288 // Validate that the 'cell' pointer is not NULL to prevent dereferencing a null pointer.
289 if (cell == NULL) {
290 LOG_ALLOW(LOCAL,LOG_ERROR, "'cell' is NULL.\n");
291 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "'cell' is NULL.");
292 }
293
294 // Validate that the 'd' pointer is not NULL to ensure there is memory allocated for distance storage.
295 if (d == NULL) {
296 LOG_ALLOW(LOCAL,LOG_ERROR, "'d' is NULL.\n");
297 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "'d' is NULL.");
298 }
299
300 // Compute the centroid of the entire cell
301 Cmpnts cell_centroid = {0.0, 0.0, 0.0};
302 for (int i = 0; i < 8; ++i) {
303 cell_centroid.x += cell->vertices[i].x;
304 cell_centroid.y += cell->vertices[i].y;
305 cell_centroid.z += cell->vertices[i].z;
306 }
307 cell_centroid.x /= 8.0;
308 cell_centroid.y /= 8.0;
309 cell_centroid.z /= 8.0;
310
311 LOG_ALLOW(LOCAL,LOG_DEBUG, "Cell Centroid: (%.6e, %.6e, %.6e)\n",
312 cell_centroid.x, cell_centroid.y, cell_centroid.z);
313
314
315 // Compute the signed distance from point 'p' to the BACK face of the cell.
316 // The BACK face is defined by vertices 0, 3, 2, and 1, with its normal vector pointing in the -z direction.
318 cell->vertices[0], // Vertex 0
319 cell->vertices[3], // Vertex 3
320 cell->vertices[2], // Vertex 2
321 cell->vertices[1], // Vertex 1
322 cell_centroid, // Cell centroid
323 p, // Target point
324 &d[BACK], // Storage location for the BACK face distance
325 threshold // Threshold for zero distance
326 ); CHKERRQ(ierr);
327
328 // Compute the signed distance from point 'p' to the FRONT face of the cell.
329 // The FRONT face is defined by vertices 4, 7, 6, and 5, with its normal vector pointing in the +z direction.
331 cell->vertices[4], // Vertex 4
332 cell->vertices[7], // Vertex 7
333 cell->vertices[6], // Vertex 6
334 cell->vertices[5], // Vertex 5
335 cell_centroid, // Cell centroid
336 p, // Target point
337 &d[FRONT], // Storage location for the FRONT face distance
338 threshold // Threshold for zero distance
339 ); CHKERRQ(ierr);
340
341 // Compute the signed distance from point 'p' to the BOTTOM face of the cell.
342 // The BOTTOM face is defined by vertices 0, 1, 6, and 7, with its normal vector pointing in the -y direction.
344 cell->vertices[0], // Vertex 0
345 cell->vertices[1], // Vertex 1
346 cell->vertices[6], // Vertex 6
347 cell->vertices[7], // Vertex 7
348 cell_centroid, // Cell centroid
349 p, // Target point
350 &d[BOTTOM], // Storage location for the BOTTOM face distance
351 threshold // Threshold for zero distance
352 ); CHKERRQ(ierr);
353
354 // Compute the signed distance from point 'p' to the TOP face of the cell.
355 // The TOP face is defined by vertices 3, 4, 5, and 2, with its normal vector pointing in the +y direction.
357 cell->vertices[3], // Vertex 3
358 cell->vertices[4], // Vertex 4
359 cell->vertices[5], // Vertex 5
360 cell->vertices[2], // Vertex 2
361 cell_centroid, // Cell centroid
362 p, // Target point
363 &d[TOP], // Storage location for the TOP face distance
364 threshold // Threshold for zero distance
365 ); CHKERRQ(ierr);
366
367 // Compute the signed distance from point 'p' to the LEFT face of the cell.
368 // The LEFT face is defined by vertices 0, 7, 4, and 3, with its normal vector pointing in the -x direction.
370 cell->vertices[0], // Vertex 0
371 cell->vertices[7], // Vertex 7
372 cell->vertices[4], // Vertex 4
373 cell->vertices[3], // Vertex 3
374 cell_centroid, // Cell centroid
375 p, // Target point
376 &d[LEFT], // Storage location for the LEFT face distance
377 threshold // Threshold for zero distance
378 ); CHKERRQ(ierr);
379
380 // Compute the signed distance from point 'p' to the RIGHT face of the cell.
381 // The RIGHT face is defined by vertices 1, 2, 5, and 6, with its normal vector pointing in the +x direction.
383 cell->vertices[1], // Vertex 1
384 cell->vertices[2], // Vertex 2
385 cell->vertices[5], // Vertex 5
386 cell->vertices[6], // Vertex 6
387 cell_centroid, // Cell centroid
388 p, // Target point
389 &d[RIGHT], // Storage location for the RIGHT face distance
390 threshold // Threshold for zero distance
391 ); CHKERRQ(ierr);
392
393 LOG_ALLOW(LOCAL, LOG_DEBUG, "Populated d: "
394 "d[LEFT=%d]=%.3e, d[RIGHT=%d]=%.3e, d[BOTTOM=%d]=%.3e, "
395 "d[TOP=%d]=%.3e, d[FRONT=%d]=%.3e, d[BACK=%d]=%.3e\n",
396 LEFT, d[LEFT], RIGHT, d[RIGHT], BOTTOM, d[BOTTOM],
397 TOP, d[TOP], FRONT, d[FRONT], BACK, d[BACK]);
398 LOG_ALLOW(LOCAL, LOG_DEBUG, "Raw d: "
399 "[%.3e, %.3e, %.3e, %.3e, %.3e, %.3e]\n",
400 d[0], d[1], d[2], d[3], d[4], d[5]);
401
403 PetscFunctionReturn(0); // Indicate successful execution of the function.
404}
@ LOG_DEBUG
Detailed debugging information.
Definition logging.h:32
@ TOP
Definition variables.h:145
@ FRONT
Definition variables.h:145
@ BOTTOM
Definition variables.h:145
@ BACK
Definition variables.h:145
@ LEFT
Definition variables.h:145
@ RIGHT
Definition variables.h:145
PetscErrorCode ComputeSignedDistanceToPlane(const Cmpnts v1, const Cmpnts v2, const Cmpnts v3, const Cmpnts v4, const Cmpnts cell_centroid, const Cmpnts p_target, PetscReal *d_signed, const PetscReal threshold)
Computes the signed distance from a point to the plane approximating a quadrilateral face.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ DeterminePointPosition()

PetscErrorCode DeterminePointPosition ( PetscReal *  d,
PetscInt *  result 
)

Classifies a point based on precomputed face distances.

Given an array of six distances d[NUM_FACES] from the point to each face of a hexahedral cell, this function determines: 0 => inside 1 => on a single face 2 => on an edge (2 faces = 0) 3 => on a corner (3+ faces = 0) -1 => outside

Parameters
[in]dArray of six face distances.
[out]resultPointer to an integer classification: {0,1,2,3,-1}.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

Definition at line 424 of file walkingsearch.c.

425{
426
427 PetscFunctionBeginUser;
429 // Validate input pointers
430 if (d == NULL) {
431 LOG_ALLOW(LOCAL,LOG_ERROR, "'d' is NULL.\n");
432 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Input parameter 'd' is NULL.");
433 }
434 if (result == NULL) {
435 LOG_ALLOW(LOCAL,LOG_ERROR, "'result' is NULL.\n");
436 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Output parameter 'result' is NULL.");
437 }
438
439 // Initialize flags
440 PetscBool isInside = PETSC_TRUE;
441 PetscBool isOnBoundary = PETSC_FALSE;
442 PetscInt IntersectionCount = 0; // Counts the number of intersections of the point with various planes of the cell.
443
444 // Analyze distances to determine position
445 for(int i = 0; i < NUM_FACES; i++) {
446 if(d[i] < 0.0) {
447 isInside = PETSC_FALSE; // Point is outside in at least one direction
448 }
449 if(d[i] == 0.0) {
450 isOnBoundary = PETSC_TRUE; // Point is exactly at least one face
451 IntersectionCount++;
452 }
453 }
454
455 // Set the result based on flags
456 if(isInside) {
457 *result = 0; // Inside the cell
458 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Particle is inside the cell.\n");
459 }
460 else if(isOnBoundary) {
461 if(IntersectionCount == 1){
462 *result = 1; // on face
463 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Particle is on a face of the cell.\n");
464 }
465 else if(IntersectionCount == 2){
466 *result = 2; // on edge
467 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Particle is on an edge of the cell.\n");
468 }
469 else if(IntersectionCount >= 3){
470 *result = 3; // on corner
471 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Particle is on a corner of the cell.\n");
472 }
473 }
474 else {
475 *result = -1; // Outside the cell
476 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Particle is outside the cell.\n");
477 }
478
480 PetscFunctionReturn(0); // Indicate successful execution
481}
@ NUM_FACES
Definition variables.h:145
Here is the caller graph for this function:

◆ GetCellVerticesFromGrid()

PetscErrorCode GetCellVerticesFromGrid ( Cmpnts ***  coor,
PetscInt  idx,
PetscInt  idy,
PetscInt  idz,
Cell cell 
)

Retrieves the coordinates of the eight vertices of a cell based on grid indices.

This function populates the cell structure with the coordinates of the eight vertices of a hexahedral cell given its grid indices (idx, idy, idz).

Parameters
[in]coorPointer to a 3D array of Cmpnts structures representing the grid coordinates.
[in]idxThe i-index of the cell.
[in]idyThe j-index of the cell.
[in]idzThe k-index of the cell.
[out]cellPointer to a Cell structure to store the cell's vertices.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

This function populates the cell structure with the coordinates of the eight vertices of a hexahedral cell given its grid indices (idx, idy, idz). The vertex numbering follows the convention depicted in the provided figure:

(i,j+1,k)
3---------------------------2(i+1,j+1,k)
/| /|
/ | / |
/ | / |
/ | / |
/ | / |
/ | / |
/ | / |
4-------|------------------5 |
| | | |
|(i,j+1,k+1) | |
| | | |
| 0-------------------|-------1 (i+1,j,k)
| / | /
| / | /
| / | /
| / | /
| / | /
| / | /
| / |/
7---------------------------6
(i,j,k+1) (i+1,j,k+1)

Vertex Numbering and Positions:

  • Vertex 0: (i, j, k) → cell->vertices[0]
  • Vertex 1: (i+1, j, k) → cell->vertices[1]
  • Vertex 2: (i+1, j+1, k) → cell->vertices[2]
  • Vertex 3: (i, j+1, k) → cell->vertices[3]
  • Vertex 4: (i, j+1, k+1) → cell->vertices[4]
  • Vertex 5: (i+1, j+1, k+1) → cell->vertices[5]
  • Vertex 6: (i+1, j, k+1) → cell->vertices[6]
  • Vertex 7: (i, j, k+1) → cell->vertices[7]
Parameters
[in]coorPointer to a 3D array of Cmpnts structures representing the grid coordinates.
[in]idxThe i-index of the cell.
[in]idyThe j-index of the cell.
[in]idzThe k-index of the cell.
[out]cellPointer to a Cell structure to store the cell's vertices.
Returns
PetscErrorCode Returns 0 to indicate successful execution. Non-zero on failure.
Note
  • It is assumed that the grid indices (idx, idy, idz) are within valid bounds.
  • The coor array should be properly initialized with accurate coordinates for all grid points.
  • The function assumes unit spacing between grid points. Modify the coordinate assignments if your grid spacing differs.
  • The boundary checks have been removed as they are performed before calling this function.

Definition at line 542 of file walkingsearch.c.

544{
545
546 PetscFunctionBeginUser;
548 // Validate input pointers
549 if (coor == NULL) {
550 LOG_ALLOW(LOCAL,LOG_ERROR, "'coor' is NULL.\n");
551 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Input array 'coor' is NULL.");
552 }
553 if (cell == NULL) {
554 LOG_ALLOW(LOCAL,LOG_ERROR, "'cell' is NULL.\n");
555 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Output parameter 'cell' is NULL.");
556 }
557
558 // Assign vertices based on grid indices matching the figure
559 // Vertex numbering follows the convention depicted in the figure
560 cell->vertices[0] = coor[idz][idy][idx]; // Vertex 0: (i, j, k)
561 cell->vertices[1] = coor[idz][idy][idx+1]; // Vertex 1: (i+1, j, k)
562 cell->vertices[2] = coor[idz][idy+1][idx+1]; // Vertex 2: (i+1, j+1, k)
563 cell->vertices[3] = coor[idz][idy+1][idx]; // Vertex 3: (i, j+1, k)
564 cell->vertices[4] = coor[idz+1][idy+1][idx]; // Vertex 4: (i, j+1, k+1)
565 cell->vertices[5] = coor[idz+1][idy+1][idx+1]; // Vertex 5: (i+1, j+1, k+1)
566 cell->vertices[6] = coor[idz+1][idy][idx+1]; // Vertex 6: (i+1, j, k+1)
567 cell->vertices[7] = coor[idz+1][idy][idx]; // Vertex 7: (i, j, k+1)
568
569 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Retrieved vertices for cell (%d, %d, %d).\n", idx, idy, idz);
570
572 PetscFunctionReturn(0); // Indicate successful execution
573}
Here is the caller graph for this function:

◆ InitializeTraversalParameters()

PetscErrorCode InitializeTraversalParameters ( UserCtx user,
Particle particle,
PetscInt *  idx,
PetscInt *  idy,
PetscInt *  idz,
PetscInt *  traversal_steps 
)

Initializes traversal parameters for locating a particle.

This function sets the initial cell indices based on the local grid boundaries and initializes the traversal step counter.

Parameters
[in]userPointer to the user-defined context containing grid information.
[in]particlePointer to the Particle structure containing its location.
[out]idxPointer to store the initial i-index of the cell.
[out]idyPointer to store the initial j-index of the cell.
[out]idzPointer to store the initial k-index of the cell.
[out]traversal_stepsPointer to store the initial traversal step count.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

This function sets the initial cell indices for the grid search.

  • If the particle has valid previous cell indices (not -1,-1,-1), the search starts from that previous cell.
  • Otherwise (e.g., first time locating the particle), the search starts from the beginning corner (xs, ys, zs) of the local grid domain owned by the current process. It also initializes the traversal step counter.
Parameters
[in]userPointer to the user-defined context containing grid information.
[in]particlePointer to the Particle structure containing its location and previous cell (if any).
[out]idxPointer to store the initial i-index of the cell.
[out]idyPointer to store the initial j-index of the cell.
[out]idzPointer to store the initial k-index of the cell.
[out]traversal_stepsPointer to store the initial traversal step count.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

Definition at line 597 of file walkingsearch.c.

598{
599 PetscErrorCode ierr;
600 DMDALocalInfo info;
601
602 PetscFunctionBeginUser;
604
605 // Validate input pointers
606 if (user == NULL || particle == NULL || idx == NULL || idy == NULL || idz == NULL || traversal_steps == NULL) {
607 LOG_ALLOW(LOCAL,LOG_ERROR, "One or more input pointers are NULL.\n");
608 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "One or more input pointers are NULL.");
609 }
610
611 // Get grid information (needed for fallback start and potentially validation)
612 ierr = DMDAGetLocalInfo(user->da, &info); CHKERRQ(ierr);
613
614 // --- Check if the particle has a valid previous cell ID ---
615 // Assuming particle->cell stores the *global* indices
616 if (particle->cell[0] >= 0 && particle->cell[1] >= 0 && particle->cell[2] >= 0) {
617
618 // It sees a valid cell ID exists.
619 // Before using it, it explicitly checks if this cell is accessible.
620 PetscBool is_handoff_cell_valid;
621 ierr = CheckCellWithinLocalGrid(user, particle->cell[0], particle->cell[1], particle->cell[2], &is_handoff_cell_valid); CHKERRQ(ierr);
622
623 if (is_handoff_cell_valid) {
624 // FAST PATH: The check passed. The cell is safe. Use it.
625 *idx = particle->cell[0];
626 *idy = particle->cell[1];
627 *idz = particle->cell[2];
628
629 LOG_ALLOW(LOCAL,LOG_DEBUG, "Particle %lld has previous cell (%d, %d, %d). Starting search there.\n",
630 (long long)particle->PID, *idx, *idy, *idz); // Cast id to long long for printing PetscInt64
631
632 } else {
633 // SAFE FALLBACK: The check failed! The handoff cell is NOT safe.
634 // Discard the unsafe handoff cell and revert to starting at the
635 // guaranteed-safe local corner.
636 *idx = info.xs;
637 *idy = info.ys;
638 *idz = info.zs;
639 }
640 } else {
641 // No valid previous cell ID (e.g., -1,-1,-1 or first time).
642 // IMPROVED: Use distance-based multi-seed approach to find the closest strategic cell.
643 // Sample the 6 face centers + volume center (7 total candidates) and start from the closest.
644
645 Cmpnts ***cent;
646
647 // Get local cell center array (lCent)
648 // For cell (i,j,k), the center is at cent[k+1][j+1][i+1]
649 ierr = DMDAVecGetArrayRead(user->fda, user->lCent, &cent); CHKERRQ(ierr);
650
651 // Define 7 strategic sampling points: 6 face centers + 1 volume center
652 PetscInt candidates[7][3];
653 const char* candidate_names[7] = {"X-min face", "X-max face", "Y-min face",
654 "Y-max face", "Z-min face", "Z-max face", "Volume center"};
655
656 // Calculate middle indices for each dimension
657 PetscInt mid_i = info.xs + (info.xm - 1) / 2;
658 PetscInt mid_j = info.ys + (info.ym - 1) / 2;
659 PetscInt mid_k = info.zs + (info.zm - 1) / 2;
660 PetscInt max_i = info.xs + info.xm - 2;
661 PetscInt max_j = info.ys + info.ym - 2;
662 PetscInt max_k = info.zs + info.zm - 2;
663
664 // Clamp max indices to ensure valid cells
665 if (max_i < info.xs) max_i = info.xs;
666 if (max_j < info.ys) max_j = info.ys;
667 if (max_k < info.zs) max_k = info.zs;
668
669 // X-min face center (i=min, j=mid, k=mid)
670 candidates[0][0] = info.xs; candidates[0][1] = mid_j; candidates[0][2] = mid_k;
671 // X-max face center (i=max, j=mid, k=mid)
672 candidates[1][0] = max_i; candidates[1][1] = mid_j; candidates[1][2] = mid_k;
673 // Y-min face center (i=mid, j=min, k=mid)
674 candidates[2][0] = mid_i; candidates[2][1] = info.ys; candidates[2][2] = mid_k;
675 // Y-max face center (i=mid, j=max, k=mid)
676 candidates[3][0] = mid_i; candidates[3][1] = max_j; candidates[3][2] = mid_k;
677 // Z-min face center (i=mid, j=mid, k=min)
678 candidates[4][0] = mid_i; candidates[4][1] = mid_j; candidates[4][2] = info.zs;
679 // Z-max face center (i=mid, j=mid, k=max)
680 candidates[5][0] = mid_i; candidates[5][1] = mid_j; candidates[5][2] = max_k;
681 // Volume center (i=mid, j=mid, k=mid)
682 candidates[6][0] = mid_i; candidates[6][1] = mid_j; candidates[6][2] = mid_k;
683
684 // Find the closest candidate cell
685 PetscReal min_distance = PETSC_MAX_REAL;
686 PetscInt best_candidate = 6; // Default to volume center
687
688 for (PetscInt i = 0; i < 7; i++) {
689 PetscInt ci = candidates[i][0];
690 PetscInt cj = candidates[i][1];
691 PetscInt ck = candidates[i][2];
692
693 // Get cell center coordinates: cent[k+1][j+1][i+1] for cell (i,j,k)
694 Cmpnts cell_center = cent[ck+1][cj+1][ci+1];
695
696 // Calculate Euclidean distance from particle to cell center
697 PetscReal dx = particle->loc.x - cell_center.x;
698 PetscReal dy = particle->loc.y - cell_center.y;
699 PetscReal dz = particle->loc.z - cell_center.z;
700 PetscReal distance = sqrt(dx*dx + dy*dy + dz*dz);
701
702 LOG_ALLOW(LOCAL, LOG_DEBUG, " Candidate %d (%s) at cell (%d,%d,%d): center=(%.3e,%.3e,%.3e), distance=%.3e\n",
703 i, candidate_names[i], ci, cj, ck, cell_center.x, cell_center.y, cell_center.z, distance);
704
705 if (distance < min_distance) {
706 min_distance = distance;
707 best_candidate = i;
708 *idx = ci;
709 *idy = cj;
710 *idz = ck;
711 }
712 }
713
714 // Restore array
715 ierr = DMDAVecRestoreArrayRead(user->fda, user->lCent, &cent); CHKERRQ(ierr);
716
717 LOG_ALLOW(LOCAL, LOG_INFO, "Particle %lld has no valid previous cell. Multi-seed search selected %s at cell (%d,%d,%d) with distance %.3e.\n",
718 (long long)particle->PID, candidate_names[best_candidate], *idx, *idy, *idz, min_distance);
719 }
720
721 // Initialize traversal step counter
722 *traversal_steps = 0;
723
724 // Log the chosen starting point
725 LOG_ALLOW(LOCAL,LOG_INFO, "Traversal for particle %lld initialized to start at cell (%d, %d, %d).\n",
726 (long long)particle->PID, *idx, *idy, *idz);
727
728 PetscFunctionReturn(0);
729}
@ LOG_INFO
Informational messages about program execution.
Definition logging.h:31
Vec lCent
Definition variables.h:776
PetscInt cell[3]
Definition variables.h:167
Cmpnts loc
Definition variables.h:168
PetscInt64 PID
Definition variables.h:166
PetscErrorCode CheckCellWithinLocalGrid(UserCtx *user, PetscInt idx, PetscInt idy, PetscInt idz, PetscBool *is_within)
Checks if the current GLOBAL CELL indices are within the LOCAL GHOSTED grid boundaries accessible by ...
Here is the call graph for this function:
Here is the caller graph for this function:

◆ CheckCellWithinLocalGrid()

PetscErrorCode CheckCellWithinLocalGrid ( UserCtx user,
PetscInt  idx,
PetscInt  idy,
PetscInt  idz,
PetscBool *  is_within 
)

Checks if the current cell indices are within the local grid boundaries.

Parameters
[in]userPointer to the user-defined context containing grid information.
[in]idxThe i-index of the current cell.
[in]idyThe j-index of the current cell.
[in]idzThe k-index of the current cell.
[out]is_withinPointer to a PetscBool that will be set to PETSC_TRUE if within bounds, else PETSC_FALSE.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

Checks if the current cell indices are within the local grid boundaries.

This function determines if the provided global cell indices (idx, idy, idz) fall within the range of cells covered by the current process's owned and ghost NODES. A cell C(i,j,k) (origin N(i,j,k)) is considered within this ghosted region if its origin node N(i,j,k) is within the rank's ghosted nodal region, AND node N(i+1,j+1,k+1) is also within it (to ensure the entire cell extent is covered by available node data). More simply, we check if the cell's origin node is within the range of nodes that can form the start of a ghosted cell.

Parameters
[in]userPointer to the user-defined context (needs user->fda for node info).
[in]idxThe global i-index of the cell's origin node.
[in]idyThe global j-index of the cell's origin node.
[in]idzThe global k-index of the cell's origin node.
[out]is_withinPointer to a PetscBool that will be set to PETSC_TRUE if within ghosted bounds, else PETSC_FALSE.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

Definition at line 753 of file walkingsearch.c.

754{
755 PetscErrorCode ierr;
756 DMDALocalInfo info_nodes; // Node information from the DMDA that defines ghost regions (user->fda)
757
758 PetscFunctionBeginUser; // Assuming this is part of your PETSc style
760
761 // Validate inputs
762 if (user == NULL || is_within == NULL) {
763 LOG_ALLOW(LOCAL, LOG_ERROR, "Input pointer is NULL.\n");
764 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Input pointer is NULL in CheckCellWithinLocalGrid.");
765 }
766 if (user->fda == NULL) {
767 LOG_ALLOW(LOCAL, LOG_ERROR, "user->fda is NULL.Cannot get ghost info.\n");
768 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "user->fda is NULL. Cannot get ghost info.");
769 }
770
771 // Get node info from user->fda (this DMDA has the ghost layer information for nodes)
772 ierr = DMDAGetLocalInfo(user->fda, &info_nodes); CHKERRQ(ierr);
773
774 // Determine the range of GLOBAL CELL INDICES that are covered by this rank's ghosted NODAL region.
775 // A cell C(i,j,k) has origin node N(i,j,k).
776 // The ghosted nodal region starts at global node index info_nodes.gxs and has info_nodes.gxm nodes.
777
778 // Global starting index of the first cell whose origin node is within the ghosted nodal region.
779 PetscInt gxs_cell_global_start = info_nodes.gxs;
780 PetscInt gys_cell_global_start = info_nodes.gys;
781 PetscInt gzs_cell_global_start = info_nodes.gzs;
782
783 // Number of cells that can be formed starting from nodes within the ghosted nodal region.
784 // If there are N ghosted nodes (info_nodes.gxm), they can be origins for N-1 cells.
785 PetscInt gxm_cell_local_count = (info_nodes.gxm > 0) ? info_nodes.gxm - 1 : 0;
786 PetscInt gym_cell_local_count = (info_nodes.gym > 0) ? info_nodes.gym - 1 : 0;
787 PetscInt gzm_cell_local_count = (info_nodes.gzm > 0) ? info_nodes.gzm - 1 : 0;
788
789 // Global exclusive end index for cells whose origins are in the ghosted nodal region.
790 PetscInt gxe_cell_global_end_exclusive = gxs_cell_global_start + gxm_cell_local_count;
791 PetscInt gye_cell_global_end_exclusive = gys_cell_global_start + gym_cell_local_count;
792 PetscInt gze_cell_global_end_exclusive = gzs_cell_global_start + gzm_cell_local_count;
793
794 // Check if the given global cell index (idx, idy, idz) falls within this range.
795 // This means the origin node of cell (idx,idy,idz) is within the rank's accessible ghosted node region,
796 // and that node can indeed serve as a cell origin (i.e., it's not the very last node in the ghosted region).
797 if (idx >= gxs_cell_global_start && idx < gxe_cell_global_end_exclusive &&
798 idy >= gys_cell_global_start && idy < gye_cell_global_end_exclusive &&
799 idz >= gzs_cell_global_start && idz < gze_cell_global_end_exclusive) {
800 *is_within = PETSC_TRUE;
801 } else {
802 *is_within = PETSC_FALSE;
803 }
804
805 LOG_ALLOW(LOCAL, LOG_DEBUG, "Cell (origin node glob idx) (%d, %d, %d) is %s the ghosted local grid (covered cell origins x:[%d..%d), y:[%d..%d), z:[%d..%d)).\n",
806 idx, idy, idz, (*is_within) ? "within" : "outside",
807 gxs_cell_global_start, gxe_cell_global_end_exclusive,
808 gys_cell_global_start, gye_cell_global_end_exclusive,
809 gzs_cell_global_start, gze_cell_global_end_exclusive);
810
812 PetscFunctionReturn(0);
813}
Here is the caller graph for this function:

◆ RetrieveCurrentCell()

PetscErrorCode RetrieveCurrentCell ( UserCtx user,
PetscInt  idx,
PetscInt  idy,
PetscInt  idz,
Cell cell 
)

Retrieves the coordinates of the eight vertices of the current cell.

Parameters
[in]userPointer to the user-defined context containing grid information.
[in]idxThe i-index of the current cell.
[in]idyThe j-index of the current cell.
[in]idzThe k-index of the current cell.
[out]cellPointer to a Cell structure to store the cell's vertices.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

Definition at line 828 of file walkingsearch.c.

829{
830 PetscErrorCode ierr;
831 Vec Coor;
832 Cmpnts ***coor_array;
833 PetscMPIInt rank;
834 DMDALocalInfo info_nodes;
835
836 PetscFunctionBeginUser;
838
839 // Validate input pointers
840 if (user == NULL || cell == NULL) {
841 LOG_ALLOW(LOCAL,LOG_ERROR, "One or more input pointers are NULL.\n");
842 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "One or more input pointers are NULL.");
843 }
844
845 // Get local coordinates
846 ierr = DMGetCoordinatesLocal(user->da, &Coor); CHKERRQ(ierr);
847 ierr = DMDAVecGetArrayRead(user->fda, Coor, &coor_array); CHKERRQ(ierr);
848
849 // Get the local grid information FOR THE GHOSTED ARRAY from user->fda
850 // This info contains the mapping between global and local indices.
851 ierr = DMDAGetLocalInfo(user->fda, &info_nodes); CHKERRQ(ierr);
852
853 PetscInt idx_local = idx; // - info_nodes.gxs;
854 PetscInt idy_local = idy; // - info_nodes.gys;
855 PetscInt idz_local = idz; // - info_nodes.gzs;
856
857 LOG_ALLOW(LOCAL,LOG_VERBOSE," [Rank %d] Getting vertex coordinates for cell: %d,%d,%d whose local coordinates are %d,%d,%d.\n",rank,idx,idy,idz,idx_local,idy_local,idz_local);
858
859 // Get the current cell's vertices
860 ierr = GetCellVerticesFromGrid(coor_array, idx_local, idy_local, idz_local, cell); CHKERRQ(ierr);
861
862 // Restore array
863 ierr = DMDAVecRestoreArrayRead(user->fda, Coor, &coor_array); CHKERRQ(ierr);
864
865 // Get MPI rank for debugging
866 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
867
868 // Debug: Print cell vertices
869 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Cell (%d, %d, %d) vertices \n", idx, idy, idz);
870 ierr = LOG_CELL_VERTICES(cell, rank); CHKERRQ(ierr);
871
873 PetscFunctionReturn(0);
874}
PetscErrorCode LOG_CELL_VERTICES(const Cell *cell, PetscMPIInt rank)
Prints the coordinates of a cell's vertices.
Definition logging.c:187
PetscErrorCode GetCellVerticesFromGrid(Cmpnts ***coor, PetscInt idx, PetscInt idy, PetscInt idz, Cell *cell)
Retrieves the coordinates of the eight vertices of a cell based on grid indices.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ EvaluateParticlePosition()

PetscErrorCode EvaluateParticlePosition ( const Cell cell,
PetscReal *  d,
const Cmpnts  p,
PetscInt *  position,
const PetscReal  threshold 
)

Determines the spatial relationship of a particle relative to a cubic cell.

This function evaluates whether a particle located at a specific point p in 3D space is positioned inside the cell, on the boundary of the cell, or outside the cell. The determination is based on the signed distances from the particle to each of the six faces of the cell. The function utilizes the computed distances to ascertain the particle's position with respect to the cell boundaries, considering a threshold to account for floating-point precision.

Parameters
[in]cellA pointer to a Cell structure that defines the cubic cell via its vertices. The cell's geometry is essential for accurately computing the distances to each face.
[in]dA pointer to an array of six PetscReal values that store the signed distances from the particle to each face of the cell. These distances are typically computed using the CalculateDistancesToCellFaces function.
[in]pThe location of the particle in 3D space, represented by the Cmpnts structure. This point is the reference for distance calculations to the cell's faces.
[out]positionA pointer to an integer that will be set based on the particle's position relative to the cell:
  • 0: The particle is inside the cell.
  • 1: The particle is on the boundary of the cell.
  • -1: The particle is outside the cell.
[in]thresholdA PetscReal value that defines the minimum distance below which a computed distance is considered to be zero. This threshold helps in mitigating inaccuracies due to floating-point arithmetic, especially when determining if the particle lies exactly on the boundary.
Returns
PetscErrorCode Returns 0 if the function executes successfully. If an error occurs, a non-zero error code is returned, indicating the type of failure.
Note
  • It is assumed that the d array has been properly allocated and contains valid distance measurements before calling this function.
  • The function relies on CalculateDistancesToCellFaces to accurately compute the signed distances to each face. Any inaccuracies in distance calculations can affect the determination of the particle's position.
  • The threshold parameter should be chosen based on the specific precision requirements of the application to balance between sensitivity and robustness against floating-point errors.
  • The function includes a debug statement that prints the face distances, which can be useful for verifying the correctness of distance computations during development or troubleshooting.

Definition at line 923 of file walkingsearch.c.

924{
925 PetscErrorCode ierr;
926 PetscReal cellSize;
927 PetscReal cellThreshold;
928 PetscFunctionBeginUser;
930
931 // Validate input pointers to ensure they are not NULL, preventing potential segmentation faults.
932 if (cell == NULL || position == NULL) {
933 LOG_ALLOW(LOCAL,LOG_ERROR, "One or more input pointers are NULL.\n");
934 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "One or more input pointers are NULL.");
935 }
936
937 // Compute a local cell size
938 ierr = GetCellCharacteristicSize(cell, &cellSize); CHKERRQ(ierr);
939
940 // scale base threshold by cell size (ex: if threshold = 1e-6 and cell size = 0.01, cellThreshold = 1e-8)
941 cellThreshold = threshold*cellSize;
942
943 // Invoke the function to calculate signed distances from the particle to each face of the cell.
944 // The distances are stored in the array pointed to by 'd'.
945 ierr = CalculateDistancesToCellFaces(p, cell, d, cellThreshold); CHKERRQ(ierr);
946 CHKERRQ(ierr); // Check for errors in distance calculation.
947
948
949 // Catch degenerate-plane error manually:
950 if (ierr == PETSC_ERR_USER) {
952 "Skipping cell due to degenerate face.\n");
953 // We can set *position = -1 here
954 *position = -1; // treat as outside
955 return 0; // not a fatal error, just skip
956 } else {
957 CHKERRQ(ierr);
958 }
959
960
961 // Debugging output: Print the computed distances to each face for verification purposes.
962 LOG_ALLOW(LOCAL,LOG_DEBUG, "Face Distances:\n");
964 CHKERRQ(ierr); // Check for errors in printing distances.
965
966 // Determine the particle's position relative to the cell based on the computed distances.
967 // The function sets the value pointed to by 'position' accordingly:
968 // 0 for inside, 1 (2 or 3) for on the boundary, and -1 for outside.
969 ierr = DeterminePointPosition(d,position);
970 CHKERRQ(ierr); // Check for errors in position determination.
971
973 PetscFunctionReturn(0); // Indicate successful execution of the function.
974}
PetscBool is_function_allowed(const char *functionName)
Checks if a given function is in the allow-list.
Definition logging.c:157
LogLevel get_log_level()
Retrieves the current logging level from the environment variable LOG_LEVEL.
Definition logging.c:39
PetscErrorCode LOG_FACE_DISTANCES(PetscReal *d)
Prints the signed distances to each face of the cell.
Definition logging.c:227
PetscErrorCode GetCellCharacteristicSize(const Cell *cell, PetscReal *cellSize)
Estimates a characteristic length of the cell for threshold scaling.
PetscErrorCode CalculateDistancesToCellFaces(const Cmpnts p, const Cell *cell, PetscReal *d, const PetscReal threshold)
Computes the signed distances from a point to each face of a cubic cell.
PetscErrorCode DeterminePointPosition(PetscReal *d, PetscInt *result)
Classifies a point based on precomputed face distances.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ LocateParticleInGrid()

PetscErrorCode LocateParticleInGrid ( UserCtx user,
Particle particle 
)

Locates the grid cell containing a particle using a walking search.

This function implements a walking search algorithm to find the specific cell (identified by global indices i, j, k) that encloses the particle's physical location (particle->loc).

The search starts from an initial guess cell (either the particle's previously known cell or the corner of the local process domain) and iteratively steps to adjacent cells based on the signed distances from the particle to the faces of the current cell.

Upon successful location (position == 0), the particle's cell field is updated with the found indices (i,j,k), and the corresponding interpolation weights are calculated and stored using the distances (d) relative to the final cell.

Handles particles exactly on cell boundaries by attempting a tie-breaker if the search gets stuck oscillating between adjacent cells due to the boundary condition.

Parameters
[in]userPointer to the UserCtx structure containing grid information (DMDA, coordinates) and domain boundaries.
[in,out]particlePointer to the Particle structure. Its loc field provides the target position. On successful return, its cell and weights fields are updated. On failure, cell is set to {-1, -1, -1} and weights to {0.0, 0.0, 0.0}.
Returns
PetscErrorCode 0 on success. Non-zero error codes may indicate issues during coordinate access, distance calculation, or other internal errors. A return code of 0 does not guarantee the particle was found; check particle->cell[0] >= 0 afterward.
Note
Relies on helper functions like InitializeTraversalParameters, CheckCellWithinLocalGrid, RetrieveCurrentCell, EvaluateParticlePosition, UpdateCellIndicesBasedOnDistances, UpdateParticleWeights, and FinalizeTraversal.
Warning
Ensure particle->loc is set correctly before calling.
The function may fail to find the particle if it lies outside the domain accessible by the current process (including ghost cells) or if MAX_TRAVERSAL steps are exceeded.

◆ FinalizeTraversal()

PetscErrorCode FinalizeTraversal ( UserCtx user,
Particle particle,
PetscInt  traversal_steps,
PetscBool  cell_found,
PetscInt  idx,
PetscInt  idy,
PetscInt  idz 
)

Finalizes the traversal by reporting the results.

This function prints the outcome of the traversal, indicating whether the particle was found within a cell or not, and updates the particle's cell indices accordingly.

Parameters
[in]userPointer to the user-defined context containing grid information.
[out]particlePointer to the Particle structure to update with cell indices.
[in]traversal_stepsThe number of traversal steps taken.
[in]cell_foundFlag indicating whether the particle was found within a cell.
[in]idxThe i-index of the found cell.
[in]idyThe j-index of the found cell.
[in]idzThe k-index of the found cell.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

Definition at line 1100 of file walkingsearch.c.

1101{
1102 PetscFunctionBeginUser;
1103 // Validate input pointers
1104 if (user == NULL || particle == NULL) {
1105 LOG_ALLOW(LOCAL,LOG_ERROR, "One or more input pointers are NULL.\n");
1106 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "One or more input pointers are NULL.");
1107 }
1108
1109 if (cell_found) {
1110 LOG_ALLOW(LOCAL,LOG_INFO, "Particle located in cell (%d, %d, %d) after %d traversal steps.\n",
1111 idx, idy, idz, traversal_steps);
1112 }
1113 else {
1114 LOG_ALLOW(LOCAL,LOG_WARNING, "Particle could not be located within the grid after %d traversal steps.\n", (PetscInt)traversal_steps);
1115 particle->cell[0] = -1;
1116 particle->cell[1] = -1;
1117 particle->cell[2] = -1;
1118 }
1119
1120 LOG_ALLOW(LOCAL, LOG_INFO, "Completed final traversal sync across all ranks.\n");
1121
1122
1124 PetscFunctionReturn(0);
1125}

◆ FindOwnerOfCell()

PetscErrorCode FindOwnerOfCell ( UserCtx user,
PetscInt  i,
PetscInt  j,
PetscInt  k,
PetscMPIInt *  owner_rank 
)

Finds the MPI rank that owns a given global cell index.

This function performs a linear search through the pre-computed decomposition map (user->RankCellInfoMap) to determine which process is responsible for the cell with global indices (i, j, k). It is the definitive method for resolving cell ownership in the "Walk and Handoff" migration algorithm.

If the provided indices are outside the range of any rank (e.g., negative or beyond the global domain), the function will not find an owner and owner_rank will be set to -1.

Parameters
[in]userPointer to the UserCtx structure, which must contain the initialized RankCellInfoMap and num_ranks.
[in]iGlobal i-index of the cell to find.
[in]jGlobal j-index of the cell to find.
[in]kGlobal k-index of the cell to find.
[out]owner_rankPointer to a PetscMPIInt where the resulting owner rank will be stored. It is set to -1 if no owner is found.
Returns
PetscErrorCode 0 on success, or a non-zero PETSc error code on failure.

Definition at line 1152 of file walkingsearch.c.

1153{
1154 PetscErrorCode ierr;
1155 PetscMPIInt size;
1156
1157 PetscFunctionBeginUser;
1158
1160
1161 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size); CHKERRQ(ierr);
1162
1163 // --- 1. Input Validation ---
1164 if (!user || !user->RankCellInfoMap) {
1165 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "UserCtx or RankCellInfoMap is not initialized in FindOwnerOfCell.");
1166 }
1167 if (!owner_rank) {
1168 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Output pointer owner_rank is NULL in FindOwnerOfCell.");
1169 }
1170
1171 // --- 2. Linear Search through the Decomposition Map ---
1172 // Initialize to a "not found" state.
1173 *owner_rank = -1;
1174
1175 // Loop through the map, which contains the ownership info for every rank 'r'.
1176 for (PetscMPIInt r = 0; r < size; ++r) {
1177 const RankCellInfo *info = &user->RankCellInfoMap[r];
1178
1179 // A rank owns a cell if the cell's index is within its start (inclusive)
1180 // and end (exclusive) range for all three dimensions.
1181 if ((i >= info->xs_cell && i < info->xs_cell + info->xm_cell) &&
1182 (j >= info->ys_cell && j < info->ys_cell + info->ym_cell) &&
1183 (k >= info->zs_cell && k < info->zs_cell + info->zm_cell))
1184 {
1185 *owner_rank = r; // We found the owner.
1186 break; // The search is over, exit the loop.
1187 }
1188 }
1189
1190 // --- 3. Logging for Diagnostics ---
1191 // This is extremely useful for debugging particle migration issues.
1192 if (*owner_rank == -1) {
1193 LOG_ALLOW(LOCAL, LOG_DEBUG, "No owner found for global cell (%d, %d, %d). It is likely outside the domain.\n", i, j, k);
1194 } else {
1195 LOG_ALLOW(LOCAL, LOG_DEBUG, "Owner of cell (%d, %d, %d) is Rank %d.\n", i, j, k, *owner_rank);
1196 }
1197
1199 PetscFunctionReturn(0);
1200}
PetscInt ys_cell
Definition variables.h:187
PetscInt xs_cell
Definition variables.h:187
PetscInt zm_cell
Definition variables.h:188
PetscInt zs_cell
Definition variables.h:187
PetscInt xm_cell
Definition variables.h:188
RankCellInfo * RankCellInfoMap
Definition variables.h:799
PetscInt ym_cell
Definition variables.h:188
A lean struct to hold the global cell ownership range for a single MPI rank.
Definition variables.h:186
Here is the caller graph for this function:

◆ LocateParticleOrFindMigrationTarget()

PetscErrorCode LocateParticleOrFindMigrationTarget ( UserCtx user,
Particle particle,
ParticleLocationStatus status_out 
)

Locates a particle's host cell or identifies its migration target using a robust walk search.

This is the core search engine. It starts from a guess cell and walks through the grid. It returns a definitive, actionable status indicating the outcome:

  • ACTIVE_AND_LOCATED: The particle was found in a cell on the current rank.
  • MIGRATING_OUT: The particle was found to belong to another rank. particle->destination_rank is set.
  • LOST: The search failed to find the particle within the global domain.

This function is globally aware and can walk across MPI rank boundaries. It contains robust checks for global domain boundaries and a tie-breaker for numerically "stuck" particles on cell faces.

Parameters
[in]userPointer to the UserCtx containing all grid and domain info.
[in,out]particlePointer to the Particle struct. Its fields are updated based on the search outcome.
[out]status_outThe final, actionable status of the particle after the search.
Returns
PetscErrorCode 0 on success, or a non-zero PETSc error code on failure.

Definition at line 1226 of file walkingsearch.c.

1229{
1230 PetscErrorCode ierr;
1231 PetscMPIInt rank;
1232
1233 // --- Search State Variables ---
1234 PetscInt idx, idy, idz; // Current search cell global indices
1235 PetscInt traversal_steps; // Counter to prevent infinite loops
1236 PetscBool search_concluded = PETSC_FALSE; // Flag to terminate the main while loop
1237
1238 // --- Oscillation/Stuck Loop Detection Variables ---
1239 PetscInt repeatedIndexCount = 0;
1240 PetscInt prevIdx = PETSC_MIN_INT, prevIdy = PETSC_MIN_INT, prevIdz = PETSC_MIN_INT;
1241 PetscInt last_position_result = -999;
1242
1243 PetscFunctionBeginUser;
1245 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
1246
1247 // --- 1. Initialize the Search ---
1248 ierr = InitializeTraversalParameters(user, particle, &idx, &idy, &idz, &traversal_steps); CHKERRQ(ierr);
1249
1250 LOG_ALLOW(LOCAL,LOG_VERBOSE, " The Threshold for considering a particle to be at a face is %.16f.\n",DISTANCE_THRESHOLD);
1251
1252 LOG_ALLOW(LOCAL,LOG_TRACE," [PID %lld]Traversal Initiated at : i = %d, j = %d, k = %d.\n",(long long)particle->PID,idx,idy,idz);
1253
1254 // --- 2. Main Walking Search Loop ---
1255 while (!search_concluded && traversal_steps < MAX_TRAVERSAL) {
1256 traversal_steps++;
1257
1258 // --- 2a. GLOBAL Domain Boundary Check ---
1259 // IMPROVED: Instead of failing immediately when hitting a boundary, clamp the indices
1260 // to the boundary and continue searching. This allows the search to explore other
1261 // directions that might lead to the particle.
1262 PetscBool hit_boundary = PETSC_FALSE;
1263 if (idx < 0) {
1264 idx = 0;
1265 hit_boundary = PETSC_TRUE;
1266 LOG_ALLOW(LOCAL, LOG_DEBUG, "[PID %lld]: Hit global i-min boundary, clamped to idx=0.\n", (long long)particle->PID);
1267 }
1268 if (idx >= (user->IM - 1)) {
1269 idx = user->IM - 2;
1270 hit_boundary = PETSC_TRUE;
1271 LOG_ALLOW(LOCAL, LOG_DEBUG, "[PID %lld]: Hit global i-max boundary, clamped to idx=%d.\n", (long long)particle->PID, idx);
1272 }
1273 if (idy < 0) {
1274 idy = 0;
1275 hit_boundary = PETSC_TRUE;
1276 LOG_ALLOW(LOCAL, LOG_DEBUG, "[PID %lld]: Hit global j-min boundary, clamped to idy=0.\n", (long long)particle->PID);
1277 }
1278 if (idy >= (user->JM - 1)) {
1279 idy = user->JM - 2;
1280 hit_boundary = PETSC_TRUE;
1281 LOG_ALLOW(LOCAL, LOG_DEBUG, "[PID %lld]: Hit global j-max boundary, clamped to idy=%d.\n", (long long)particle->PID, idy);
1282 }
1283 if (idz < 0) {
1284 idz = 0;
1285 hit_boundary = PETSC_TRUE;
1286 LOG_ALLOW(LOCAL, LOG_DEBUG, "[PID %lld]: Hit global k-min boundary, clamped to idz=0.\n", (long long)particle->PID);
1287 }
1288 if (idz >= (user->KM - 1)) {
1289 idz = user->KM - 2;
1290 hit_boundary = PETSC_TRUE;
1291 LOG_ALLOW(LOCAL, LOG_DEBUG, "[PID %lld]: Hit global k-max boundary, clamped to idz=%d.\n", (long long)particle->PID, idz);
1292 }
1293
1294 if (hit_boundary) {
1295 LOG_ALLOW(LOCAL, LOG_INFO, "[PID %lld]: Hit GLOBAL domain boundary. Clamped to boundary cell (%d,%d,%d) and continuing search.\n",
1296 (long long)particle->PID, idx, idy, idz);
1297 }
1298
1299 // --- 2b. LOCAL GHOST REGION CHECK (PREVENTS SEGV) ---
1300 // Before trying to access cell data, check if we have it.
1301 PetscBool is_cell_local;
1302 ierr = CheckCellWithinLocalGrid(user, idx, idy, idz, &is_cell_local); CHKERRQ(ierr);
1303 if (!is_cell_local) {
1304 // We have walked outside the local rank's ghost region.
1305 // This definitively means the particle belongs to another rank.
1306 // Conclude the search here; the current (idx,idy,idz) is the handoff cell.
1307 LOG_ALLOW(LOCAL, LOG_INFO, "[PID %lld]: Walked outside local ghost region to cell (%d,%d,%d). Concluding search for handoff.\n",
1308 (long long)particle->PID, idx, idy, idz);
1309 search_concluded = PETSC_TRUE;
1310 continue; // Skip the rest of the loop; proceed to finalization.
1311 }
1312
1313 // --- 2c. Stuck Loop Detection & Enhanced Tie-Breaker ---
1314 if (idx == prevIdx && idy == prevIdy && idz == prevIdz) {
1315 repeatedIndexCount++;
1316 if (repeatedIndexCount > REPEAT_COUNT_THRESHOLD) {
1317 // Only apply tie-breaker if we are stuck for the right reason (on a boundary)
1318 if (last_position_result >= 1) {
1319 LOG_ALLOW(LOCAL, LOG_WARNING, "[PID %lld]: Stuck on boundary of cell (%d,%d,%d) for %d steps. Applying enhanced tie-breaker.\n",
1320 (long long)particle->PID, idx, idy, idz, repeatedIndexCount);
1321
1322 // Re-evaluate at the stuck cell to get definitive weights
1323 Cell final_cell;
1324 PetscReal final_d[NUM_FACES];
1325 PetscInt final_position; // Dummy variable
1326
1327 ierr = RetrieveCurrentCell(user, idx, idy, idz, &final_cell); CHKERRQ(ierr);
1328 ierr = EvaluateParticlePosition(&final_cell, final_d, particle->loc, &final_position, DISTANCE_THRESHOLD);
1329
1330 if (ierr == 0) { // If evaluation succeeded
1331 ierr = UpdateParticleWeights(final_d, particle); CHKERRQ(ierr);
1332 search_concluded = PETSC_TRUE; // Conclude search, accepting this cell.
1333 } else { // Evaluation failed (e.g., degenerate cell)
1334 LOG_ALLOW(LOCAL, LOG_ERROR, "[PID %lld]: Tie-breaker failed during final evaluation at cell (%d,%d,%d). Search fails.\n",
1335 (long long)particle->PID, idx, idy, idz);
1336 idx = -1; // Invalidate result
1337 search_concluded = PETSC_TRUE;
1338 }
1339 } else { // Stuck for the wrong reason (not on a boundary)
1340 LOG_ALLOW(LOCAL, LOG_ERROR, "[PID %lld]: Search is stuck at cell (%d,%d,%d) but not on a boundary. This indicates a logic error. Failing search.\n",
1341 (long long)particle->PID, idx, idy, idz);
1342 idx = -1; // Invalidate result
1343 search_concluded = PETSC_TRUE;
1344 }
1345 if(search_concluded) continue;
1346 }
1347 } else {
1348 repeatedIndexCount = 0;
1349 }
1350 prevIdx = idx; prevIdy = idy; prevIdz = idz;
1351
1352 // --- 2d. Geometric Evaluation ---
1353 Cell current_cell;
1354 PetscReal distances[NUM_FACES];
1355 PetscInt position_in_cell;
1356
1357 ierr = RetrieveCurrentCell(user, idx, idy, idz, &current_cell); CHKERRQ(ierr);
1358 ierr = EvaluateParticlePosition(&current_cell, distances, particle->loc, &position_in_cell, DISTANCE_THRESHOLD); CHKERRQ(ierr);
1359 last_position_result = position_in_cell;
1360
1361 // --- 2e. Decision Making ---
1362 if (position_in_cell >= 0) { // Particle is INSIDE or ON THE BOUNDARY
1363 search_concluded = PETSC_TRUE;
1364 ierr = UpdateParticleWeights(distances, particle); CHKERRQ(ierr);
1365 } else { // Particle is OUTSIDE
1366 ierr = UpdateCellIndicesBasedOnDistances(distances, &idx, &idy, &idz); CHKERRQ(ierr);
1367 }
1368 }
1369
1370 // --- 3. Finalize and Determine Actionable Status ---
1371 if (idx == -1 || (!search_concluded && traversal_steps >= MAX_TRAVERSAL)) {
1372 if (idx != -1) {
1373 LOG_ALLOW(LOCAL, LOG_ERROR, "[PID %lld]: Search FAILED, exceeded MAX_TRAVERSAL limit of %d.\n",
1374 (long long)particle->PID, MAX_TRAVERSAL);
1375 }
1376 *status_out = LOST;
1377 particle->cell[0] = -1; particle->cell[1] = -1; particle->cell[2] = -1;
1378 } else {
1379 // Search succeeded in finding a candidate cell, now determine its owner.
1380 PetscMPIInt owner_rank;
1381 ierr = FindOwnerOfCell(user, idx, idy, idz, &owner_rank); CHKERRQ(ierr);
1382
1383 LOG_ALLOW(LOCAL,LOG_TRACE," [PID %ld] Owner rank : %d.\n",particle->PID,owner_rank);
1384
1385 // Always update the particle's cell index. It's a good guess for the receiving rank.
1386 particle->cell[0] = idx; particle->cell[1] = idy; particle->cell[2] = idz;
1387
1388 if (owner_rank == rank) {
1389 *status_out = ACTIVE_AND_LOCATED;
1390 } else if (owner_rank != -1) {
1391 // Particle belongs to another rank. Return the direct, actionable status.
1392 *status_out = MIGRATING_OUT;
1393 particle->destination_rank = owner_rank;
1394 } else { // Found a valid index, but no owner in the map.
1395 *status_out = LOST;
1396 particle->cell[0] = -1; particle->cell[1] = -1; particle->cell[2] = -1;
1397 }
1398 }
1399
1400 LOG_ALLOW(LOCAL,LOG_DEBUG,"[Rank %d][PID %ld] Search complete.\n",rank,particle->PID);
1401
1402 // --- 4. Report the Final Outcome ---
1403 ierr = ReportSearchOutcome(particle, *status_out, traversal_steps); CHKERRQ(ierr);
1404
1406 PetscFunctionReturn(0);
1407}
PetscErrorCode UpdateParticleWeights(PetscReal *d, Particle *particle)
Updates a particle's interpolation weights based on distances to cell faces.
@ LOG_TRACE
Very fine-grained tracing information for in-depth debugging.
Definition logging.h:33
@ LOST
Definition variables.h:139
@ ACTIVE_AND_LOCATED
Definition variables.h:137
@ MIGRATING_OUT
Definition variables.h:138
PetscInt KM
Definition variables.h:737
PetscMPIInt destination_rank
Definition variables.h:172
PetscInt JM
Definition variables.h:737
PetscInt IM
Definition variables.h:737
Defines the vertices of a single hexahedral grid cell.
Definition variables.h:160
#define DISTANCE_THRESHOLD
#define MAX_TRAVERSAL
PetscErrorCode EvaluateParticlePosition(const Cell *cell, PetscReal *d, const Cmpnts p, PetscInt *position, const PetscReal threshold)
Determines the spatial relationship of a particle relative to a cubic cell.
PetscErrorCode ReportSearchOutcome(const Particle *particle, ParticleLocationStatus status, PetscInt traversal_steps)
Logs the final outcome of the particle location search.
#define REPEAT_COUNT_THRESHOLD
PetscErrorCode UpdateCellIndicesBasedOnDistances(PetscReal d[NUM_FACES], PetscInt *idx, PetscInt *idy, PetscInt *idz)
Updates the cell indices based on the signed distances to each face.
PetscErrorCode FindOwnerOfCell(UserCtx *user, PetscInt i, PetscInt j, PetscInt k, PetscMPIInt *owner_rank)
Finds the MPI rank that owns a given global cell index.
PetscErrorCode RetrieveCurrentCell(UserCtx *user, PetscInt idx, PetscInt idy, PetscInt idz, Cell *cell)
Retrieves the coordinates of the eight vertices of the current cell.
PetscErrorCode InitializeTraversalParameters(UserCtx *user, Particle *particle, PetscInt *idx, PetscInt *idy, PetscInt *idz, PetscInt *traversal_steps)
Initializes traversal parameters for locating a particle.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ReportSearchOutcome()

PetscErrorCode ReportSearchOutcome ( const Particle particle,
ParticleLocationStatus  status,
PetscInt  traversal_steps 
)

Logs the final outcome of the particle location search.

Definition at line 1415 of file walkingsearch.c.

1418{
1419 PetscFunctionBeginUser;
1421 switch (status) {
1422 case ACTIVE_AND_LOCATED:
1423 LOG_ALLOW(LOCAL, LOG_INFO, "Search SUCCESS [PID %lld]: Located in global cell (%d, %d, %d) after %d steps.\n",
1424 (long long)particle->PID, particle->cell[0], particle->cell[1], particle->cell[2], traversal_steps);
1425 break;
1426 case MIGRATING_OUT:
1427 LOG_ALLOW(LOCAL, LOG_INFO, "Search SUCCESS [PID %lld]: Identified for migration to Rank %d. Handoff cell is (%d, %d, %d) after %d steps.\n",
1428 (long long)particle->PID, particle->destination_rank, particle->cell[0], particle->cell[1], particle->cell[2], traversal_steps);
1429 break;
1430 case LOST:
1431 LOG_ALLOW(LOCAL, LOG_WARNING, "Search FAILED [PID %lld]: Particle is LOST after %d steps.\n",
1432 (long long)particle->PID, traversal_steps);
1433 break;
1434 default:
1435 LOG_ALLOW(LOCAL, LOG_WARNING, "Search ended with unexpected status %d for PID %lld.\n", status, (long long)particle->PID);
1436 break;
1437 }
1439 PetscFunctionReturn(0);
1440}
Here is the caller graph for this function:

◆ UpdateCellIndicesBasedOnDistances()

PetscErrorCode UpdateCellIndicesBasedOnDistances ( PetscReal  d[NUM_FACES],
PetscInt *  idx,
PetscInt *  idy,
PetscInt *  idz 
)

Updates the cell indices based on the signed distances to each face.

This function modifies the cell indices (idx, idy, idz) to move towards the direction where the particle is likely to be located, based on positive distances indicating that the particle is outside in that particular direction.

Parameters
[in]dAn array of six PetscReal values representing the signed distances to each face:
  • d[LEFT]: Left Face
  • d[RIGHT]: Right Face
  • d[BOTTOM]: Bottom Face
  • d[TOP]: Top Face
  • d[FRONT]: Front Face
  • d[BACK]: Back Face
[out]idxPointer to the i-index of the cell to be updated.
[out]idyPointer to the j-index of the cell to be updated.
[out]idzPointer to the k-index of the cell to be updated.
[in]infoDMDALocalInfo structure that holds local & global domain bounds.
Returns
PetscErrorCode Returns 0 on success, non-zero on failure.

Definition at line 999 of file walkingsearch.c.

1000{
1001
1002 PetscFunctionBeginUser;
1004 /*
1005 PetscInt cxm,cxs; // maximum & minimum cell ID in x
1006 PetscInt cym,cys; // maximum & minimum cell ID in y
1007 PetscInt czm,czs; // maximum & minimum cell ID in z
1008
1009 cxs = info->xs; cxm = cxs + info->xm - 2;
1010 cys = info->ys; cym = cys + info->ym - 2;
1011 czs = info->zs; czm = czs + info->zm - 2;
1012 */
1013
1014 // LOG_ALLOW(LOCAL, LOG_DEBUG, "Received d: "
1015 // "d[LEFT=%d]=%.3e, d[RIGHT=%d]=%.3e, d[BOTTOM=%d]=%.3e, "
1016 // "d[TOP=%d]=%.3e, d[FRONT=%d]=%.3e, d[BACK=%d]=%.3e\n",
1017 // LEFT, d[LEFT], RIGHT, d[RIGHT], BOTTOM, d[BOTTOM],
1018 // TOP, d[TOP], FRONT, d[FRONT], BACK, d[BACK]);
1019 // LOG_ALLOW(LOCAL, LOG_DEBUG, "Raw d: "
1020 // "[%.3e, %.3e, %.3e, %.3e, %.3e, %.3e]\n",
1021 // d[0], d[1], d[2], d[3], d[4], d[5]);
1022
1023 // Validate input pointers
1024 if (d == NULL) {
1025 LOG_ALLOW(LOCAL,LOG_ERROR, "'d' is NULL.\n");
1026 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Input array 'd' is NULL.");
1027 }
1028 if (idx == NULL || idy == NULL || idz == NULL) {
1029 LOG_ALLOW(LOCAL,LOG_ERROR, "One or more index pointers are NULL.\n");
1030 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "One or more index pointers are NULL.");
1031 }
1032
1033 // Debug: Print current face distances
1034 LOG_ALLOW(LOCAL,LOG_DEBUG, "Current Face Distances:\n");
1036
1037 // Update k-direction based on FRONT and BACK distances
1038 if (d[FRONT] < 0.0) {
1039 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Condition met: d[FRONT] < 0.0, incrementing idz.\n");
1040 (*idz) += 1;
1041 }
1042 else if(d[BACK] < 0.0){
1043 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Condition met: d[BACK] < 0.0, decrementing idz.\n");
1044 (*idz) -= 1;
1045 }
1046
1047 // Update i-direction based on LEFT and RIGHT distances
1048 if (d[LEFT] < 0.0) {
1049 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Condition met: d[LEFT] < 0.0, decrementing idx.\n");
1050 (*idx) -= 1;
1051 }
1052 else if (d[RIGHT] < 0.0) {
1053 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Condition met: d[RIGHT] < 0.0, incrementing idx.\n");
1054 (*idx) += 1;
1055 }
1056
1057 // Update j-direction based on BOTTOM and TOP distances
1058 if (d[BOTTOM] < 0.0) {
1059 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Condition met: d[BOTTOM] < 0.0, decrementing idy.\n");
1060 (*idy) -= 1;
1061 }
1062 else if (d[TOP] < 0.0) {
1063 LOG_ALLOW(LOCAL,LOG_VERBOSE, "Condition met: d[TOP] < 0.0, incrementing idy.\n");
1064 (*idy) += 1;
1065 }
1066
1067 /*
1068 // The 'cell' corners you can reference go from [xs .. xs+xm-1], but
1069 // to form a valid cell in x, you need (idx+1) in range, so max is (xs+xm-2).
1070 *idx = PetscMax(cxs, PetscMin(*idx, cxm));
1071 *idy = PetscMax(cys, PetscMin(*idy, cym));
1072 *idz = PetscMax(czs, PetscMin(*idz, czm));
1073 */
1074
1075 LOG_ALLOW(LOCAL,LOG_DEBUG, "Updated Indices - idx, idy, idz: %d, %d, %d\n", *idx, *idy, *idz);
1076
1078 PetscFunctionReturn(0); // Indicate successful execution
1079}
Here is the call graph for this function:
Here is the caller graph for this function: