MESA_PD_KernelBenchmark.cpp 21.2 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
//======================================================================================================================
//
//  This file is part of waLBerla. waLBerla is free software: you can
//  redistribute it and/or modify it under the terms of the GNU General Public
//  License as published by the Free Software Foundation, either version 3 of
//  the License, or (at your option) any later version.
//
//  waLBerla 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 General Public License
//  for more details.
//
//  You should have received a copy of the GNU General Public License along
//  with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
//
//! \file   MESA_PD_KernelBenchmark.cpp
//! \author Sebastian Eibl <sebastian.eibl@fau.de>
//
//======================================================================================================================

Sebastian Eibl's avatar
Sebastian Eibl committed
21
22
23
24
25
26
27
28
29
#include "check.h"
#include "Contact.h"
#include "CreateParticles.h"
#include "NodeTimings.h"
#include "Parameters.h"
#include "SelectProperty.h"
#include "sortParticleStorage.h"
#include "SQLProperties.h"

30
31
32
33
#include <mesa_pd/vtk/ParticleVtkOutput.h>

#include <mesa_pd/collision_detection/AnalyticContactDetection.h>
#include <mesa_pd/data/LinkedCells.h>
34
#include <mesa_pd/data/ParticleAccessorWithShape.h>
35
36
37
#include <mesa_pd/data/ParticleStorage.h>
#include <mesa_pd/data/ShapeStorage.h>
#include <mesa_pd/domain/BlockForestDomain.h>
38
#include <mesa_pd/kernel/AssocToBlock.h>
39
#include <mesa_pd/kernel/DoubleCast.h>
40
#include <mesa_pd/kernel/ExplicitEuler.h>
41
42
43
#include <mesa_pd/kernel/InsertParticleIntoLinkedCells.h>
#include <mesa_pd/kernel/ParticleSelector.h>
#include <mesa_pd/kernel/SpringDashpot.h>
44
#include <mesa_pd/kernel/SpringDashpotSpring.h>
45
#include <mesa_pd/mpi/ContactFilter.h>
46
#include <mesa_pd/mpi/ReduceContactHistory.h>
47
48
#include <mesa_pd/mpi/ReduceProperty.h>
#include <mesa_pd/mpi/SyncNextNeighbors.h>
49
#include <mesa_pd/mpi/SyncNextNeighborsBlockForest.h>
50
#include <mesa_pd/mpi/notifications/ForceTorqueNotification.h>
Sebastian Eibl's avatar
Sebastian Eibl committed
51
52
#include <mesa_pd/sorting/HilbertCompareFunctor.h>
#include <mesa_pd/sorting/LinearizedCompareFunctor.h>
53
54
55
56
57

#include <blockforest/BlockForest.h>
#include <blockforest/Initialization.h>
#include <core/Abort.h>
#include <core/Environment.h>
Sebastian Eibl's avatar
Sebastian Eibl committed
58
#include <core/Hostname.h>
59
#include <core/math/Random.h>
60
#include <core/math/Sample.h>
Sebastian Eibl's avatar
Sebastian Eibl committed
61
62
#include <core/mpi/Gatherv.h>
#include <core/mpi/RecvBuffer.h>
63
#include <core/mpi/Reduce.h>
Sebastian Eibl's avatar
Sebastian Eibl committed
64
#include <core/mpi/SendBuffer.h>
65
66
67
#include <core/grid_generator/SCIterator.h>
#include <core/logging/Logging.h>
#include <core/OpenMP.h>
68
#include <core/perf_analysis/extern/likwid.h>
69
70
71
#include <core/timing/Timer.h>
#include <core/timing/TimingPool.h>
#include <core/waLBerlaBuildInfo.h>
72
#include <sqlite/SQLite.h>
73
74
75
76
77
78
79
80
81
82
#include <vtk/VTKOutput.h>

#include <functional>
#include <memory>
#include <string>
#include <type_traits>

namespace walberla {
namespace mesa_pd {

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
class ContactDetection
{
public:
   ContactDetection(const std::shared_ptr<domain::BlockForestDomain>& domain)
   : domain_(domain)
   {
      contacts_.reserve(4000000);
   }

   template <typename T>
   inline
   void operator()(const size_t idx1, const size_t idx2, T& ac)
   {
      ++contactsChecked_;
      if (double_cast_(idx1, idx2, ac, acd_, ac))
      {
         ++contactsDetected_;
         if (contact_filter_(acd_.getIdx1(), acd_.getIdx2(), ac, acd_.getContactPoint(),
                             *domain_))
         {
            ++contactsTreated_;
            contacts_.emplace_back(acd_.getIdx1(), acd_.getIdx2(), acd_.getContactPoint(),
                 acd_.getContactNormal(), acd_.getPenetrationDepth());
         }
      }
   }

   inline const auto& getContacts() const {return contacts_;}

   inline
   void resetCounters()
   {
      contactsChecked_ = 0;
      contactsDetected_ = 0;
      contactsTreated_ = 0;
      contacts_.clear();
   }

   inline
   int64_t getContactsChecked() const
   {
      return contactsChecked_;
   }

   inline
   int64_t getContactsDetected() const
   {
      return contactsDetected_;
   }

   inline
   int64_t getContactsTreated() const
   {
      return contactsTreated_;
   }

private:
   kernel::DoubleCast double_cast_;
   mpi::ContactFilter contact_filter_;
   std::shared_ptr<domain::BlockForestDomain> domain_;
   collision_detection::AnalyticContactDetection acd_;
   std::vector<Contact> contacts_;
   int64_t contactsChecked_ = 0;
   int64_t contactsDetected_ = 0;
   int64_t contactsTreated_ = 0;
};

template <typename T>
void reportOverRanks(const std::string& info, const T& value)
{
   math::Sample sample;
   sample.insert(value);
   sample.mpiGatherRoot();
   WALBERLA_LOG_INFO_ON_ROOT(info);
   WALBERLA_LOG_INFO_ON_ROOT(sample.format());
}

160
161
int main( int argc, char ** argv )
{
162
163
   LIKWID_MARKER_INIT;

164
165
166
167
168
169
170
171
   using namespace walberla::timing;

   Environment env(argc, argv);
   auto mpiManager = walberla::mpi::MPIManager::instance();
   mpiManager->useWorldComm();

   logging::Logging::instance()->setStreamLogLevel(logging::Logging::INFO);
   logging::Logging::instance()->setFileLogLevel(logging::Logging::INFO);
172
//   logging::Logging::instance()->includeLoggingToFile("KernelBenchmark");
173
174
175
176
177
178
179
180
181
182

   WALBERLA_LOG_INFO_ON_ROOT( "config file: " << argv[1] );
   WALBERLA_LOG_INFO_ON_ROOT( "waLBerla Revision: " << WALBERLA_GIT_SHA1 );

   math::seedRandomGenerator( static_cast<unsigned int>(1337 * mpiManager->worldRank()) );

   WALBERLA_LOG_INFO_ON_ROOT("*** READING CONFIG FILE ***");
   auto cfg = env.config();
   if (cfg == nullptr) WALBERLA_ABORT("No config specified!");
   const Config::BlockHandle mainConf  = cfg->getBlock( "GranularGas" );
Sebastian Eibl's avatar
Sebastian Eibl committed
183
184
   Parameters params;
   loadFromConfig(params, mainConf);
185
186
187
188
189
190
191
192
193

   WALBERLA_LOG_INFO_ON_ROOT("*** BLOCKFOREST ***");
   // create forest
   auto forest = blockforest::createBlockForestFromConfig( mainConf );
   if (!forest)
   {
      WALBERLA_LOG_INFO_ON_ROOT( "No BlockForest created ... exiting!");
      return EXIT_SUCCESS;
   }
194
195
196
   auto domain = std::make_shared<domain::BlockForestDomain>(forest);

   LIKWID_MARKER_THREADINIT;
197
198
199
200
201
202
203

   auto simulationDomain = forest->getDomain();
   auto localDomain = forest->begin()->getAABB();
   for (auto& blk : *forest)
   {
      localDomain.merge(blk.getAABB());
   }
204
//   WALBERLA_LOG_DEVEL_VAR(localDomain);
205
206
207
208
209

   WALBERLA_LOG_INFO_ON_ROOT("*** SETUP - START ***");

   //init data structures
   auto ss = std::make_shared<data::ShapeStorage>();
210
   auto ps = std::make_shared<data::ParticleStorage>(100);
211
   data::ParticleAccessorWithShape accessor(ps, ss);
Sebastian Eibl's avatar
Sebastian Eibl committed
212
   data::LinkedCells     lc(localDomain.getExtended(params.spacing), params.spacing );
213

Sebastian Eibl's avatar
Sebastian Eibl committed
214
   auto  smallSphere = ss->create<data::Sphere>( params.radius );
215
216
217
   ss->shapes[smallSphere]->updateMassAndInertia(real_t(2707));
   for (auto& iBlk : *forest)
   {
Sebastian Eibl's avatar
Sebastian Eibl committed
218
      for (auto pt : grid_generator::SCGrid(iBlk.getAABB(),
219
                                            Vector3<real_t>(params.spacing) * real_c(0.5) + params.shift,
Sebastian Eibl's avatar
Sebastian Eibl committed
220
                                            params.spacing))
221
222
      {
         WALBERLA_CHECK(iBlk.getAABB().contains(pt));
Sebastian Eibl's avatar
Sebastian Eibl committed
223
         createSphere(*ps, pt, params.radius, smallSphere);
224
//         WALBERLA_LOG_DEVEL_VAR_ON_ROOT(pt);
225
226
227
228
229
230
      }
   }
   int64_t numParticles = int64_c(ps->size());
   walberla::mpi::reduceInplace(numParticles, walberla::mpi::SUM);
   WALBERLA_LOG_INFO_ON_ROOT("#particles created: " << numParticles);

Sebastian Eibl's avatar
Sebastian Eibl committed
231
   auto shift = (params.spacing - params.radius - params.radius) * real_t(0.5);
232
233
234
235
   auto confiningDomain = simulationDomain.getExtended(shift);

   if (!forest->isPeriodic(0))
   {
236
237
      createPlane(*ps, *ss, confiningDomain.minCorner() + params.shift, Vec3(+1,0,0));
      createPlane(*ps, *ss, confiningDomain.maxCorner() + params.shift, Vec3(-1,0,0));
238
239
240
241
   }

   if (!forest->isPeriodic(1))
   {
242
243
      createPlane(*ps, *ss, confiningDomain.minCorner() + params.shift, Vec3(0,+1,0));
      createPlane(*ps, *ss, confiningDomain.maxCorner() + params.shift, Vec3(0,-1,0));
244
245
246
247
   }

   if (!forest->isPeriodic(2))
   {
248
249
      createPlane(*ps, *ss, confiningDomain.minCorner() + params.shift, Vec3(0,0,+1));
      createPlane(*ps, *ss, confiningDomain.maxCorner() + params.shift, Vec3(0,0,-1));
250
251
252
253
254
   }

   WALBERLA_LOG_INFO_ON_ROOT("*** SETUP - END ***");

   WALBERLA_LOG_INFO_ON_ROOT("*** VTK ***");
Sebastian Eibl's avatar
Sebastian Eibl committed
255
   auto vtkDomainOutput = walberla::vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, params.vtk_out, "simulation_step" );
256
   auto vtkOutput       = make_shared<mesa_pd::vtk::ParticleVtkOutput>(ps) ;
Sebastian Eibl's avatar
Sebastian Eibl committed
257
   auto vtkWriter       = walberla::vtk::createVTKOutput_PointData(vtkOutput, "Bodies", 1, params.vtk_out, "simulation_step", false, false);
258
259
   vtkOutput->addOutput<SelectRank>("rank");
   vtkOutput->addOutput<data::SelectParticleOwner>("owner");
Sebastian Eibl's avatar
Sebastian Eibl committed
260
   vtkOutput->addOutput<SelectIdx>("idx");
261
//   vtkDomainOutput->write();
262
263
264

   WALBERLA_LOG_INFO_ON_ROOT("*** SIMULATION - START ***");
   // Init kernels
265
   kernel::ExplicitEuler                 explicitEuler( params.dt );
266
   kernel::AssocToBlock                  assoc(forest);
267
   kernel::InsertParticleIntoLinkedCells ipilc;
268
269
270
271
272
273
274
275
276
   kernel::SpringDashpot                 sd(1);
   sd.setDampingT (0, 0, real_t(0));
   sd.setFriction (0, 0, real_t(0));
   sd.setParametersFromCOR(0, 0, real_t(0.9), params.dt*real_t(20), ss->shapes[smallSphere]->getMass() * real_t(0.5));
   kernel::SpringDashpotSpring           sds(1);
   sds.setParametersFromCOR(0, 0, real_t(0.9), params.dt*real_t(20), ss->shapes[smallSphere]->getMass() * real_t(0.5));
   sds.setCoefficientOfFriction(0,0,real_t(0.4));
   sds.setStiffnessT(0,0,real_t(0.9) * sds.getStiffnessN(0,0));

277
278
   mpi::ReduceProperty                   RP;
   mpi::SyncNextNeighbors                SNN;
279
   mpi::ReduceContactHistory             RCH;
280
   ContactDetection                      CD(domain);
281
282

   // initial sync
283
284
   ps->forEachParticle(false, kernel::SelectLocal(), accessor, assoc, accessor);
   SNN(*ps, *domain);
Sebastian Eibl's avatar
Sebastian Eibl committed
285
286
   sortParticleStorage(*ps, params.sorting, lc.domain_, uint_c(lc.numCellsPerDim_[0]));
//   vtkWriter->write();
287

Sebastian Eibl's avatar
Sebastian Eibl committed
288
   for (int64_t outerIteration = 0; outerIteration < params.numOuterIterations; ++outerIteration)
289
290
291
292
293
   {
      WALBERLA_LOG_INFO_ON_ROOT("*** RUNNING OUTER ITERATION " << outerIteration << " ***");

      WcTimingPool tp;

294
295
296
297
298
299
300
301
302
303
304
      LIKWID_MARKER_REGISTER("SNN");
      WALBERLA_MPI_BARRIER();
      LIKWID_MARKER_START("SNN");
      tp["SNN"].start();
      for (int64_t i=0; i < params.simulationSteps; ++i)
      {
         SNN(*ps, *domain);
      }
      tp["SNN"].end();
      LIKWID_MARKER_STOP("SNN");

305
      LIKWID_MARKER_REGISTER("AssocToBlock");
306
      WALBERLA_MPI_BARRIER();
307
308
309
310
311
312
313
314
315
316
317
318
      LIKWID_MARKER_START("AssocToBlock");
      tp["AssocToBlock"].start();
      for (int64_t i=0; i < params.simulationSteps; ++i)
      {
         ps->forEachParticle(false, kernel::SelectLocal(), accessor, assoc, accessor);
      }
      tp["AssocToBlock"].end();
      LIKWID_MARKER_STOP("AssocToBlock");

      LIKWID_MARKER_REGISTER("GenerateLinkedCells");
      WALBERLA_MPI_BARRIER();
      LIKWID_MARKER_START("GenerateLinkedCells");
319
      tp["GenerateLinkedCells"].start();
Sebastian Eibl's avatar
Sebastian Eibl committed
320
      for (int64_t i=0; i < params.simulationSteps; ++i)
321
322
323
324
325
      {
         lc.clear();
         ps->forEachParticle(true, kernel::SelectAll(), accessor, ipilc, accessor, lc);
      }
      tp["GenerateLinkedCells"].end();
326
      LIKWID_MARKER_STOP("GenerateLinkedCells");
327

328
      LIKWID_MARKER_REGISTER("ContactDetection");
329
      WALBERLA_MPI_BARRIER();
330
      LIKWID_MARKER_START("ContactDetection");
331
      tp["ContactDetection"].start();
Sebastian Eibl's avatar
Sebastian Eibl committed
332
      for (int64_t i=0; i < params.simulationSteps; ++i)
333
      {
334
         CD.resetCounters();
335
336
337
         lc.forEachParticlePairHalf(true,
                                    kernel::SelectAll(),
                                    accessor,
338
339
                                    CD,
                                    accessor);
340
341
      }
      tp["ContactDetection"].end();
342
      LIKWID_MARKER_STOP("ContactDetection");
343

344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
      LIKWID_MARKER_REGISTER("SpringDashpot");
      WALBERLA_MPI_BARRIER();
      LIKWID_MARKER_START("SpringDashpot");
      tp["SpringDashpot"].start();
      for (int64_t i=0; i < params.simulationSteps; ++i)
      {
         for (auto& c : CD.getContacts())
         {
            sd(c.idx1_, c.idx2_, accessor, c.contactPoint_, c.contactNormal_, c.penetrationDepth_);
         }
      }
      tp["SpringDashpot"].end();
      LIKWID_MARKER_STOP("SpringDashpot");

      LIKWID_MARKER_REGISTER("SpringDashpotSpring");
359
      WALBERLA_MPI_BARRIER();
360
361
      LIKWID_MARKER_START("SpringDashpotSpring");
      tp["SpringDashpotSpring"].start();
Sebastian Eibl's avatar
Sebastian Eibl committed
362
      for (int64_t i=0; i < params.simulationSteps; ++i)
363
      {
364
         for (auto& c : CD.getContacts())
365
         {
366
            sds(c.idx1_, c.idx2_, accessor, c.contactPoint_, c.contactNormal_, c.penetrationDepth_, params.dt);
367
368
         }
      }
369
370
      tp["SpringDashpotSpring"].end();
      LIKWID_MARKER_STOP("SpringDashpotSpring");
371

372
      LIKWID_MARKER_REGISTER("ReduceForce");
373
      WALBERLA_MPI_BARRIER();
374
      LIKWID_MARKER_START("ReduceForce");
375
      tp["ReduceForce"].start();
Sebastian Eibl's avatar
Sebastian Eibl committed
376
      for (int64_t i=0; i < params.simulationSteps; ++i)
377
378
379
380
      {
         RP.operator()<ForceTorqueNotification>(*ps);
      }
      tp["ReduceForce"].end();
381
      LIKWID_MARKER_STOP("ReduceForce");
382

383
384
385
386
387
388
389
390
391
392
393
      LIKWID_MARKER_REGISTER("ReduceContactHistory");
      WALBERLA_MPI_BARRIER();
      LIKWID_MARKER_START("ReduceContactHistory");
      tp["ReduceContactHistory"].start();
      for (int64_t i=0; i < params.simulationSteps; ++i)
      {
         RCH(*ps);
      }
      tp["ReduceContactHistory"].end();
      LIKWID_MARKER_STOP("ReduceContactHistory");

394
      LIKWID_MARKER_REGISTER("Euler");
395
      WALBERLA_MPI_BARRIER();
396
      LIKWID_MARKER_START("Euler");
397
      tp["Euler"].start();
Sebastian Eibl's avatar
Sebastian Eibl committed
398
      for (int64_t i=0; i < params.simulationSteps; ++i)
399
      {
400
         ps->forEachParticle(true, kernel::SelectLocal(), accessor, explicitEuler, accessor);
401
402
      }
      tp["Euler"].end();
403
      LIKWID_MARKER_STOP("Euler");
404
405
406

      WALBERLA_LOG_INFO_ON_ROOT("*** SIMULATION - END ***");

Sebastian Eibl's avatar
Sebastian Eibl committed
407
      if (params.checkSimulation)
408
      {
409
410
         //if you want to activate checking you have to deactivate sorting
         check(*ps, *forest, params.spacing, params.shift);
411
412
413
      }

      WALBERLA_LOG_INFO_ON_ROOT("*** SQL OUTPUT - START ***");
414
415
416
417
418
419
420
421
      int64_t contactsChecked = CD.getContactsChecked();
      int64_t contactsDetected = CD.getContactsDetected();
      int64_t contactsTreated = CD.getContactsTreated();

      reportOverRanks("Contacts Checked:", real_c(contactsChecked));
      reportOverRanks("Contacts Detected:", real_c(contactsDetected));
      reportOverRanks("Contacts Treated:", real_c(contactsTreated));

422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
      auto SNNBytesSent     = SNN.getBytesSent();
      auto SNNBytesReceived = SNN.getBytesReceived();
      auto SNNSends         = SNN.getNumberOfSends();
      auto SNNReceives      = SNN.getNumberOfReceives();
      auto RPBytesSent      = RP.getBytesSent();
      auto RPBytesReceived  = RP.getBytesReceived();
      auto RPSends          = RP.getNumberOfSends();
      auto RPReceives       = RP.getNumberOfReceives();
      walberla::mpi::reduceInplace(SNNBytesSent, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(SNNBytesReceived, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(SNNSends, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(SNNReceives, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(RPBytesSent, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(RPBytesReceived, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(RPSends, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(RPReceives, walberla::mpi::SUM);
      auto cC = walberla::mpi::reduce(contactsChecked, walberla::mpi::SUM);
      auto cD = walberla::mpi::reduce(contactsDetected, walberla::mpi::SUM);
      auto cT = walberla::mpi::reduce(contactsTreated, walberla::mpi::SUM);
      WALBERLA_LOG_DEVEL_ON_ROOT( "SNN bytes communicated:   " << SNNBytesSent << " / " << SNNBytesReceived );
      WALBERLA_LOG_DEVEL_ON_ROOT( "SNN communication partners: " << SNNSends << " / " << SNNReceives );
      WALBERLA_LOG_DEVEL_ON_ROOT( "RP bytes communicated:  " << RPBytesSent << " / " << RPBytesReceived );
      WALBERLA_LOG_DEVEL_ON_ROOT( "RP communication partners: " << RPSends << " / " << RPReceives );
      WALBERLA_LOG_DEVEL_ON_ROOT( "contacts checked/detected/treated: " << cC << " / " << cD << " / " << cT );

      auto tp_reduced = tp.getReduced();
      WALBERLA_LOG_INFO_ON_ROOT(*tp_reduced);

      numParticles = 0;
      int64_t numGhostParticles = 0;
      ps->forEachParticle(false,
                          kernel::SelectAll(),
                          accessor,
                          [&numParticles, &numGhostParticles](const size_t idx, auto& ac)
      {
         if (data::particle_flags::isSet( ac.getFlagsRef(idx), data::particle_flags::GHOST))
         {
            ++numGhostParticles;
         } else
         {
            ++numParticles;
         }
      },
      accessor);
466
467
468
469
470
471
472
473
474
475

      reportOverRanks("Total Particles:", real_c(ps->size()));
      reportOverRanks("Number of Particles:", real_c(numParticles));
      reportOverRanks("Number of Ghost Particles:", real_c(numGhostParticles));
      auto estGhostParticles = (localDomain.xSize()+2) * (localDomain.ySize()+2) * (localDomain.zSize()+2);
      estGhostParticles -= localDomain.xSize() * localDomain.ySize() * localDomain.zSize();
      estGhostParticles -= 4*localDomain.xSize() + 4*localDomain.ySize() + 4*localDomain.zSize();
      estGhostParticles -= 8;
      WALBERLA_LOG_DEVEL_ON_ROOT("Estimation: " << estGhostParticles);

476
477
478
479
480
481
482
483
484
      walberla::mpi::reduceInplace(numParticles, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(numGhostParticles, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(contactsChecked, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(contactsDetected, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(contactsTreated, walberla::mpi::SUM);
      double linkedCellsVolume = lc.domain_.volume();
      walberla::mpi::reduceInplace(linkedCellsVolume, walberla::mpi::SUM);
      size_t numLinkedCells = lc.cells_.size();
      walberla::mpi::reduceInplace(numLinkedCells, walberla::mpi::SUM);
485
486
487
      size_t local_aabbs         = domain->getNumLocalAABBs();
      size_t neighbor_subdomains = domain->getNumNeighborSubdomains();
      size_t neighbor_processes  = domain->getNumNeighborProcesses();
488
489
490
      walberla::mpi::reduceInplace(local_aabbs, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(neighbor_subdomains, walberla::mpi::SUM);
      walberla::mpi::reduceInplace(neighbor_processes, walberla::mpi::SUM);
Sebastian Eibl's avatar
Sebastian Eibl committed
491
492

      uint_t runId = uint_c(-1);
493
494
495
496
497
498
499
500
501
502
      WALBERLA_ROOT_SECTION()
      {
         std::map< std::string, walberla::int64_t > integerProperties;
         std::map< std::string, double >            realProperties;
         std::map< std::string, std::string >       stringProperties;

         stringProperties["walberla_git"]         = WALBERLA_GIT_SHA1;
         stringProperties["tag"]                  = "mesa_pd";
         integerProperties["mpi_num_processes"]   = mpiManager->numProcesses();
         integerProperties["omp_max_threads"]     = omp_get_max_threads();
Sebastian Eibl's avatar
Sebastian Eibl committed
503
         integerProperties["outerIteration"]      = int64_c(outerIteration);
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
         integerProperties["num_particles"]       = numParticles;
         integerProperties["num_ghost_particles"] = numGhostParticles;
         integerProperties["contacts_checked"]    = contactsChecked;
         integerProperties["contacts_detected"]   = contactsDetected;
         integerProperties["contacts_treated"]    = contactsTreated;
         integerProperties["local_aabbs"]         = int64_c(local_aabbs);
         integerProperties["neighbor_subdomains"] = int64_c(neighbor_subdomains);
         integerProperties["neighbor_processes"]  = int64_c(neighbor_processes);
         integerProperties["SNNBytesSent"]        = SNNBytesSent;
         integerProperties["SNNBytesReceived"]    = SNNBytesReceived;
         integerProperties["SNNSends"]            = SNNSends;
         integerProperties["SNNReceives"]         = SNNReceives;
         integerProperties["RPBytesSent"]         = RPBytesSent;
         integerProperties["RPBytesReceived"]     = RPBytesReceived;
         integerProperties["RPSends"]             = RPSends;
         integerProperties["RPReceives"]          = RPReceives;
         realProperties["linkedCellsVolume"]      = linkedCellsVolume;
         integerProperties["numLinkedCells"]      = int64_c(numLinkedCells);
Sebastian Eibl's avatar
Sebastian Eibl committed
522

523
         addBuildInfoToSQL( integerProperties, realProperties, stringProperties );
Sebastian Eibl's avatar
Sebastian Eibl committed
524
525
526
527
         saveToSQL(params, integerProperties, realProperties, stringProperties );
         addDomainPropertiesToSQL(*forest, integerProperties, realProperties, stringProperties);
         addSlurmPropertiesToSQL(integerProperties, realProperties, stringProperties);

528
529
         runId = sqlite::storeRunInSqliteDB( params.sqlFile, integerProperties, stringProperties, realProperties );
         sqlite::storeTimingPoolInSqliteDB( params.sqlFile, runId, *tp_reduced, "Timeloop" );
Sebastian Eibl's avatar
Sebastian Eibl committed
530
531
532
533
      }
      if (params.storeNodeTimings)
      {
         storeNodeTimings(runId, params.sqlFile, "NodeTiming", tp);
534
535
536
537
      }
      WALBERLA_LOG_INFO_ON_ROOT("*** SQL OUTPUT - END ***");
   }

538
539
   LIKWID_MARKER_CLOSE;

540
541
542
543
544
545
546
547
548
549
   return EXIT_SUCCESS;
}

} // namespace mesa_pd
} // namespace walberla

int main( int argc, char* argv[] )
{
   return walberla::mesa_pd::main( argc, argv );
}