Skip to content

Commit

Permalink
Allowing the parallel for a more general 6*(2n)^2 blocks (#176)
Browse files Browse the repository at this point in the history
Allowing a more general 6*(2n)^2 blocks parallel
Currently still a backdoor solution:

Need to add #define USE_NBLOCKS and NBLOCKS=600 in configure.hpp
Change configure.hpp before making
Originally allowed configurations are 6*2^(2n), not changed if the backdoor change is not taken
  • Loading branch information
cshsgy authored Aug 29, 2024
1 parent 7024e25 commit 7d8e2b5
Show file tree
Hide file tree
Showing 5 changed files with 87 additions and 68 deletions.
27 changes: 25 additions & 2 deletions src/exo3/cs_find_block_id.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -9,8 +9,9 @@
#include "cubed_sphere.hpp"

int CubedSphere::FindBlockID(LogicalLocation const& loc) {
int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);

int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);

// Determine the block number
int block_id;
Expand All @@ -27,6 +28,7 @@ int CubedSphere::FindBlockID(LogicalLocation const& loc) {
std::stringstream msg;
msg << "Error: something wrong, check the geometry setup of the "
"cubed sphere. \n";
msg << "lv2_lx2: " << lv2_lx2 << " lv2_lx3: " << lv2_lx3 << "bound_lim: " << bound_lim << "loc.lx2: " << loc.lx2 << "loc.lx3: " << loc.lx3 << std::endl;
msg << "----------------------------------" << std::endl;
ATHENA_ERROR(msg);
}
Expand Down Expand Up @@ -73,3 +75,24 @@ int CubedSphere::FindBlockID(LogicalLocation const& loc) {

return block_id;
}

void CubedSphere::GetLocalIndex(int *lv2_lx2, int *lv2_lx3, int *local_lx2, int *local_lx3, int *bound_lim, LogicalLocation const& loc) {
#ifdef USE_NBLOCKS
// Updated method, need to manually setup in configure.hpp, allow 6*n^2 blocks
*bound_lim = (int)(sqrt(NBLOCKS / 6) - 0.5);
// Find relative location within block

*lv2_lx2 = loc.lx2 / (*bound_lim + 1);
*lv2_lx3 = loc.lx3 / (*bound_lim + 1);
*local_lx2 = loc.lx2 - (*lv2_lx2 * (*bound_lim + 1));
*local_lx3 = loc.lx3 - (*lv2_lx3 * (*bound_lim + 1));
#else
// Old method, suitable for 6*4^n blocks
*bound_lim = (1 << (loc.level - 2)) - 1;
// Find relative location within block
*lv2_lx2 = loc.lx2 >> (loc.level - 2);
*lv2_lx3 = loc.lx3 >> (loc.level - 2);
*local_lx2 = loc.lx2 - (*lv2_lx2 << (loc.level - 2));
*local_lx3 = loc.lx3 - (*lv2_lx3 << (loc.level - 2));
#endif
}
18 changes: 10 additions & 8 deletions src/exo3/cs_transform_ox.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -12,16 +12,12 @@ void CubedSphere::TransformOX(int *ox2, int *ox3, int *tox2, int *tox3,
LogicalLocation const &loc) {
// Find the block ID
int block_id = FindBlockID(loc);
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);

// Find relative location within block
int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int local_lx2 = loc.lx2 - (lv2_lx2 << (loc.level - 2));
int local_lx3 = loc.lx3 - (lv2_lx3 << (loc.level - 2));
int bound_lim = (1 << (loc.level - 2)) - 1;

// Hard code the cases...
// No need to consider the corner cases, abandon in reading buffers.
// Hard code the cases...
// The corner cases are not used, abandoned after reading buffers.
int target_block = -1; // Block id of target
int target_loc_2, target_loc_3; // local x2 and x3 in target block

Expand Down Expand Up @@ -168,6 +164,7 @@ void CubedSphere::TransformOX(int *ox2, int *ox3, int *tox2, int *tox3,
std::stringstream msg;
msg << "Error: something wrong, check the geometry setup of the cubed "
"sphere. \n";
msg << "Block ID: " << block_id << std::endl;
msg << "----------------------------------" << std::endl;
ATHENA_ERROR(msg);
}
Expand Down Expand Up @@ -210,8 +207,13 @@ void CubedSphere::TransformOX(int *ox2, int *ox3, int *tox2, int *tox3,
msg << "----------------------------------" << std::endl;
ATHENA_ERROR(msg);
}
#ifdef USE_NBLOCKS
lx3_0 = (lx3_0 * (bound_lim + 1));
lx2_0 = (lx2_0 * (bound_lim + 1));
#else
lx3_0 = (lx3_0 << (loc.level - 2));
lx2_0 = (lx2_0 << (loc.level - 2));
#endif
// Add up first block and local positions
int lx3_t = lx3_0 + target_loc_3;
int lx2_t = lx2_0 + target_loc_2;
Expand Down
85 changes: 41 additions & 44 deletions src/exo3/cubed_sphere.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -39,9 +39,10 @@ CubedSphere::CubedSphere(MeshBlock *pmb) : pmy_block_(pmb) {

Real CubedSphere::GenerateMeshX2(Real x, LogicalLocation const &loc) {
Real x_l, x_u;
int lx2_lv2 = loc.lx2 >> (loc.level - 2);
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);

switch (lx2_lv2) {
switch (lv2_lx2) {
case 0:
x_l = -0.5;
x_u = 0.0;
Expand All @@ -56,9 +57,10 @@ Real CubedSphere::GenerateMeshX2(Real x, LogicalLocation const &loc) {

Real CubedSphere::GenerateMeshX3(Real x, LogicalLocation const &loc) {
Real x_l, x_u;
int lx3_lv2 = loc.lx3 >> (loc.level - 2);
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);

switch (lx3_lv2) {
switch (lv2_lx3) {
case 0:
x_l = -0.5;
x_u = -1.0 / 6.0;
Expand All @@ -81,10 +83,8 @@ Real CubedSphere::GenerateMeshX3(Real x, LogicalLocation const &loc) {
void CubedSphere::GetLatLon(Real *lat, Real *lon, int k, int j, int i) const {
auto pcoord = pmy_block_->pcoord;
auto &loc = pmy_block_->loc;
int blockID = FindBlockID(loc);

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
// Calculate the needed parameters
Real dX = tan(pcoord->x2v(j));
Real dY = tan(pcoord->x3v(k));
Expand All @@ -99,10 +99,7 @@ void CubedSphere::GetLatLonFace2(Real *lat, Real *lon, int k, int j,
int i) const {
auto pcoord = pmy_block_->pcoord;
auto &loc = pmy_block_->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
int blockID = FindBlockID(loc);

// Calculate the needed parameters
Real dX = tan(pcoord->x2f(j));
Expand All @@ -118,10 +115,7 @@ void CubedSphere::GetLatLonFace3(Real *lat, Real *lon, int k, int j,
int i) const {
auto pcoord = pmy_block_->pcoord;
auto &loc = pcoord->pmy_block->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
int blockID = FindBlockID(loc);

// Calculate the needed parameters
Real dX = tan(pcoord->x2v(j));
Expand All @@ -136,10 +130,7 @@ void CubedSphere::GetUV(Real *U, Real *V, Real V2, Real V3, int k, int j,
int i) const {
auto pcoord = pmy_block_->pcoord;
auto &loc = pmy_block_->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
int blockID = FindBlockID(loc);

Real X = tan(pcoord->x2v(j));
Real Y = tan(pcoord->x3v(k));
Expand All @@ -153,10 +144,7 @@ void CubedSphere::GetVyVz(Real *V2, Real *V3, Real U, Real V, int k, int j,
int i) const {
auto pcoord = pmy_block_->pcoord;
auto &loc = pmy_block_->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
int blockID = FindBlockID(loc);

// Calculate the needed parameters
Real X = tan(pcoord->x2v(j));
Expand All @@ -169,10 +157,7 @@ void CubedSphere::CalculateCoriolisForce2(int i2, int i3, Real v2, Real v3,
Real *cF3) const {
auto pcoord = pmy_block_->pcoord;
auto &loc = pmy_block_->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
int blockID = FindBlockID(loc);

Real x = tan(pcoord->x2v(i2));
Real y = tan(pcoord->x3v(i3));
Expand Down Expand Up @@ -203,10 +188,7 @@ void CubedSphere::CalculateCoriolisForce3(int i2, int i3, Real v1, Real v2,
Real *cF3) const {
auto pcoord = pmy_block_->pcoord;
auto &loc = pmy_block_->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
int blockID = FindBlockID(loc);

Real x = tan(pcoord->x2v(i2));
Real y = tan(pcoord->x3v(i3));
Expand Down Expand Up @@ -333,14 +315,10 @@ void CubedSphere::sendNeighborBlocks(int ox2, int ox3, int tg_rank,
int tg_gid) {
MeshBlock *pmb = pmy_block_;
auto &loc = pmb->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int blockID = lv2_lx2 + lv2_lx3 * 2 + 1;
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);
int blockID = FindBlockID(loc);
// Calculate local ID
int local_lx2 = loc.lx2 - (lv2_lx2 << (loc.level - 2));
int local_lx3 = loc.lx3 - (lv2_lx3 << (loc.level - 2));
int bound_lim = (1 << (loc.level - 2)) - 1;
int tox2, tox3;
int DirTag, DirNum, ownTag; // Tag for target, and the numbering of axis
int ox2_bkp = ox2;
Expand Down Expand Up @@ -560,6 +538,16 @@ void CubedSphere::sendNeighborBlocks(int ox2, int ox3, int tg_rank,
}

// Calculate the tag of destination
#ifdef USE_NBLOCKS
if (tox2 == -1)
DirTag = 0 + 4 * pmb->gid + 24 * (bound_lim + 1) * tg_gid;
if (tox2 == 1)
DirTag = 1 + 4 * pmb->gid + 24 * (bound_lim + 1) * tg_gid;
if (tox3 == -1)
DirTag = 2 + 4 * pmb->gid + 24 * (bound_lim + 1) * tg_gid;
if (tox3 == 1)
DirTag = 3 + 4 * pmb->gid + 24 * (bound_lim + 1) * tg_gid;
#else
if (tox2 == -1)
DirTag = 0 + 4 * pmb->gid + 24 * (1 << (loc.level - 2)) * tg_gid;
if (tox2 == 1)
Expand All @@ -568,6 +556,7 @@ void CubedSphere::sendNeighborBlocks(int ox2, int ox3, int tg_rank,
DirTag = 2 + 4 * pmb->gid + 24 * (1 << (loc.level - 2)) * tg_gid;
if (tox3 == 1)
DirTag = 3 + 4 * pmb->gid + 24 * (1 << (loc.level - 2)) * tg_gid;
#endif
// Send by MPI: we don't care whether it is in the same process for now
if (ox2 == -1) ownTag = 0;
if (ox2 == 1) ownTag = 1;
Expand All @@ -586,13 +575,10 @@ void CubedSphere::recvNeighborBlocks(int ox2, int ox3, int tg_rank,
MeshBlock *pmb = pmy_block_;
auto &loc = pmb->loc;

int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);
int blockID = FindBlockID(loc);
// Calculate local ID
int local_lx2 = loc.lx2 - (lv2_lx2 << (loc.level - 2));
int local_lx3 = loc.lx3 - (lv2_lx3 << (loc.level - 2));
int bound_lim = (1 << (loc.level - 2)) - 1;

int tox2, tox3;
int DirTag, DirNum, ownTag; // Tag for receiving, and the numbering of axis
// to place the values
Expand Down Expand Up @@ -658,6 +644,16 @@ void CubedSphere::recvNeighborBlocks(int ox2, int ox3, int tg_rank,
int dsize = ((kb2 - kb1 + 1) * (jb2 - jb1 + 1) * (ib2 - ib1 + 1) * NWAVE);
Real *data = new Real[dsize];
// Calculate the tag for receiving
#ifdef USE_NBLOCKS
if (ox2 == -1)
DirTag = 0 + 4 * tg_gid + 24 * (bound_lim + 1) * pmb->gid;
if (ox2 == 1)
DirTag = 1 + 4 * tg_gid + 24 * (bound_lim + 1) * pmb->gid;
if (ox3 == -1)
DirTag = 2 + 4 * tg_gid + 24 * (bound_lim + 1) * pmb->gid;
if (ox3 == 1)
DirTag = 3 + 4 * tg_gid + 24 * (bound_lim + 1) * pmb->gid;
#else
if (ox2 == -1)
DirTag = 0 + 4 * tg_gid + 24 * (1 << (loc.level - 2)) * pmb->gid;
if (ox2 == 1)
Expand All @@ -666,6 +662,7 @@ void CubedSphere::recvNeighborBlocks(int ox2, int ox3, int tg_rank,
DirTag = 2 + 4 * tg_gid + 24 * (1 << (loc.level - 2)) * pmb->gid;
if (ox3 == 1)
DirTag = 3 + 4 * tg_gid + 24 * (1 << (loc.level - 2)) * pmb->gid;
#endif

#ifdef MPI_PARALLEL
MPI_Recv(data, dsize, MPI_DOUBLE, tg_rank, DirTag, MPI_COMM_WORLD,
Expand Down
4 changes: 3 additions & 1 deletion src/exo3/cubed_sphere.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,9 @@ class CubedSphere {
~CubedSphere() {}

static int FindBlockID(LogicalLocation const &loc);
static void GetLocalIndex(int *lv2_lx2, int *lv2_lx3, int *local_lx2,
int *local_lx3, int *bound_lim,
LogicalLocation const &loc);
static void TransformOX(int *ox2, int *ox3, int *tox2, int *tox3,
LogicalLocation const &loc);

Expand Down Expand Up @@ -66,7 +69,6 @@ class CubedSphere {
#endif

std::vector<Real> LRDataBuffer[4];

MeshBlock *pmy_block_;
};

Expand Down
21 changes: 8 additions & 13 deletions src/exo3/cubed_sphere_utility.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -102,10 +102,9 @@ void InteprolateX2(const AthenaArray<Real> &src, AthenaArray<Real> &tgt,
int sk, int ek, int TypeFlag) {
// Interpolation along X2 (j) axis, used before sending data to X3 (k) axis
// Get the local indices
int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int local_lx2 = loc.lx2 - (lv2_lx2 << (loc.level - 2));
int bound_lim = (1 << (loc.level - 2)) - 1;
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
CubedSphere::GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);

int meshblock_size = ej - sj + 1;
int N_blk = meshblock_size *
(bound_lim + 1); // N in X2 direction for each panel. This value
Expand Down Expand Up @@ -233,10 +232,9 @@ void InteprolateX3(const AthenaArray<Real> &src, AthenaArray<Real> &tgt,
int sk, int ek, int TypeFlag) {
// Interpolation along X3 (k) axis, used before sending data to ghost zone in
// X2 (j) direction Get the local indices
int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int local_lx3 = loc.lx3 - (lv2_lx3 << (loc.level - 2));
int bound_lim = (1 << (loc.level - 2)) - 1;
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
CubedSphere::GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);

int meshblock_size = ek - sk + 1;
int N_blk = meshblock_size *
(bound_lim + 1); // N in X2 direction for each panel. This value
Expand Down Expand Up @@ -403,11 +401,8 @@ void PackData(const AthenaArray<Real> &src, Real *buf, int sn, int en, int si,
}

// Get the local indices
int lv2_lx2 = loc.lx2 >> (loc.level - 2);
int lv2_lx3 = loc.lx3 >> (loc.level - 2);
int local_lx2 = loc.lx2 - (lv2_lx2 << (loc.level - 2));
int local_lx3 = loc.lx3 - (lv2_lx3 << (loc.level - 2));
int bound_lim = (1 << (loc.level - 2)) - 1;
int lv2_lx2, lv2_lx3, local_lx2, local_lx3, bound_lim;
CubedSphere::GetLocalIndex(&lv2_lx2, &lv2_lx3, &local_lx2, &local_lx3, &bound_lim, loc);

// Work on interpolation
AthenaArray<Real> interpolatedSrc;
Expand Down

0 comments on commit 7d8e2b5

Please sign in to comment.