grid.h File Reference

Public interface for grid, solver, and metric setup routines. More...

#include "variables.h"
#include "logging.h"
#include "io.h"
#include "setup.h"

Include dependency graph for grid.h:

This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Functions
PetscErrorCode	DefineAllGridDimensions (SimCtx *simCtx)
	Orchestrates the parsing and setting of grid dimensions for all blocks.
PetscErrorCode	InitializeAllGridDMs (SimCtx *simCtx)
	Orchestrates the creation of DMDA objects for every block and multigrid level.
PetscErrorCode	AssignAllGridCoordinates (SimCtx *simCtx)
	Orchestrates the assignment of physical coordinates to all DMDA objects.
PetscErrorCode	ComputeLocalBoundingBox (UserCtx *user, BoundingBox *localBBox)
	Computes the local bounding box of the grid on the current process.
PetscErrorCode	GatherAllBoundingBoxes (UserCtx *user, BoundingBox **allBBoxes)
	Gathers local bounding boxes from all MPI processes to rank 0.
PetscErrorCode	BroadcastAllBoundingBoxes (UserCtx *user, BoundingBox **bboxlist)
	Broadcasts the bounding box information collected on rank 0 to all other ranks.

Detailed Description

Public interface for grid, solver, and metric setup routines.

Definition in file grid.h.

Function Documentation

DefineAllGridDimensions()

PetscErrorCode DefineAllGridDimensions

(

SimCtx *

simCtx

)

Orchestrates the parsing and setting of grid dimensions for all blocks.

This function serves as the high-level entry point for defining the geometric properties of each grid block in the simulation. It iterates through every block defined by simCtx->block_number.

For each block, it performs two key actions:

1. It explicitly sets the block's index (_this) in the corresponding UserCtx struct for the finest multigrid level. This makes the context "self-aware".
2. It calls a helper function (ParseAndSetGridInputs) to handle the detailed work of parsing options or files to populate the rest of the geometric properties for that specific block (e.g., IM, Min_X, rx).

Parameters

simCtx

The master SimCtx, which contains the number of blocks and the UserCtx hierarchy to be configured.

Returns

PetscErrorCode 0 on success, or a PETSc error code on failure.

Definition at line 57 of file grid.c.

{
    PetscErrorCode ierr;
    PetscInt       nblk = simCtx->block_number;
    UserCtx        *fine_users;
 
    PetscFunctionBeginUser;
 
    if (simCtx->usermg.mglevels == 0) {
        SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "MG levels not set. Cannot get finest_users.");
    }
    // Get the UserCtx array for the finest grid level
    fine_users = simCtx->usermg.mgctx[simCtx->usermg.mglevels - 1].user;
 
    LOG_ALLOW(GLOBAL, LOG_INFO, "Defining grid dimensions for %d blocks...\n", nblk);
 
    // Loop over each block to configure its grid dimensions and geometry.
    for (PetscInt bi = 0; bi < nblk; bi++) {
        LOG_ALLOW(LOCAL, LOG_DEBUG, "Rank %d: --- Configuring Geometry for Block %d ---\n", simCtx->rank, bi);
 
        // Before calling any helpers, set the block index in the context.
        // This makes the UserCtx self-aware of which block it represents.
        fine_users[bi]._this = bi;
 
        // Call the helper function for this specific block. It can now derive
        // all necessary information from the UserCtx pointer it receives.
        ierr = ParseAndSetGridInputs(&fine_users[bi]); CHKERRQ(ierr);
    }
 
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

InitializeAllGridDMs()

PetscErrorCode InitializeAllGridDMs

(

SimCtx *

simCtx

)

Orchestrates the creation of DMDA objects for every block and multigrid level.

This function systematically builds the entire DMDA hierarchy. It first calculates the dimensions (IM, JM, KM) for all coarse grids based on the finest grid's dimensions and the semi-coarsening flags. It then iterates from the coarsest to the finest level, calling a powerful helper function (InitializeSingleGridDM) to create the DMs for each block, ensuring that finer grids are properly aligned with their coarser parents for multigrid efficiency.

Parameters

simCtx

The master SimCtx, containing the configured UserCtx hierarchy.

Returns

PetscErrorCode 0 on success, or a PETSc error code on failure.

Definition at line 230 of file grid.c.

{
    PetscErrorCode ierr;
    UserMG         *usermg = &simCtx->usermg;
    MGCtx          *mgctx = usermg->mgctx;
    PetscInt       nblk = simCtx->block_number;
 
    PetscFunctionBeginUser;
    LOG_ALLOW(GLOBAL, LOG_INFO, "Creating DMDA objects for all levels and blocks...\n");
 
    // --- Part 1: Calculate Coarse Grid Dimensions & VALIDATE ---
    LOG_ALLOW(GLOBAL, LOG_DEBUG, "Calculating and validating coarse grid dimensions...\n");
    for (PetscInt level = usermg->mglevels - 2; level >= 0; level--) {
        for (PetscInt bi = 0; bi < nblk; bi++) {
            UserCtx *user_coarse = &mgctx[level].user[bi];
            UserCtx *user_fine   = &mgctx[level + 1].user[bi];
 
            user_coarse->IM = user_fine->isc ? user_fine->IM : (user_fine->IM + 1) / 2;
            user_coarse->JM = user_fine->jsc ? user_fine->JM : (user_fine->JM + 1) / 2;
            user_coarse->KM = user_fine->ksc ? user_fine->KM : (user_fine->KM + 1) / 2;
 
            LOG_ALLOW_SYNC(LOCAL, LOG_DEBUG, "Rank %d: Block %d, Level %d dims calculated: %d x %d x %d\n",
                      simCtx->rank, bi, level, user_coarse->IM, user_coarse->JM, user_coarse->KM);
 
            // Validation check from legacy MGDACreate to ensure coarsening is possible
            PetscInt check_i = user_coarse->IM * (2 - user_coarse->isc) - (user_fine->IM + 1 - user_coarse->isc);
            PetscInt check_j = user_coarse->JM * (2 - user_coarse->jsc) - (user_fine->JM + 1 - user_coarse->jsc);
            PetscInt check_k = user_coarse->KM * (2 - user_coarse->ksc) - (user_fine->KM + 1 - user_coarse->ksc);
 
            if (check_i + check_j + check_k != 0) {
          //      SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG,
          //           "Grid at level %d, block %d cannot be coarsened from %dx%dx%d to %dx%dx%d with the given semi-coarsening flags. Check grid dimensions.",
              //          level, bi, user_fine->IM, user_fine->JM, user_fine->KM, user_coarse->IM, user_coarse->JM, user_coarse->KM);
          LOG(GLOBAL,LOG_WARNING,"WARNING: Grid at level %d, block %d can't be consistently coarsened further.\n", level, bi);
            }
        }
    }
 
    // --- Part 2: Create DMs from Coarse to Fine for each Block ---
    for (PetscInt bi = 0; bi < nblk; bi++) {
        LOG_ALLOW_SYNC(GLOBAL, LOG_DEBUG, "--- Creating DMs for Block %d ---\n", bi);
 
        // Create the coarsest level DM first (passing NULL for the coarse_user)
        ierr = InitializeSingleGridDM(&mgctx[0].user[bi], NULL); CHKERRQ(ierr);
 
        // Create finer level DMs, passing the next-coarser context for alignment
        for (PetscInt level = 1; level < usermg->mglevels; level++) {
            ierr = InitializeSingleGridDM(&mgctx[level].user[bi], &mgctx[level-1].user[bi]); CHKERRQ(ierr);
        }
    }
 
    // --- Optional: View the finest DM for debugging verification ---
    if (get_log_level() >= LOG_DEBUG) {
        LOG_ALLOW_SYNC(GLOBAL, LOG_INFO, "--- Viewing Finest DMDA (Level %d, Block 0) ---\n", usermg->mglevels - 1);
        ierr = DMView(mgctx[usermg->mglevels - 1].user[0].da, PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
    }
 
    LOG_ALLOW(GLOBAL, LOG_INFO, "DMDA object creation complete.\n");
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

AssignAllGridCoordinates()

PetscErrorCode AssignAllGridCoordinates

(

SimCtx *

simCtx

)

Orchestrates the assignment of physical coordinates to all DMDA objects.

This function manages the entire process of populating the coordinate vectors for every DMDA across all multigrid levels and blocks. It follows a two-part strategy that is essential for multigrid methods:

1. Populate Finest Level: It first loops through each block and calls a helper (SetFinestLevelCoordinates) to set the physical coordinates for the highest-resolution grid (the finest multigrid level).
2. Restrict to Coarser Levels: It then iterates downwards from the finest level, calling a helper (RestrictCoordinates) to copy the coordinate values from the fine grid nodes to their corresponding parent nodes on the coarser grids. This ensures all levels represent the exact same geometry.

Parameters

simCtx

The master SimCtx, containing the configured UserCtx hierarchy.

Returns

PetscErrorCode 0 on success, or a PETSc error code on failure.

Definition at line 317 of file grid.c.

{
    PetscErrorCode ierr;
    UserMG         *usermg = &simCtx->usermg;
    PetscInt       nblk = simCtx->block_number;
 
    PetscFunctionBeginUser;
    LOG_ALLOW(GLOBAL, LOG_INFO, "Assigning physical coordinates to all grid DMs...\n");
 
    // --- Part 1: Populate the Finest Grid Level ---
    LOG_ALLOW(GLOBAL, LOG_DEBUG, "Setting coordinates for the finest grid level (%d)...\n", usermg->mglevels - 1);
    for (PetscInt bi = 0; bi < nblk; bi++) {
        UserCtx *fine_user = &usermg->mgctx[usermg->mglevels - 1].user[bi];
        ierr = SetFinestLevelCoordinates(fine_user); CHKERRQ(ierr);
    }
    LOG_ALLOW(GLOBAL, LOG_DEBUG, "Finest level coordinates have been set for all blocks.\n");
 
 
    // --- Part 2: Restrict Coordinates to Coarser Levels ---
    LOG_ALLOW(GLOBAL, LOG_DEBUG, "Restricting coordinates to coarser grid levels...\n");
    for (PetscInt level = usermg->mglevels - 2; level >= 0; level--) {
        for (PetscInt bi = 0; bi < nblk; bi++) {
            UserCtx *coarse_user = &usermg->mgctx[level].user[bi];
            UserCtx *fine_user   = &usermg->mgctx[level + 1].user[bi];
            ierr = RestrictCoordinates(coarse_user, fine_user); CHKERRQ(ierr);
        }
    }
 
    LOG_ALLOW(GLOBAL, LOG_INFO, "Physical coordinates assigned to all grid levels and blocks.\n");
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

ComputeLocalBoundingBox()

PetscErrorCode ComputeLocalBoundingBox	(	UserCtx *	user,
		BoundingBox *	localBBox
	)

Computes the local bounding box of the grid on the current process.

This function calculates the minimum and maximum coordinates of the local grid points owned by the current MPI process and stores the computed bounding box in the provided structure.

Parameters

[in]	user	Pointer to the user-defined context containing grid information.
[out]	localBBox	Pointer to the BoundingBox structure to store the computed bounding box.

Returns

PetscErrorCode Returns 0 on success, non-zero on failure.

This function calculates the minimum and maximum coordinates (x, y, z) of the local grid points owned by the current MPI process. It iterates over the local portion of the grid, examines each grid point's coordinates, and updates the minimum and maximum values accordingly.

The computed bounding box is stored in the provided localBBox structure, and the user->bbox field is also updated with this bounding box for consistency within the user context.

Parameters

[in]	user	Pointer to the user-defined context containing grid information. This context must be properly initialized before calling this function.
[out]	localBBox	Pointer to the BoundingBox structure where the computed local bounding box will be stored. The structure should be allocated by the caller.

Returns

PetscErrorCode Returns 0 on success, non-zero on failure.

Definition at line 702 of file grid.c.

{
    PetscErrorCode ierr;
    PetscInt i, j, k;
    PetscMPIInt rank;
    PetscInt xs, ys, zs, xe, ye, ze;
    DMDALocalInfo info;
    Vec coordinates;
    Cmpnts ***coordArray;
    Cmpnts minCoords, maxCoords;
 
    // Start of function execution
    LOG_ALLOW(GLOBAL, LOG_INFO, "Entering the function.\n");
 
    // Validate input Pointers
    if (!user) {
        LOG_ALLOW(LOCAL, LOG_ERROR, "Input 'user' Pointer is NULL.\n");
        return PETSC_ERR_ARG_NULL;
    }
    if (!localBBox) {
        LOG_ALLOW(LOCAL, LOG_ERROR, "Output 'localBBox' Pointer is NULL.\n");
        return PETSC_ERR_ARG_NULL;
    }
 
    // Get MPI rank
    ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr);
 
    // Get the local coordinates vector from the DMDA
    ierr = DMGetCoordinatesLocal(user->da, &coordinates);
    if (ierr) {
        LOG_ALLOW(LOCAL, LOG_ERROR, "Error getting local coordinates vector.\n");
        return ierr;
    }
 
    if (!coordinates) {
        LOG_ALLOW(LOCAL, LOG_ERROR, "Coordinates vector is NULL.\n");
        return PETSC_ERR_ARG_NULL;
    }
 
    // Access the coordinate array for reading
    ierr = DMDAVecGetArrayRead(user->fda, coordinates, &coordArray);
    if (ierr) {
        LOG_ALLOW(LOCAL, LOG_ERROR, "Error accessing coordinate array.\n");
        return ierr;
    }
 
    // Get the local grid information (indices and sizes)
    ierr = DMDAGetLocalInfo(user->da, &info);
    if (ierr) {
        LOG_ALLOW(LOCAL, LOG_ERROR, "Error getting DMDA local info.\n");
        return ierr;
    }
 
    
    xs = info.gxs; xe = xs + info.gxm;
    ys = info.gys; ye = ys + info.gym;
    zs = info.gzs; ze = zs + info.gzm;
    
    /*
    xs = info.xs; xe = xs + info.xm;
    ys = info.ys; ye = ys + info.ym;
    zs = info.zs; ze = zs + info.zm;
    */
    
    // Initialize min and max coordinates with extreme values
    minCoords.x = minCoords.y = minCoords.z = PETSC_MAX_REAL;
    maxCoords.x = maxCoords.y = maxCoords.z = PETSC_MIN_REAL;
 
    LOG_ALLOW(LOCAL, LOG_DEBUG, "[Rank %d] Grid indices (Including Ghosts): xs=%d, xe=%d, ys=%d, ye=%d, zs=%d, ze=%d.\n",rank, xs, xe, ys, ye, zs, ze);
 
    // Iterate over the local grid to find min and max coordinates
    for (k = zs; k < ze; k++) {
        for (j = ys; j < ye; j++) {
            for (i = xs; i < xe; i++) {
                Cmpnts coord = coordArray[k][j][i];
 
                // Update min and max coordinates
                if (coord.x < minCoords.x) minCoords.x = coord.x;
                if (coord.y < minCoords.y) minCoords.y = coord.y;
                if (coord.z < minCoords.z) minCoords.z = coord.z;
 
                if (coord.x > maxCoords.x) maxCoords.x = coord.x;
                if (coord.y > maxCoords.y) maxCoords.y = coord.y;
                if (coord.z > maxCoords.z) maxCoords.z = coord.z;
            }
        }
    }
 
 
    // Add tolerance to bboxes.
    minCoords.x =  minCoords.x - BBOX_TOLERANCE;
    minCoords.y =  minCoords.y - BBOX_TOLERANCE;
    minCoords.z =  minCoords.z - BBOX_TOLERANCE;
 
    maxCoords.x =  maxCoords.x + BBOX_TOLERANCE;
    maxCoords.y =  maxCoords.y + BBOX_TOLERANCE;
    maxCoords.z =  maxCoords.z + BBOX_TOLERANCE;
 
    LOG_ALLOW(LOCAL,LOG_INFO," Tolerance added to the limits: %.8e .\n",(PetscReal)BBOX_TOLERANCE);
       
    // Log the computed min and max coordinates
     LOG_ALLOW(LOCAL, LOG_INFO,"[Rank %d]minCoords=(%.6f, %.6f, %.6f), maxCoords=(%.6f, %.6f, %.6f).\n",rank,minCoords.x, minCoords.y, minCoords.z, maxCoords.x, maxCoords.y, maxCoords.z);
 
 
    
    // Restore the coordinate array
    ierr = DMDAVecRestoreArrayRead(user->fda, coordinates, &coordArray);
    if (ierr) {
        LOG_ALLOW(LOCAL, LOG_ERROR, "Error restoring coordinate array.\n");
        return ierr;
    }
 
    // Set the local bounding box
    localBBox->min_coords = minCoords;
    localBBox->max_coords = maxCoords;
 
    // Update the bounding box inside the UserCtx for consistency
    user->bbox = *localBBox;
 
    LOG_ALLOW(GLOBAL, LOG_INFO, "Exiting the function successfully.\n");
    return 0;
}

Here is the caller graph for this function:

GatherAllBoundingBoxes()

PetscErrorCode GatherAllBoundingBoxes	(	UserCtx *	user,
		BoundingBox **	allBBoxes
	)

Gathers local bounding boxes from all MPI processes to rank 0.

This function computes the local bounding box on each process, then collects all local bounding boxes on the root process (rank 0) using MPI. The result is stored in an array of BoundingBox structures on rank 0.

Parameters

[in]	user	Pointer to the user-defined context containing grid information.
[out]	allBBoxes	Pointer to a pointer where the array of gathered bounding boxes will be stored on rank 0. The caller on rank 0 must free this array.

Returns

PetscErrorCode Returns 0 on success, non-zero on failure.

Each rank computes its local bounding box, then all ranks participate in an MPI_Gather to send their BoundingBox to rank 0. Rank 0 allocates the result array and returns it via allBBoxes.

Parameters

[in]	user	Pointer to UserCtx (must be non-NULL).
[out]	allBBoxes	On rank 0, receives malloc’d array of size size. On other ranks, set to NULL.

Returns

PetscErrorCode

Definition at line 837 of file grid.c.

{
  PetscErrorCode ierr;
  PetscMPIInt    rank, size;
  BoundingBox    *bboxArray = NULL;
  BoundingBox     localBBox;
 
  /* Validate */
  if (!user || !allBBoxes) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
                                    "GatherAllBoundingBoxes: NULL pointer");
 
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRMPI(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size); CHKERRMPI(ierr);
 
  /* Compute local bbox */
  ierr = ComputeLocalBoundingBox(user, &localBBox); CHKERRQ(ierr);
 
  /* Ensure everyone is synchronized before the gather */
  MPI_Barrier(PETSC_COMM_WORLD);
  LOG_ALLOW_SYNC(GLOBAL, LOG_INFO,
    "Rank %d: about to MPI_Gather(localBBox)\n", rank);
 
  /* Allocate on root */
  if (rank == 0) {
    bboxArray = (BoundingBox*)malloc(size * sizeof(BoundingBox));
    if (!bboxArray) SETERRABORT(PETSC_COMM_WORLD, PETSC_ERR_MEM,
                                "GatherAllBoundingBoxes: malloc failed");
  }
 
  /* Collective: every rank must call */
  ierr = MPI_Gather(&localBBox, sizeof(BoundingBox), MPI_BYTE,
                    bboxArray, sizeof(BoundingBox), MPI_BYTE,
                    0, PETSC_COMM_WORLD);
  CHKERRMPI(ierr);
 
  MPI_Barrier(PETSC_COMM_WORLD);
  LOG_ALLOW_SYNC(GLOBAL, LOG_INFO,
    "Rank %d: completed MPI_Gather(localBBox)\n", rank);
 
  /* Return result */
  if (rank == 0) {
    *allBBoxes = bboxArray;
  } else {
    *allBBoxes = NULL;
  }
 
  return 0;
}

Here is the call graph for this function:

Here is the caller graph for this function:

BroadcastAllBoundingBoxes()

PetscErrorCode BroadcastAllBoundingBoxes	(	UserCtx *	user,
		BoundingBox **	bboxlist
	)

Broadcasts the bounding box information collected on rank 0 to all other ranks.

This function assumes that GatherAllBoundingBoxes() was previously called, so bboxlist is allocated and populated on rank 0. All other ranks will allocate memory for bboxlist, and this function will use MPI_Bcast to distribute the bounding box data to them.

Parameters

[in]	user	Pointer to the UserCtx structure. (Currently unused in this function, but kept for consistency.)
[in,out]	bboxlist	Pointer to the array of BoundingBoxes. On rank 0, this should point to a valid array of size 'size' (where size is the number of MPI ranks). On non-root ranks, this function will allocate memory for bboxlist.

Returns

PetscErrorCode Returns 0 on success, non-zero on MPI or PETSc-related errors.

Broadcasts the bounding box information collected on rank 0 to all other ranks.

After GatherAllBoundingBoxes, rank 0 has an array of size boxes. This routine makes sure every rank ends up with its own malloc’d copy.

Parameters

[in]	user	Pointer to UserCtx (unused here, but kept for signature).
[in,out]	bboxlist	On entry: rank 0’s array; on exit: every rank’s array.

Returns

PetscErrorCode

Definition at line 896 of file grid.c.

{
  PetscErrorCode ierr;
  PetscMPIInt    rank, size;
 
  if (!bboxlist) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
                          "BroadcastAllBoundingBoxes: NULL pointer");
 
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRMPI(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size); CHKERRMPI(ierr);
 
  /* Non-root ranks must allocate before the Bcast */
  if (rank != 0) {
    *bboxlist = (BoundingBox*)malloc(size * sizeof(BoundingBox));
    if (!*bboxlist) SETERRABORT(PETSC_COMM_WORLD, PETSC_ERR_MEM,
                                "BroadcastAllBoundingBoxes: malloc failed");
  }
 
  MPI_Barrier(PETSC_COMM_WORLD);
  LOG_ALLOW_SYNC(GLOBAL, LOG_INFO,
    "Rank %d: about to MPI_Bcast(%d boxes)\n", rank, size);
 
  /* Collective: every rank must call */
  ierr = MPI_Bcast(*bboxlist, size * sizeof(BoundingBox), MPI_BYTE,
                   0, PETSC_COMM_WORLD);
  CHKERRMPI(ierr);
 
  MPI_Barrier(PETSC_COMM_WORLD);
  LOG_ALLOW_SYNC(GLOBAL, LOG_INFO,
    "Rank %d: completed MPI_Bcast(%d boxes)\n", rank, size);
 
  return 0;
}

Here is the caller graph for this function: