Skip to content

Modifying field values for neighboring cells inside a postprocessor

OpenLB – Open Source Lattice Boltzmann Code Forums on OpenLB General Topics Modifying field values for neighboring cells inside a postprocessor

Viewing 15 posts - 31 through 45 (of 54 total)
  • Author
    Posts
  • #9322
    Danial.Khazaeipoul
    Participant

    Thank you, Mathias, for the explanation. The code I’ve implemented works well in a single partition setup where localCellID == globalCellID. However, in configurations with multiple partitions or blocks, localCellID != globalCellID, leading to non-unique IDs.

    For local cell IDs, which are of type CellID or std::uint32_t, there’s an advantage in CCL algorithms as these IDs are unique and remain consistent throughout the simulation. Because of this, I’m looking for a way to generate global cell IDs with similar characteristics. For example, I am thinking about a method like the following which creates a unique global CellID for each cell:

    std::uint32_t globalCellID = blockID × cellsPerBlock + localCellID

    #9326
    Adrian
    Keymaster

    Yes, the cell IDs are block-local. Computing a global ID in the way you posted will work if you use the maximum global block size for the multiplier.

    • This reply was modified 4 months, 1 week ago by Adrian.
    #9333
    Danial.Khazaeipoul
    Participant

    Just to confirm my understanding, by

    if you use the maximum global block size

    Are you referring to the fact that cuboids may not have an equal number of associated cells, and that using the maximum block size can prevent overlapping global cell IDs?

    #9334
    Adrian
    Keymaster

    Yes, the cuboids commonly do not have exactly the same number of cells. You can also see this in the cuboid geometry output.

    #9362
    Danial.Khazaeipoul
    Participant

    Using the cell interface, is there any method that can be used to identify the global cuboid number to which the current cell belongs?

    #9363
    Adrian
    Keymaster

    No, this is not intended as all individual threads of a kernel always only process cells of the same block (all LBM processing is split into per block operators). You should be able to handle the cuboid number / cell ID mapping at the same place where you perform the inter block communication of your algorithm.

    However, in principle nothing prevents you from passing the cuboid number as a block-specific parameter to your operator (with application scope PerCellWithParameters)

    #9474
    Danial.Khazaeipoul
    Participant

    Hello Adrian,

    Is my understanding correct regarding the use of a block-specific parameter passed to an operator? Suppose I have an offset vector, synchronized across all processors, which records the number of cells in each preceding block, shown below.

    struct CELL_COUNT : public descriptors::TYPED_FIELD_BASE<std::size_t,1> { };
    std::vector<std::size_t> offset(_sLattice.getCuboidGeometry().getNc(), 0);

    The size of the offset vector matches the number of global cuboids. For a domain with 4 blocks, the resulting offset vector might look like {0, 431433, 855297, 1279161}, as an example. Given this setup, will the following method yield the expected results for the CELL_COUNT, as shown in the image?

    for (int iC = 0; iC < _sLattice.getLoadBalancer().size(); ++iC)
    {
      auto& block = _sLattice.getBlock(iC);
      block.template setParameter<FreeSurface::CELL_COUNT>(offset[_sLattice.getLoadBalancer().glob(iC)]);
    }

    Please note that the image is only for demonstration purposes. When an operator processes a cell, will calling the parameter use the correct value of CELL_COUNT based on the block in which the cell resides?

    #9478
    Adrian
    Keymaster

    Yes, this should work as you convert the local cuboid index iC to the global one for accessing the offset array.

    #9541
    Danial.Khazaeipoul
    Participant

    Dear Adrian,

    I believe I’ve identified the issue with the CCL implementation when MPI is enabled. Thank you for your insights so far in working toward a robust solution. Since you’ve already reviewed the code, you’re familiar with how cell IDs are used to assign bubble IDs and to navigate between cells in each kernel to access their bubble IDs.

    I’ll do my best to clearly explain the problem. The CCL algorithm assumes that cell IDs follow a consistent order in each direction within a block, which works fine for a single block. However, when MPI is enabled and the domain is divided into multiple blocks, this consistency breaks down for neighboring cells located in different blocks. As a result, the bubble ID assignment becomes incorrect for blocks other than Block 0. In the worst-case scenario, this can lead to illegal memory access, especially since blocks may have an unequal number of cells, see the below picture.

    When the orange cell in Block #2 is being processed, it accesses the bubble ID of the red cell, as the red cell is a valid neighbor. However, since the red cell belongs to a different block, its cell ID does not follow the block-local order of cell IDs in Block #2. When the orange cell adopts the bubble ID of the red cell, it may incorrectly reference a cell within Block #2, as indicated by the blue arrow. This can result in either assigning an incorrect bubble ID to the orange cell or, in the worst-case scenario, causing illegal memory access. This happens because Block #2 may not contain a cell with the cell ID from Block #0.

    Is it possible to defer the automatic communication between blocks for a specific global data structure until it is explicitly triggered when needed?

    Regards,
    Danial

    #9555
    Adrian
    Keymaster

    Sure, the overlap communication only happens at discrete and configurable stages. i.e. you can move all communication related to the CCL implementation into a distinct (or more than one) communicator stage fully within your control. See e.g. User Guide Section 2.2.3.

    You can create a new communicator stage for CCL simply by requesting it:

    
    namespace stage { struct CCL { }; }
    
    {
      auto& comm = sLattice.getCommunicator(stage::CCL{});
      comm.template requestField<descriptors::POPULATION>();
      comm.requestOverlap(1);
      comm.exchangeRequests();
    }
    

    this stage will only be triggered if you explicitly call sLattice.getCommunicator(stage::CCL{}).communicate(). It is also possible to request more than one field or other subsets of the overlap than just a constant width if necessary.

    #9604
    Danial.Khazaeipoul
    Participant

    Hello Adrian,

    Can you confirm that the “apply” operators used in the post-processors operate only on the internal nodes? As I understand it, these operators cannot act directly on overlap nodes. The values of overlap nodes are updated only when communication is triggered using the communicate() method. However, it is still possible to access the values of an overlap node, such as when an internal node has an overlap node as its neighbor.

    Is it possible to determine whether a neighbor of an internal node is an overlap node? I have found a temporary workaround to address my issue, as shown in the picture below. However, a cleaner and safer approach would be preferable as I need to avoid using the overlap nodes values until ids are finalized.

    After the temporary fix is applied, as you can see, the connected regions on each block have their own unique ID.

    #9605
    Adrian
    Keymaster

    Yes, by default all post processors are only applied to the non-overlap nodes and attempts to explicitly set them on overlap nodes are filtered out. In the next release this will change as I also needed more control on this for my FSI implementations. Right now the most straight forward way would be to mark the overlap cells with a dedicated field.

    If you need more control on the post processor placement right now I can probably also extract you a patch.

    #9606
    Danial.Khazaeipoul
    Participant

    Thank you for the assistance Adrian.

    Right now the most straight forward way would be to mark the overlap cells with a dedicated field.

    Are you suggesting that I define a new dedicated field, where for instance overlap cells are assigned a value of “1” and all other cells are assigned a value of “0”? This field could then be used to distinguish overlap cells from regular cells.

    If you need more control on the post processor placement right now I can probably also extract you a patch.

    That would be greatly appreciated, as it would help me gain a better understanding of the code. Additionally, it could provide me with more control over the process.

    #9619
    Adrian
    Keymaster

    Are you suggesting that I define a new dedicated field, where for instance overlap cells are assigned a value of “1” and all other cells are assigned a value of “0”? This field could then be used to distinguish overlap cells from regular cells.

    Yes, exactly (this is of course quite ugly but it will work).

    If you want to add post processors to the overlap you can drop the isPadding check in line 618 of core/blockLattice.h and add them per-cell.

    #9620
    Danial.Khazaeipoul
    Participant

    Yes, exactly (this is of course quite ugly but it will work).

    FYI, I implemented a cleaner method to replace this and now one can check whether a cell is a padding or not directly from the Cell interface. However, this required me to modify the src/core files to introduce a new function inside the BlockStructureD class that retrieves a local LatticR<D> from a linearized CellID.

    FYI, I implemented a more streamlined method to replace this, allowing direct checks on whether a cell is padding via the Cell interface. To achieve this, I updated the src/core files to include a new function in the BlockStructureD class, which retrieves a local LatticeR<D> from a linearized CellID.

      /// Get D-dimensional LatticeR from 1D cell ID
      LatticeR<D> getLatticeR(CellID iCell) const
      {
        if constexpr (D == 3) {
          std::int32_t iX = iCell / _projection[0];
          std::int32_t remainder = iCell % _projection[0];
          std::int32_t iY = remainder / _projection[1];
          std::int32_t iZ = remainder % _projection[1];
          return LatticeR<D>(iX - _padding, iY - _padding, iZ - _padding);
        }
        else {
          std::int32_t iX = iCell / _projection[0];
          std::int32_t iY = iCell % _projection[0];
          return LatticeR<D>(iX - _padding, iY - _padding);
        }
      };
      /// Return whether iCell is valid
      bool isPadding(CellID iCell) const
      {
        return isPadding(getLatticeR(iCell));
      };

    For more details, please refer to the “feature-CCL” branch in my repository, specifically commits 2ce2d4fc and b28e5305.

Viewing 15 posts - 31 through 45 (of 54 total)
  • You must be logged in to reply to this topic.