Line data Source code
1 : // This file is part of_outputCells_Localconditions of distribution
2 : // and use, please see the copyright notice in Mamico's main folder, or at
3 : // www5.in.tum.de/mamico
4 :
5 : // Implementation of SequentialFilter.h
6 :
7 : // CONSTRUCTOR
8 : template <unsigned int dim>
9 0 : coupling::filtering::SequentialFilter<dim>::SequentialFilter(coupling::filtering::FilterInterface<dim>* filter, const MPI_Comm comm)
10 0 : : coupling::filtering::FilterInterface<dim>("SEQUENTIALIZED_FILTER"), _filter(filter), _comm(comm), _processingRank(0), // TODO: Make this customizable.
11 0 : _cellbuf() {
12 0 : MPI_Comm_rank(_comm, &_myRank);
13 0 : MPI_Comm_size(_comm, &_commSize);
14 :
15 : #ifdef DEBUG_SEQ_FILTER
16 : std::cout << " SEQFILTER (" << _myRank << "): Now initializing sequentialized Filter of type: " << _filter->getType() << std::endl;
17 : #endif
18 :
19 0 : if (_commSize > 1) {
20 : // allocate global cell data structures for processing rank only
21 0 : if (_processingRank == _myRank) {
22 : // filter's input cell pointers will be overwritten, but we need this
23 : // information in contribute()
24 :
25 0 : _inputCells_Local = _filter->getInputCells();
26 0 : _outputCells_Local = _filter->getOutputCells();
27 :
28 : // allocate global cell datastructures
29 0 : unsigned int numcells = coupling::indexing::CellIndex<dim,
30 : /*global*/ coupling::indexing::IndexTrait::md2macro,
31 : coupling::indexing::IndexTrait::noGhost>::linearNumberCellsInDomain;
32 0 : for (unsigned int c = 0; c < numcells; c++) {
33 0 : _inputCells_Global.push_back(new coupling::datastructures::CouplingCell<dim>());
34 0 : _outputCells_Global.push_back(new coupling::datastructures::CouplingCell<dim>());
35 : }
36 0 : _cellRanks.resize(numcells);
37 : } // end: if rank processing rank
38 :
39 : // Resize if buffers are too small to contain cell data
40 0 : if (_cellbuf.size() < 4 + 3 * dim) {
41 0 : _cellbuf.resize(4 + 3 * dim);
42 : #ifdef DEBUG_SEQ_FILTER
43 : std::cout << " SEQFILTER (" << _myRank << "): Cell buffer resized to: " << _cellbuf.size() << std::endl;
44 : #endif
45 : }
46 : } // end: if _commSize > 1
47 :
48 : // case: _commSize = 1, i.e. single process
49 : else {
50 0 : _commSize = 1;
51 : #ifdef DEBUG_SEQ_FILTER
52 : std::cout << " SEQFILTER (" << _myRank << "): Warning: SequentialFilter operating on one rank only." << std::endl;
53 : #endif
54 : }
55 0 : }
56 :
57 : // MEMBER FUNCTIONS
58 0 : template <unsigned int dim> void coupling::filtering::SequentialFilter<dim>::operator()() {
59 0 : if (_commSize > 1) {
60 0 : if (_processingRank == _myRank) {
61 0 : contribute();
62 0 : process(FILTER_SEQUENTIAL);
63 : } else
64 0 : contribute();
65 : } else
66 0 : process(FILTER_PARALLEL);
67 0 : }
68 :
69 : // PRIVATE FUNCTIONS
70 0 : template <unsigned int dim> void coupling::filtering::SequentialFilter<dim>::contribute() {
71 : using coupling::indexing::CellIndex;
72 : using coupling::indexing::IndexTrait;
73 :
74 : // buffer used locally for raw indices
75 0 : unsigned int indexbuf{};
76 :
77 : #ifdef DEBUG_SEQ_FILTER
78 : std::cout << " SEQFILTER (" << _myRank << "): Sending cells and indices." << std::endl;
79 : #endif
80 :
81 0 : if (_processingRank == _myRank) {
82 : // count cells that we update. This should be equal to
83 : // _inputCells_Global.size().
84 : #ifdef DEBUG_SEQ_FILTER
85 : unsigned int cellsUpdated = 0;
86 : #endif
87 :
88 : // Iterate over all (contributing) ranks...
89 0 : for (int rank = 0; rank < _commSize; rank++) {
90 :
91 : // Case: The contributing rank is the processing rank. In this case, we
92 : // only need to copy cell data.
93 0 : if (rank == _myRank) {
94 : #ifdef DEBUG_SEQ_FILTER
95 : std::cout << " SEQFILTER (" << _myRank << "): Copying " << _inputCells_Local.size() << " cells from local on processing rank " << rank
96 : << std::endl;
97 : #endif
98 : // copy all cell data to global data structures. TODO: use std list and
99 : // concat
100 0 : for (const auto& index : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>()) {
101 : // determine corresponding global index
102 0 : indexbuf = CellIndex<dim, /*global*/ IndexTrait::md2macro,
103 0 : IndexTrait::noGhost>(index).get(); // then convert back to unsigned int
104 :
105 0 : *(_inputCells_Global[indexbuf]) = *(_inputCells_Local[index.get()]);
106 : #ifdef DEBUG_SEQ_FILTER
107 : cellsUpdated++;
108 : #endif
109 : }
110 :
111 : #ifdef DEBUG_SEQ_FILTER
112 : std::cout << " SEQFILTER (" << _myRank << "): Copied a total of " << cellsUpdated << " cells." << std::endl;
113 : #endif
114 : }
115 :
116 : // Case: the contributing rank is not the processing rank. We need to
117 : // communicate cell data using MPI.
118 : else {
119 : // receive number of cells+indices that will be sent from that rank
120 : unsigned int numberOfCellsFromRank;
121 0 : MPI_Recv(&numberOfCellsFromRank, 1, MPI_UNSIGNED, rank, 0, _comm, MPI_STATUS_IGNORE);
122 : #ifdef DEBUG_SEQ_FILTER
123 : std::cout << " SEQFILTER (" << _myRank << "): Expecting a total of " << numberOfCellsFromRank << " cells and indices from rank " << rank
124 : << std::endl;
125 : #endif
126 :
127 : // Receive that many cells
128 0 : for (unsigned int c = 0; c < numberOfCellsFromRank; c++) { // TODO: receive all cells at once
129 : // receive cell
130 0 : MPI_Recv(_cellbuf.data(), 4 + 3 * dim, MPI_DOUBLE, rank, 0, _comm, MPI_STATUS_IGNORE);
131 : // receive global index
132 0 : MPI_Recv(&indexbuf, 1, MPI_UNSIGNED, rank, 0, _comm, MPI_STATUS_IGNORE);
133 :
134 : // insert data received into global cell data structures TODO: use std
135 : // list and concat
136 0 : bufferToCouplingCell(_cellbuf, _inputCells_Global[indexbuf]);
137 :
138 : // save rank for this cell to be able to send it back later
139 0 : _cellRanks[indexbuf] = rank;
140 :
141 : #ifdef DEBUG_SEQ_FILTER
142 : cellsUpdated++;
143 : #endif
144 :
145 : #ifdef DEBUG_SEQ_FILTER_VERBOSE
146 : std::cout << " SEQFILTER (" << _myRank << "): Received cell from rank " << rank << std::endl;
147 : #endif
148 : }
149 : }
150 :
151 : #ifdef DEBUG_SEQ_FILTER
152 : std::cout << " SEQFILTER (" << _myRank << "): Received cells and indices from rank " << rank << std::endl;
153 : #endif
154 : }
155 :
156 : // This rank now processes...
157 : } else /* contributing rank */ {
158 :
159 : // communicate to processing rank how many cells+indices will be send
160 0 : unsigned int numberOfCellsFromThisRank = _filter->getInputCells().size();
161 0 : MPI_Send(&numberOfCellsFromThisRank, 1, MPI_UNSIGNED, _processingRank, 0, _comm);
162 :
163 : // send that many cells (and indices)
164 0 : for (auto& i : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro,
165 : IndexTrait::noGhost>()) { // TODO: send all cells at once
166 : // fill buffer
167 0 : couplingCellToBuffer(_cellbuf, _filter->getInputCells()[i.get()]);
168 :
169 : // send cell
170 0 : MPI_Send(_cellbuf.data(), 4 + 3 * dim, MPI_DOUBLE, _processingRank, 0, _comm);
171 : // send global index
172 0 : indexbuf = CellIndex<dim, /*global*/ IndexTrait::md2macro, IndexTrait::noGhost>(i).get();
173 0 : MPI_Send(&indexbuf, 1, MPI_UNSIGNED, _processingRank, 0, _comm);
174 : }
175 :
176 : // the processing rank now processes...
177 :
178 : // receive output data
179 0 : for (auto& i [[maybe_unused]] : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro,
180 : IndexTrait::noGhost>()) { // TODO: receive all cells at once
181 : // receive cell
182 0 : MPI_Recv(_cellbuf.data(), 4 + 3 * dim, MPI_DOUBLE, _processingRank, 0, _comm, MPI_STATUS_IGNORE);
183 : // receive global index
184 0 : MPI_Recv(&indexbuf, 1, MPI_UNSIGNED, _processingRank, 0, _comm, MPI_STATUS_IGNORE);
185 :
186 : // convert received index
187 0 : indexbuf = CellIndex<dim, IndexTrait::local, IndexTrait::md2macro,
188 : IndexTrait::noGhost>( // get local index
189 0 : CellIndex<dim, /*global*/ IndexTrait::md2macro,
190 : IndexTrait::noGhost>(indexbuf) // by constructing and converting global index
191 : )
192 0 : .get(); // then convert back to unsigned int
193 :
194 : // find cell in filter's domain and apply
195 0 : bufferToCouplingCell(_cellbuf, _filter->getOutputCells()[indexbuf]);
196 : }
197 : }
198 0 : }
199 :
200 0 : template <unsigned int dim> void coupling::filtering::SequentialFilter<dim>::process(bool sequential) {
201 0 : if (sequential) {
202 0 : _filter->updateCellData(_inputCells_Global, _outputCells_Global);
203 : }
204 :
205 : #ifdef DEBUG_SEQ_FILTER
206 : std::cout << " SEQFILTER (" << _myRank << "): Now applying filter..." << std::endl;
207 : #endif
208 :
209 : // Apply _filter
210 0 : (*_filter)();
211 :
212 : #ifdef DEBUG_SEQ_FILTER
213 : std::cout << " SEQFILTER (" << _myRank << "): ...done applying filter." << std::endl;
214 : #endif
215 :
216 0 : if (sequential) {
217 : // Now ready to distribute data back to the data's original ranks...
218 0 : for (unsigned int c = 0, localCellsUpdated = 0; c < _outputCells_Global.size(); c++) {
219 : // fill buffers with cell data
220 0 : couplingCellToBuffer(_cellbuf, _outputCells_Global[c]);
221 :
222 : // get origin rank of cell
223 0 : int target_rank = _cellRanks[c];
224 :
225 : #ifdef DEBUG_SEQ_FILTER_VERBOSE
226 : std::cout << " SEQFILTER (" << _myRank << "): Sending cell " << c << " to rank " << target_rank << std::endl;
227 : #endif
228 :
229 : // Case: cell belongs to this rank, no communication needed
230 0 : if (target_rank == _processingRank) {
231 0 : bufferToCouplingCell(_cellbuf, _outputCells_Local[localCellsUpdated]);
232 0 : localCellsUpdated++;
233 :
234 : #ifdef DEBUG_SEQ_FILTER_VERBOSE
235 : std::cout << " SEQFILTER (" << _myRank << "): Wrote cell " << c << " to local cell." << std::endl;
236 : #endif
237 : }
238 : // Case: other rank, we need to send cell and index
239 : else {
240 : // send cell
241 0 : MPI_Send(_cellbuf.data(), 4 + 3 * dim, MPI_DOUBLE, target_rank, 0, _comm);
242 : // send corresponding global index
243 0 : MPI_Send(&c, 1, MPI_UNSIGNED, target_rank, 0, _comm);
244 : }
245 : }
246 :
247 : // Re-update _filter's internal member variables in order to ensure proper
248 : // functionality of filters after this one, i.e. make getOutputVector() work
249 : // as expected.
250 0 : _filter->updateCellData(_inputCells_Local, _outputCells_Local);
251 : }
252 0 : }
253 :
254 : template <unsigned int dim>
255 0 : void coupling::filtering::SequentialFilter<dim>::couplingCellToBuffer(std::vector<double>& buf, const coupling::datastructures::CouplingCell<dim>* cell) {
256 : // Resize if buffer is too small to contain cell data
257 0 : if (buf.capacity() < 4 + 3 * dim) {
258 0 : buf.resize(4 + 3 * dim);
259 : #ifdef DEBUG_SEQ_FILTER
260 : std::cout << " SEQFILTER (" << _myRank << "): Cell buffer resized to: " << buf.size() << std::endl;
261 : #endif
262 : }
263 :
264 : // copy cell data to buffer
265 0 : unsigned int i = 0;
266 0 : /*micro mass*/ buf[i] = cell->getMicroscopicMass();
267 0 : i++; // std::cout << buf.back() << " ";
268 0 : /*micro momentum*/ for (unsigned d = 0; d < dim; d++) {
269 0 : buf[i] = cell->getMicroscopicMomentum()[d];
270 0 : i++; /*std::cout << buf.back() << " ";*/
271 : }
272 0 : /*macro mass*/ buf[i] = cell->getMacroscopicMass();
273 0 : i++; // std::cout << buf.back() << " ";
274 0 : /*macro momentum*/ for (unsigned d = 0; d < dim; d++) {
275 0 : buf[i] = cell->getMacroscopicMomentum()[d];
276 0 : i++; /*std::cout << buf.back() << " ";*/
277 : }
278 0 : /*pot energy*/ buf[i] = cell->getPotentialEnergy();
279 0 : i++; // std::cout << buf.back() << " ";
280 0 : /*velocity*/ for (unsigned d = 0; d < dim; d++) {
281 0 : buf[i] = cell->getCurrentVelocity()[d];
282 0 : i++; /*std::cout << buf.back() << " ";*/
283 : }
284 0 : /*temperature*/ buf[i] = cell->getTemperature();
285 0 : i++; // std::cout << buf.back() << " ";
286 : // std::cout << std::endl;
287 0 : }
288 :
289 : template <unsigned int dim>
290 0 : void coupling::filtering::SequentialFilter<dim>::bufferToCouplingCell(const std::vector<double>& buf, coupling::datastructures::CouplingCell<dim>* cell) {
291 0 : if (buf.capacity() < 4 + 3 * dim)
292 0 : throw std::runtime_error("Buffer too small for cell data!");
293 :
294 0 : tarch::la::Vector<dim, double> mvec_buf;
295 0 : unsigned int i = 0;
296 : // copy buffer data to cell
297 0 : /*micro mass*/ cell->setMicroscopicMass(buf[i]);
298 0 : i++;
299 0 : /*micro momentum*/ for (unsigned int d = 0; d < dim; d++) {
300 0 : mvec_buf[d] = buf[i];
301 0 : i++;
302 : }
303 0 : cell->setMicroscopicMomentum(mvec_buf);
304 0 : /*macro mass*/ cell->setMacroscopicMass(buf[i]);
305 0 : i++;
306 0 : /*macro momentum*/ for (unsigned int d = 0; d < dim; d++) {
307 0 : mvec_buf[d] = buf[i];
308 0 : i++;
309 : }
310 0 : cell->setMacroscopicMomentum(mvec_buf);
311 0 : /*pot energy*/ cell->setPotentialEnergy(buf[i]);
312 0 : i++;
313 0 : /*velocity*/ for (unsigned int d = 0; d < dim; d++) {
314 0 : mvec_buf[d] = buf[i];
315 0 : i++;
316 : }
317 0 : cell->setCurrentVelocity(mvec_buf);
318 0 : /*temperature*/ cell->setTemperature(buf[i]);
319 0 : i++;
320 : // std::cout << cell << std::endl;
321 0 : }
|