1 //////////////////////////////////////////////////////////////////////////////////////
2 // This file is distributed under the University of Illinois/NCSA Open Source License.
3 // See LICENSE file in top directory for details.
4 //
5 // Copyright (c) 2020 QMCPACK developers.
6 //
7 // File developed by: Ken Esler, kpesler@gmail.com, University of Illinois at Urbana-Champaign
8 // Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign
9 // Anouar Benali, benali@anl.gov, Argonne National Laboratory
10 // Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign
11 // Cynthia Gu, zg1@ornl.gov, Oak Ridge National Laboratory
12 // Mark Dewing, markdewing@gmail.com, University of Illinois at Urbana-Champaign
13 // Peter Doak, doakpw@ornl.gov, Oak Ridge National Laboratory
14 //
15 // File created by: Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign
16 //////////////////////////////////////////////////////////////////////////////////////
17
18
19 #include "Communicate.h"
20 #include <iostream>
21 #include <sstream>
22 #include <string>
23 #include <cstdio>
24 #include <fstream>
25 #include "config.h"
26 #include "Platforms/Host/sysutil.h"
27 #include "Utilities/FairDivide.h"
28
29 #ifdef HAVE_MPI
30 #include "mpi3/shared_communicator.hpp"
31 #endif
32
33
34 //Global Communicator is created without initialization
35 Communicate* OHMMS::Controller = new Communicate;
36
37 //default constructor: ready for a serial execution
Communicate()38 Communicate::Communicate()
39 : myMPI(MPI_COMM_NULL), d_mycontext(0), d_ncontexts(1), d_groupid(0), d_ngroups(1), GroupLeaderComm(nullptr)
40 {}
41
42 #ifdef HAVE_MPI
Communicate(const mpi3::environment & env)43 Communicate::Communicate(const mpi3::environment& env) : GroupLeaderComm(nullptr) { initialize(env); }
44 #endif
45
~Communicate()46 Communicate::~Communicate()
47 {
48 if (GroupLeaderComm != nullptr)
49 delete GroupLeaderComm;
50 }
51
52 //exclusive: MPI or Serial
53 #ifdef HAVE_MPI
54
Communicate(const mpi3::communicator & in_comm)55 Communicate::Communicate(const mpi3::communicator& in_comm) : d_groupid(0), d_ngroups(1), GroupLeaderComm(nullptr)
56 {
57 // const_cast is used to enable non-const call to duplicate()
58 comm = const_cast<mpi3::communicator&>(in_comm).duplicate();
59 myMPI = comm.get();
60 d_mycontext = comm.rank();
61 d_ncontexts = comm.size();
62 }
63
64
Communicate(const Communicate & in_comm,int nparts)65 Communicate::Communicate(const Communicate& in_comm, int nparts)
66 {
67 std::vector<int> nplist(nparts + 1);
68 int p = FairDivideLow(in_comm.rank(), in_comm.size(), nparts, nplist); //group
69 // const_cast is used to enable non-const call to split()
70 comm = const_cast<mpi3::communicator&>(in_comm.comm).split(p, in_comm.rank());
71 myMPI = comm.get();
72 // TODO: mpi3 needs to define comm
73 d_mycontext = comm.rank();
74 d_ncontexts = comm.size();
75 d_groupid = p;
76 d_ngroups = nparts;
77 // create a communicator among group leaders.
78
79 nplist.pop_back();
80 mpi3::communicator leader_comm = in_comm.comm.subcomm(nplist);
81 if (isGroupLeader())
82 GroupLeaderComm = new Communicate(leader_comm);
83 else
84 GroupLeaderComm = nullptr;
85 }
86
87
initialize(const mpi3::environment & env)88 void Communicate::initialize(const mpi3::environment& env)
89 {
90 comm = env.get_world_instance();
91 myMPI = comm.get();
92 d_mycontext = comm.rank();
93 d_ncontexts = comm.size();
94 d_groupid = 0;
95 d_ngroups = 1;
96 #ifdef __linux__
97 for (int proc = 0; proc < OHMMS::Controller->size(); proc++)
98 {
99 if (OHMMS::Controller->rank() == proc)
100 {
101 fprintf(stderr, "Rank = %4d Free Memory = %5zu MB\n", proc, (freemem() >> 20));
102 }
103 comm.barrier();
104 }
105 comm.barrier();
106 #endif
107 std::string when = "qmc." + getDateAndTime("%Y%m%d_%H%M");
108 }
109
110 // For unit tests until they can be changed and this will be removed.
initialize(int argc,char ** argv)111 void Communicate::initialize(int argc, char** argv) {}
112
initializeAsNodeComm(const Communicate & parent)113 void Communicate::initializeAsNodeComm(const Communicate& parent)
114 {
115 // const_cast is used to enable non-const call to split_shared()
116 comm = const_cast<mpi3::communicator&>(parent.comm).split_shared();
117 myMPI = comm.get();
118 d_mycontext = comm.rank();
119 d_ncontexts = comm.size();
120 }
121
finalize()122 void Communicate::finalize()
123 {
124 static bool has_finalized = false;
125
126 if (!has_finalized)
127 {
128 has_finalized = true;
129 }
130 }
131
cleanupMessage(void *)132 void Communicate::cleanupMessage(void*) {}
133
abort() const134 void Communicate::abort() const { comm.abort(1); }
135
barrier() const136 void Communicate::barrier() const { comm.barrier(); }
137 #else
138
initialize(int argc,char ** argv)139 void Communicate::initialize(int argc, char** argv) { std::string when = "qmc." + getDateAndTime("%Y%m%d_%H%M"); }
140
initializeAsNodeComm(const Communicate & parent)141 void Communicate::initializeAsNodeComm(const Communicate& parent) {}
142
finalize()143 void Communicate::finalize() {}
144
abort() const145 void Communicate::abort() const { std::abort(); }
146
barrier() const147 void Communicate::barrier() const {}
148
cleanupMessage(void *)149 void Communicate::cleanupMessage(void*) {}
150
Communicate(const Communicate & in_comm,int nparts)151 Communicate::Communicate(const Communicate& in_comm, int nparts)
152 : myMPI(MPI_COMM_NULL), d_mycontext(0), d_ncontexts(1), d_groupid(0)
153 {
154 GroupLeaderComm = new Communicate();
155 }
156 #endif // !HAVE_MPI
157
barrier_and_abort(const std::string & msg) const158 void Communicate::barrier_and_abort(const std::string& msg) const
159 {
160 if (!rank())
161 std::cerr << "Fatal Error. Aborting at " << msg << std::endl;
162 Communicate::barrier();
163 Communicate::abort();
164 }
165