-
Notifications
You must be signed in to change notification settings - Fork 3
/
asagi.cpp
274 lines (226 loc) · 6.83 KB
/
asagi.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
/**
* @file
* This file is part of ASAGI.
*
* ASAGI is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as
* published by the Free Software Foundation, either version 3 of
* the License, or (at your option) any later version.
*
* ASAGI is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with ASAGI. If not, see <http://www.gnu.org/licenses/>.
*
* Diese Datei ist Teil von ASAGI.
*
* ASAGI ist Freie Software: Sie koennen es unter den Bedingungen
* der GNU Lesser General Public License, wie von der Free Software
* Foundation, Version 3 der Lizenz oder (nach Ihrer Option) jeder
* spaeteren veroeffentlichten Version, weiterverbreiten und/oder
* modifizieren.
*
* ASAGI wird in der Hoffnung, dass es nuetzlich sein wird, aber
* OHNE JEDE GEWAEHELEISTUNG, bereitgestellt; sogar ohne die implizite
* Gewaehrleistung der MARKTFAEHIGKEIT oder EIGNUNG FUER EINEN BESTIMMTEN
* ZWECK. Siehe die GNU Lesser General Public License fuer weitere Details.
*
* Sie sollten eine Kopie der GNU Lesser General Public License zusammen
* mit diesem Programm erhalten haben. Wenn nicht, siehe
* <http://www.gnu.org/licenses/>.
*
* @copyright 2012-2015 Sebastian Rettenberger <[email protected]>
*/
#include "asagi.h"
#include <cassert>
#include "utils/logger.h"
#include "grid/grid.h"
#ifndef ASAGI_NOMPI
#include "mpi/commthread.h"
#endif // ASAGI_NOMPI
// Static c++ functions
asagi::Grid* asagi::Grid::create(Type type)
{
return new grid::Grid(type);
}
asagi::Grid* asagi::Grid::createArray(Type type)
{
return new grid::Grid(type, true);
}
asagi::Grid* asagi::Grid::createStruct(unsigned int count,
unsigned int blockLength[], unsigned long displacements[], Type types[])
{
return new grid::Grid(count, blockLength, displacements, types);
}
#ifndef ASAGI_NOMPI
asagi::Grid::Error asagi::Grid::startCommThread(int schedCPU, MPI_Comm comm)
{
return mpi::CommThread::commThread.init(schedCPU, comm);
}
void asagi::Grid::stopCommThread()
{
mpi::CommThread::commThread.finialize();
}
int asagi::Grid::nodeLocalRank(MPI_Comm comm)
{
// The main idea for this function is taken from:
// https://blogs.fau.de/wittmann/2013/02/mpi-node-local-rank-determination/
// http://git.rrze.uni-erlangen.de/gitweb/?p=apsm.git;a=blob;f=MpiNodeRank.cpp;hb=HEAD
int mpiResult; NDBG_UNUSED(mpiResult);
typedef char procName_t[MPI_MAX_PROCESSOR_NAME+1];
// Get the processor name
procName_t procName;
int procNameLength;
MPI_Get_processor_name(procName, &procNameLength);
assert(procNameLength <= MPI_MAX_PROCESSOR_NAME);
procName[procNameLength] = '\0';
// Compute Adler32 hash
const uint8_t* buffer = reinterpret_cast<const uint8_t*>(procName);
uint32_t s1 = 1;
uint32_t s2 = 0;
for (int i = 0; i < procNameLength; i++) {
s1 = (s1 + buffer[i]) % 65521;
s2 = (s2 + s1) % 65521;
}
uint32_t hash = (s2 << 16) | s1;
int rank;
mpiResult = MPI_Comm_rank(comm, &rank);
assert(mpiResult == MPI_SUCCESS);
MPI_Comm nodeComm;
mpiResult = MPI_Comm_split(comm, hash, rank, &nodeComm);
assert(mpiResult == MPI_SUCCESS);
// Gather all proc names of this node to detect Adler32 collisions
int nodeSize;
mpiResult = MPI_Comm_size(nodeComm, &nodeSize);
assert(mpiResult == MPI_SUCCESS);
procName_t* procNames = new procName_t[nodeSize];
mpiResult = MPI_Allgather(procName, MPI_MAX_PROCESSOR_NAME+1, MPI_CHAR,
procNames, MPI_MAX_PROCESSOR_NAME+1, MPI_CHAR, nodeComm);
assert(mpiResult == MPI_SUCCESS);
// recv contains now an array of hostnames from all MPI ranks of
// this communicator. They are sorted ascending by the MPI rank.
int nodeRank, realNodeRank = 0;
mpiResult = MPI_Comm_rank(nodeComm, &nodeRank);
assert(mpiResult == MPI_SUCCESS);
for (int i = 0; i < nodeRank; i++) {
if (strcmp(procName, procNames[i]) == 0)
// Detect false hash collisions
realNodeRank++;
}
mpiResult = MPI_Comm_free(&nodeComm);
assert(mpiResult == MPI_SUCCESS);
delete [] procNames;
return realNodeRank;
}
#endif // ASAGI_NOMPI
// C interfae
// Init functions
asagi_grid* asagi_grid_create(asagi_type type)
{
return asagi::Grid::create(type);
}
asagi_grid* asagi_grid_create_array(asagi_type basic_type)
{
return asagi::Grid::createArray(basic_type);
}
asagi_grid* asagi_grid_create_struct(unsigned int count,
unsigned int blockLength[],
unsigned long displacements[],
asagi_type types[])
{
return asagi::Grid::createStruct(count, blockLength, displacements, types);
}
#ifndef ASAGI_NOMPI
void asagi_grid_set_comm(asagi_grid* handle, MPI_Comm comm)
{
handle->setComm(comm);
}
#endif // ASAIG_NOMPI
void asagi_grid_set_threads(asagi_grid* handle, unsigned int threads)
{
handle->setThreads(threads);
}
void asagi_grid_set_param(asagi_grid* handle, const char* name,
const char* value, unsigned int level)
{
handle->setParam(name, value, level);
}
asagi_error asagi_grid_open(asagi_grid* handle, const char* filename,
unsigned int level)
{
return handle->open(filename, level);
}
unsigned int asagi_grid_dimensions(asagi_grid* handle)
{
return handle->getDimensions();
}
double asagi_grid_min(asagi_grid* handle, unsigned int n)
{
return handle->getMin(n);
}
double asagi_grid_max(asagi_grid* handle, unsigned int n)
{
return handle->getMax(n);
}
double asagi_grid_delta(asagi_grid* handle, unsigned int n,
unsigned int level)
{
return handle->getDelta(n, level);
}
unsigned int asagi_grid_var_size(asagi_grid* handle)
{
return handle->getVarSize();
}
// Getters
unsigned char asagi_grid_get_byte(asagi_grid* handle, const double* pos,
unsigned int level)
{
return handle->getByte(pos, level);
}
int asagi_grid_get_int(asagi_grid* handle, const double* pos,
unsigned int level)
{
return handle->getInt(pos, level);
}
long asagi_grid_get_long(asagi_grid* handle, const double* pos,
unsigned int level)
{
return handle->getLong(pos, level);
}
float asagi_grid_get_float(asagi_grid* handle, const double* pos,
unsigned int level)
{
return handle->getFloat(pos, level);
}
double asagi_grid_get_double(asagi_grid* handle, const double* pos,
unsigned int level)
{
return handle->getDouble(pos, level);
}
void asagi_grid_get_buf(asagi_grid* handle, void* buf, const double* pos,
unsigned int level)
{
handle->getBuf(buf, pos, level);
}
// destructor
void asagi_grid_close(asagi_grid* handle)
{
delete handle;
}
#ifndef ASAGI_NOMPI
asagi_error asagi_start_comm_thread(int sched_cpu, MPI_Comm comm)
{
return asagi::Grid::startCommThread(sched_cpu, comm);
}
void asagi_stop_comm_thread()
{
asagi::Grid::stopCommThread();
}
int asagi_node_local_rank(MPI_Comm comm)
{
return asagi::Grid::nodeLocalRank(comm);
}
#endif // ASAGI_NOMPI