Merge branch '71-bugs-of-blacs-context-when-sub-communicator-is-used' into 'master'

Resolve "Bugs of BLACS context when sub-communicator is used"

Closes #71

See merge request chase/chase-library/ChASE!98

Author: brunowu

external/scalapackpp/scalapackpp.hpp

@@ -11,6 +11,8 @@
 #include "algorithm/types.hpp"
 
 #ifdef HAS_SCALAPACK
+#include <mpi.h>
+
 /**
  * @page scalapackpp_namespace chase::linalg::scalapackpp Namespace
  * @brief A templated C++ interface to the ScaLAPACK library.
@@ -46,6 +48,12 @@ extern "C" void blacs_gridexit_(int*);
 extern "C" void blacs_gridmap_(int*, int*, int*, int*, int*);
 extern "C" int numroc_(std::size_t*, std::size_t*, int*, int*, int*);
 
+// C BLACS interface (used to bind BLACS contexts to an MPI_Comm instead of MPI_COMM_WORLD)
+extern "C" int Csys2blacs_handle(MPI_Comm);
+extern "C" void Cblacs_gridinit(int*, char*, int, int);
+extern "C" void Cblacs_gridmap(int*, int*, int, int, int);
+extern "C" void Cblacs_gridexit(int);
+
 extern "C" void pdgeqrf_(int*, int*, double*, int*, int*, int*, double*,
                          double*, int*, int*);
 extern "C" void psgeqrf_(int*, int*, float*, int*, int*, int*, float*, float*,

grid/mpiGrid2D.hpp

@@ -414,10 +414,10 @@ class MpiGrid2D : public MpiGrid2DBase
 #endif
 
 #ifdef HAS_SCALAPACK
-        int zero = 0;
         int one = 1;
-        int ictxt;
-        chase::linalg::scalapackpp::blacs_get_(&zero, &zero, &ictxt);
+        // Create BLACS contexts bound to the user communicator `comm_`
+        // instead of relying on the global BLACS system context.
+        int ictxt = chase::linalg::scalapackpp::Csys2blacs_handle(comm_);
         colComm1D_ctxt_ = ictxt;
         int userMap[dims_[0]];
         if (MajorOrder == GridMajor::ColMajor)
@@ -435,23 +435,24 @@ class MpiGrid2D : public MpiGrid2DBase
             }
         }
 
-        chase::linalg::scalapackpp::blacs_gridmap_(&colComm1D_ctxt_, userMap,
-                                                   &dims_[0], &dims_[0], &one);
-
-        int ictxt_2;
-        chase::linalg::scalapackpp::blacs_get_(&zero, &zero, &ictxt_2);
+        // 1D column communicator context (dims_[0] x 1 grid) using C BLACS
+        chase::linalg::scalapackpp::Cblacs_gridmap(
+            &colComm1D_ctxt_, userMap, dims_[0], dims_[0], 1);
+        
+        // 2D BLACS grid context over the same communicator    
+        int ictxt_2 = chase::linalg::scalapackpp::Csys2blacs_handle(comm_);
         comm2D_ctxt_ = ictxt_2;
         if (MajorOrder == GridMajor::RowMajor)
         {
             char major = 'R';
-            chase::linalg::scalapackpp::blacs_gridinit_(&comm2D_ctxt_, &major,
-                                                        &dims_[0], &dims_[1]);
+            chase::linalg::scalapackpp::Cblacs_gridinit(&comm2D_ctxt_, &major,
+                                                        dims_[0], dims_[1]);
         }
         else
         {
             char major = 'C';
-            chase::linalg::scalapackpp::blacs_gridinit_(&comm2D_ctxt_, &major,
-                                                        &dims_[0], &dims_[1]);
+            chase::linalg::scalapackpp::Cblacs_gridinit(&comm2D_ctxt_, &major,
+                                                        dims_[0], dims_[1]);
         }
 #endif
     }

linalg/distMatrix/distMatrix.hpp

@@ -2112,8 +2112,6 @@ class BlockBlockMatrix
         int one = 1;
         int info;
         int comm2D_ctxt = mpi_grid_.get()->get_blacs_comm2D_ctxt();
-        int grank;
-        MPI_Comm_rank(MPI_COMM_WORLD, &grank);
         if constexpr (std::is_same<Platform, chase::platform::GPU>::value)
         {
             if (local_matrix_.cpu_data() == nullptr)