apps/benchmarks/CNT/01_cnt_film.cpp

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#include "Config.h"
+#include "InitializeCNTs.h"
+#include "SphericalSegmentAccessor.h"
+#include "SQLProperties.h"
+#include "Statistics.h"
+#include "TerminalColors.h"
+
+#include "mesa_pd/domain/BlockForestDomain.h"
+#include "mesa_pd/data/Flags.h"
+#include "mesa_pd/data/LinkedCells.h"
+#include "mesa_pd/data/ParticleStorage.h"
+#include "mesa_pd/data/SparseLinkedCells.h"
+#include <mesa_pd/kernel/AssocToBlock.h>
+#include "mesa_pd/kernel/InsertParticleIntoLinkedCells.h"
+#include "mesa_pd/kernel/InsertParticleIntoSparseLinkedCells.h"
+#include "mesa_pd/kernel/ParticleSelector.h"
+#include "mesa_pd/kernel/cnt/AnisotropicVDWContact.h"
+#include "mesa_pd/kernel/cnt/VBondContact.h"
+#include "mesa_pd/kernel/cnt/ViscousDamping.h"
+#include "mesa_pd/kernel/cnt/Parameters.h"
+#include "mesa_pd/kernel/PFCDamping.h"
+#include "mesa_pd/kernel/VelocityVerlet.h"
+#include "mesa_pd/mpi/ContactFilter.h"
+#include "mesa_pd/mpi/notifications/ForceTorqueNotification.h"
+#include "mesa_pd/mpi/ReduceProperty.h"
+#include "mesa_pd/mpi/SyncNextNeighbors.h"
+#include "mesa_pd/sorting/LinearizedCompareFunctor.h"
+#include "mesa_pd/vtk/ParticleVtkOutput.h"
+
+#include "blockforest/Initialization.h"
+#include "core/Environment.h"
+#include "core/math/Constants.h"
+#include "core/timing/TimingPool.h"
+#include "core/timing/TimingTree.h"
+#include "sqlite/SQLite.h"
+#include "vtk/VTKOutput.h"
+
+namespace walberla {
+namespace mesa_pd {
+
+int main(int argc, char **argv)
+{
+   Environment env(argc, argv);
+   walberla::mpi::MPIManager::instance()->useWorldComm();
+
+   logging::Logging::instance()->setStreamLogLevel(logging::Logging::INFO);
+   logging::Logging::instance()->setFileLogLevel(logging::Logging::INFO);
+
+   WcTimingPool globalTP;
+   WcTimingPool loopTP;
+   WcTimingTree tt;
+
+   WALBERLA_LOG_INFO_ON_ROOT("loading configuration parameters");
+   globalTP["load config"].start();
+   //---------------Generation of CNT configuration-----------------------------
+   FilmSpecimen filmSpecimen; // Data structure for CNT film
+   auto cfg = env.config();
+   if (cfg == nullptr)
+   WALBERLA_ABORT("No config specified!");
+   const Config::BlockHandle cntConf = cfg->getBlock("CNT");
+   loadFromConfig(filmSpecimen, cntConf);
+   //===========================================================================
+   globalTP["load config"].end();
+
+   WALBERLA_LOG_INFO_ON_ROOT("creating domain setup");
+   globalTP["domain setup"].start();
+   Vector3<bool> periodicity = Vector3<bool>(true, true, false);
+   if (filmSpecimen.oopp) periodicity[2] = true;
+   real_t maxZ = filmSpecimen.sizeZ;
+   real_t minZ = 0;
+   if (!filmSpecimen.oopp)
+   {
+      maxZ += real_c(filmSpecimen.numSegs) * filmSpecimen.spacing;
+      minZ -= real_c(filmSpecimen.numSegs) * filmSpecimen.spacing;
+   }
+
+   auto forest = blockforest::createBlockForest(
+         AABB(0, 0, minZ, filmSpecimen.sizeX, filmSpecimen.sizeY, maxZ), // simulation domain
+         Vector3<uint_t>(filmSpecimen.numBlocksX, filmSpecimen.numBlocksY, filmSpecimen.numBlocksZ), // blocks in each direction
+         periodicity // periodicity
+   );
+   if (!forest) WALBERLA_ABORT("No BlockForest created ... exiting!");
+
+   auto domain = std::make_shared<domain::BlockForestDomain> (forest);
+
+   auto simulationDomain = forest->getDomain();
+   WALBERLA_UNUSED(simulationDomain);
+   auto localDomain = forest->begin()->getAABB();
+   for (auto& blk : *forest)
+   {
+      localDomain.merge(blk.getAABB());
+   }
+   globalTP["domain setup"].end();
+
+   WALBERLA_LOG_INFO_ON_ROOT("creating initial particle setup");
+   globalTP["generate CNT"].start();
+   auto ps = std::make_shared<data::ParticleStorage>(10);
+   auto ac = SphericalSegmentAccessor(ps);
+   auto lc = data::SparseLinkedCells(localDomain.getExtended(4_r * kernel::cnt::outer_radius), 4_r * kernel::cnt::outer_radius );
+
+   int64_t numParticles = 0;
+   if (filmSpecimen.initialConfigurationFile.empty())
+   {
+      numParticles = generateCNTs(filmSpecimen, ps, *domain);
+   } else
+   {
+      numParticles = loadCNTs(filmSpecimen.initialConfigurationFile, ps);
+   }
+   walberla::mpi::reduceInplace(numParticles, walberla::mpi::SUM);
+   WALBERLA_LOG_INFO_ON_ROOT("#particles created: " << numParticles);
+   globalTP["generate CNT"].end();
+
+   WALBERLA_LOG_INFO_ON_ROOT("setting up VTK output");
+   auto vtkOutput = make_shared<mesa_pd::vtk::ParticleVtkOutput>(ps);
+   vtkOutput->addOutput<data::SelectParticlePosition>("position");
+   auto vtkWriter = walberla::vtk::createVTKOutput_PointData(vtkOutput,
+                                                             "cnt",
+                                                             1,
+                                                             filmSpecimen.vtkFolder,
+                                                             "simulation_step",
+                                                             true,
+                                                             true,
+                                                             true,
+                                                             filmSpecimen.useMPIIO);
+
+   WALBERLA_LOG_INFO_ON_ROOT("setting up interaction models");
+   kernel::InsertParticleIntoSparseLinkedCells ipilc;
+   kernel::cnt::AnisotropicVDWContact anisotropicVdwContact;
+   kernel::cnt::VBondContact vBondContact;
+   kernel::cnt::ViscousDamping viscous_damping(filmSpecimen.viscousDamping * 1052.0_r,
+                                               filmSpecimen.viscousDamping * 1052.0_r);
+   kernel::PFCDamping pfcDamping(filmSpecimen.localDamping);
+   kernel::VelocityVerletPreForceUpdate  vv_pre(kernel::cnt::dT);
+   kernel::VelocityVerletPostForceUpdate vv_post(kernel::cnt::dT);
+   mpi::ContactFilter                    contactFilter;
+   mpi::ReduceProperty                   RP;
+   mpi::SyncNextNeighbors                SNN;
+
+   WALBERLA_LOG_INFO_ON_ROOT("warmup");
+   SNN(*ps, *domain);
+   sorting::LinearizedCompareFunctor linear(lc.domain_, Vector3<uint_t>(lc.numCellsPerDim_[0]));
+   ps->sort(linear);
+   vtkWriter->write(true);
+
+   WALBERLA_LOG_INFO_ON_ROOT("running simulation");
+   int64_t contactsChecked  = 0;
+   int64_t contactsDetected = 0;
+   int64_t contactsTreated  = 0;
+
+   std::vector<Statistics> statistics;
+
+   globalTP["sim loop"].start();
+   for (auto i = 0; i < filmSpecimen.simulationSteps; ++i)
+   {
+      Statistics stats;
+
+      loopTP["VV pre"].start();
+      ps->forEachParticle(false,
+                          kernel::SelectAll(),
+                          ac,
+                          vv_pre,
+                          ac);
+      loopTP["VV pre"].end();
+
+      loopTP["SNN"].start();
+      SNN(*ps, *domain);
+      loopTP["SNN"].end();
+
+      loopTP["GenerateSparseLinkedCells"].start();
+      lc.clear();
+      ps->forEachParticle(false, kernel::SelectAll(), ac, ipilc, ac, lc);
+      loopTP["GenerateSparseLinkedCells"].end();
+
+      loopTP["interaction"].start();
+      contactsChecked  = 0;
+      contactsDetected = 0;
+      contactsTreated  = 0;
+      lc.forEachParticlePairHalf(false,
+                                 kernel::SelectAll(),
+                                 ac,
+                                 [&](size_t p_idx1, size_t p_idx2)
+                                 {
+                                    //ensure ordering
+                                    if (ac.getUid(p_idx1) > ac.getUid(p_idx2))
+                                       std::swap(ac.getUidRef(p_idx1), ac.getUidRef(p_idx1));
+                                    
+                                    ++contactsChecked;
+                                    auto dist = ac.getPosition(p_idx1) - ac.getPosition(p_idx2);
+                                    auto dist2 = math::dot(dist, dist);
+                                    auto cutoff = ac.getInteractionRadius(p_idx1) + ac.getInteractionRadius(p_idx2);
+                                    if (dist2 < cutoff * cutoff)
+                                    {
+                                       auto contactPoint = ac.getPosition(p_idx2) + 0.5_r * dist;
+                                       ++contactsDetected;
+                                       if (contactFilter(p_idx1, p_idx2, ac, contactPoint, *domain))
+                                       {
+                                          ++contactsTreated;
+                                          constexpr bool period_X = false;
+                                          constexpr bool period_Y = false;
+                                          constexpr int max_segsX = 0;
+                                          constexpr int max_segsY = 0;
+                                          auto sameCluster = ac.getClusterID(p_idx2) == ac.getClusterID(p_idx1);
+                                          auto segmentDistance = std::abs(ac.getSegmentID(p_idx2) - ac.getSegmentID(p_idx1));
+                                          if ((sameCluster) && ((segmentDistance == 1) ||
+                                                                ((period_X) && (segmentDistance == max_segsX - 1)) ||
+                                                                ((period_Y) && (segmentDistance == max_segsY - 1))
+                                          )
+                                                )
+                                          {
+                                             vBondContact(p_idx1, p_idx2, ac);
+                                             stats.bendingEnergy += vBondContact.getLastBendingEnergy();
+                                             stats.shearEnergy += vBondContact.getLastShearEnergy();
+                                             stats.tensileEnergy += vBondContact.getLastTensileEnergy();
+                                             stats.twistingEnergy += vBondContact.getLastTwistingEnergy();
+                                             stats.stretchEnergy += vBondContact.getLastBendingEnergy() +
+                                                   vBondContact.getLastShearEnergy() +
+                                                   vBondContact.getLastTensileEnergy() +
+                                                   vBondContact.getLastTwistingEnergy();
+                                          }
+                                          else if ((sameCluster) && (segmentDistance < 20))
+                                          {
+                                             //do nothing
+                                          } else
+                                          {
+                                             anisotropicVdwContact(p_idx1, p_idx2, ac);
+                                             stats.vdwEnergy += anisotropicVdwContact.getLastEnergy();
+                                             stats.numAlignedPairs += anisotropicVdwContact.isParallel() ? 1 : 0;
+                                          }
+                                          viscous_damping(p_idx1, p_idx2, ac);
+                                       }
+                                    }
+                                 });
+      loopTP["interaction"].end();
+
+      loopTP["ReduceForce"].start();
+      RP.operator()<ForceTorqueNotification>(*ps);
+      loopTP["ReduceForce"].end();
+
+      loopTP["PFDDamping"].start();
+      ps->forEachParticle(false, kernel::SelectAll(), ac, pfcDamping, ac);
+      loopTP["PFDDamping"].end();
+
+      loopTP["VV post"].start();
+      ps->forEachParticle(false,
+                          kernel::SelectAll(),
+                          ac,
+                          vv_post,
+                          ac);
+      loopTP["VV post"].end();
+
+      if ((i + 1) % filmSpecimen.saveVTKEveryNthStep == 0)
+      {
+         loopTP["save vtk"].start();
+         WALBERLA_LOG_INFO_ON_ROOT(CYAN << "            SAVING VTK                  " << RESET);
+         vtkWriter->write(true);
+         loopTP["save vtk"].end();
+      }
+      if ((i + 1) % filmSpecimen.saveEnergyEveryNthStep == 0)
+      {
+         loopTP["save energy"].start();
+         WALBERLA_LOG_INFO_ON_ROOT(CYAN << "            SAVING ENERGY                  " << RESET);
+         ps->forEachParticle(false,
+                             kernel::SelectLocal(),
+                             ac,
+                             [&](const size_t p_idx)
+                             {
+                                stats.kineticEnergy += 0.5_r * ac.getMass(p_idx) * ac.getLinearVelocity(p_idx) * ac.getLinearVelocity(p_idx);
+                                stats.kineticEnergy += 0.5_r * ac.getAngularVelocity(p_idx) * ac.getInertiaBF(p_idx) * ac.getAngularVelocity(p_idx);
+                             });
+         walberla::mpi::reduceInplace(stats, walberla::mpi::SUM);
+         WALBERLA_ROOT_SECTION()
+         {
+            statistics.push_back(stats);
+            saveStatistics(filmSpecimen.energyFolder + "/statistics.txt", statistics);
+         }
+         loopTP["save energy"].end();
+      }
+      if ((i + 1) % filmSpecimen.saveConfEveryNthStep == 0)
+      {
+         loopTP["save config"].start();
+         WALBERLA_LOG_INFO_ON_ROOT(CYAN << "            SAVING CONFIG                  " << RESET);
+         saveConfig(ps,
+                    filmSpecimen.confFolder + "/conf_" + std::to_string((i + 1) * 0.00002) + "ns",
+                    IO_MODE::USE_MPI_IO);
+         loopTP["save config"].end();
+      }
+   }
+   globalTP["sim loop"].end();
+
+   auto loopTPReduced = loopTP.getReduced();
+   WALBERLA_LOG_DEVEL_ON_ROOT(*loopTPReduced);
+   auto globalTPReduced = globalTP.getReduced();
+   WALBERLA_LOG_DEVEL_ON_ROOT(*globalTPReduced);
+   auto ttReduced = tt.getReduced();
+   WALBERLA_LOG_DEVEL_ON_ROOT(ttReduced);
+
+   numParticles = 0;
+   int64_t numGhostParticles = 0;
+   ps->forEachParticle(false,
+                       kernel::SelectAll(),
+                       ac,
+                       [&](const size_t idx)
+                       {
+                          if (data::particle_flags::isSet( ac.getFlagsRef(idx), data::particle_flags::GHOST))
+                          {
+                             ++numGhostParticles;
+                          } else
+                          {
+                             ++numParticles;
+                          }
+                       });
+
+   walberla::mpi::reduceInplace(numParticles, walberla::mpi::SUM);
+   walberla::mpi::reduceInplace(numGhostParticles, walberla::mpi::SUM);
+   walberla::mpi::reduceInplace(contactsDetected, walberla::mpi::SUM);
+   walberla::mpi::reduceInplace(contactsTreated, walberla::mpi::SUM);
+
+
+   WALBERLA_ROOT_SECTION()
+   {
+      std::map<std::string, walberla::int64_t> integerProperties;
+      std::map<std::string, double> realProperties;
+      std::map<std::string, std::string> stringProperties;
+
+      addBuildInfoToSQL(integerProperties, realProperties, stringProperties);
+      addDomainPropertiesToSQL(*forest, integerProperties, realProperties, stringProperties);
+      addSlurmPropertiesToSQL(integerProperties, realProperties, stringProperties);
+
+      integerProperties["num_particles"] = numParticles;
+      integerProperties["num_ghost_particles"] = numGhostParticles;
+
+      integerProperties["contacts_detected"] = contactsDetected;
+      integerProperties["contacts_treated"] = contactsTreated;
+
+      auto runId = sqlite::storeRunInSqliteDB(filmSpecimen.sqlFile,
+                                              integerProperties,
+                                              stringProperties,
+                                              realProperties);
+      sqlite::storeTimingPoolInSqliteDB(filmSpecimen.sqlFile, runId, *loopTPReduced, "loop");
+      sqlite::storeTimingPoolInSqliteDB(filmSpecimen.sqlFile, runId, *globalTPReduced, "global");
+      sqlite::storeTimingTreeInSqliteDB(filmSpecimen.sqlFile, runId, tt, "tt");
+   }
+
+   return EXIT_SUCCESS;
+}
+
+} // namespace mesa_pd
+} // namespace walberla
+
+int main(int argc, char *argv[])
+{
+   return walberla::mesa_pd::main(argc, argv);
+}
+

apps/benchmarks/CNT/CMakeLists.txt

+waLBerla_link_files_to_builddir( *.dat )
+
+waLBerla_add_executable( NAME 01_cnt_film
+   FILES 01_cnt_film.cpp
+        Config.cpp
+        FilmSpecimen.cpp
+        InitializeCNTs.cpp
+        SQLProperties.cpp
+        Statistics.cpp
+      DEPENDS walberla::blockforest walberla::core walberla::mesa_pd walberla::sqlite walberla::vtk )

apps/benchmarks/CNT/Config.cpp

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#include "Config.h"
+#include "TerminalColors.h"
+
+#include "core/mpi/Gatherv.h"
+#include "core/mpi/MPIIO.h"
+#include "core/mpi/RecvBuffer.h"
+#include "core/mpi/SendBuffer.h"
+
+namespace walberla {
+namespace mesa_pd {
+
+void saveConfig(const shared_ptr<data::ParticleStorage>& ps,
+                const std::string& filename,
+                const IO_MODE io_mode)
+{
+   mpi::SendBuffer sb;
+
+   for (const auto& p : *ps)
+   {
+      using namespace walberla::mesa_pd::data::particle_flags;
+      if (isSet(p.getFlags(), GHOST)) continue; //skip ghosts
+
+      sb << static_cast<int64_t> (p.getSegmentID());
+      sb << static_cast<int64_t> (p.getClusterID());
+//      sb << static_cast<int64_t> (p.getGroupID());
+
+      sb << static_cast<double> (p.getPosition()[0]);
+      sb << static_cast<double> (p.getPosition()[1]);
+      sb << static_cast<double> (p.getPosition()[2]);
+
+      Vec3 orient = p.getRotation().getMatrix() * Vec3(1,0,0);
+      // Now orient has the shape (sin(th)cos(ph), sin(th)sin(ph), cos(th))
+      sb << static_cast<double> (std::acos(orient[2]));
+      sb << static_cast<double> (std::atan2(orient[1], orient[0]));
+
+      sb << static_cast<double> (p.getLinearVelocity()[0]);
+      sb << static_cast<double> (p.getLinearVelocity()[1]);
+      sb << static_cast<double> (p.getLinearVelocity()[2]);
+
+      sb << static_cast<double> (p.getAngularVelocity()[0]);
+      sb << static_cast<double> (p.getAngularVelocity()[1]);
+      sb << static_cast<double> (p.getAngularVelocity()[2]);
+   }
+
+   switch (io_mode)
+   {
+      case IO_MODE::REDUCE_TO_ROOT:
+      {
+         WALBERLA_LOG_INFO_ON_ROOT(CYAN << "Saving configuration (reduce to root): " << filename << RESET);
+         walberla::mpi::RecvBuffer rb;
+         walberla::mpi::gathervBuffer(sb, rb);
+         WALBERLA_ROOT_SECTION()
+         {
+            uint_t dataSize = sizeof(mpi::RecvBuffer::ElementType) * rb.size();
+            std::ofstream binfile;
+            binfile.open(filename + ".sav",
+                         std::ios::out | std::ios::binary);
+            binfile.write(reinterpret_cast< const char * >(rb.ptr()), numeric_cast<std::streamsize>(dataSize));
+            binfile.close();
+         }
+         break;
+      }
+      case IO_MODE::USE_MPI_IO:
+      {
+         WALBERLA_LOG_INFO_ON_ROOT(CYAN << "Saving configuration (MPIIO): " << filename << RESET);
+         walberla::mpi::writeMPIIO(filename + ".mpiio", sb);
+         break;
+      }
+      case IO_MODE::ONE_FILE_PER_RANK:
+      {
+         WALBERLA_LOG_INFO_ON_ROOT(CYAN << "Saving configuration (one file per process): " << filename << RESET);
+         uint_t dataSize = sizeof(mpi::SendBuffer::ElementType) * sb.size();
+         std::ofstream binfile;
+         binfile.open(filename + "_" + std::to_string(mpi::MPIManager::instance()->rank()) + ".sav",
+                      std::ios::out | std::ios::binary);
+         binfile.write(reinterpret_cast< const char * >(sb.ptr()), numeric_cast<std::streamsize>(dataSize));
+         binfile.close();
+      }
+         break;
+      default:
+      WALBERLA_ABORT("Unknown IO_MODE");
+   }
+}
+
+} //namespace mesa_pd
+} //namespace walberla

apps/benchmarks/CNT/Config.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "mesa_pd/data/ParticleStorage.h"
+
+namespace walberla {
+namespace mesa_pd {
+
+enum class IO_MODE {REDUCE_TO_ROOT, USE_MPI_IO, ONE_FILE_PER_RANK};
+
+void saveConfig(const shared_ptr<data::ParticleStorage>& ps,
+                const std::string& filename,
+                const IO_MODE io_mode);
+
+} //namespace mesa_pd
+} //namespace walberla
\ No newline at end of file

apps/benchmarks/CNT/FilmSpecimen.cpp

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+//======================================================================================================================
+//
+//  THIS FILE IS GENERATED - PLEASE CHANGE THE TEMPLATE !!!
+//
+//======================================================================================================================
+
+#include "FilmSpecimen.h"
+#include "TerminalColors.h"
+
+#include <core/logging/Logging.h>
+
+namespace walberla {
+namespace mesa_pd {
+
+void loadFromConfig(FilmSpecimen& params, const Config::BlockHandle& cfg)
+{
+   WALBERLA_LOG_INFO_ON_ROOT(GREEN << "Loading CNT film setup: " << RESET)
+   params.sizeX = cfg.getParameter<real_t>("sizeX", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "sizeX: " << GREEN << params.sizeX << RESET);
+   
+   params.sizeY = cfg.getParameter<real_t>("sizeY", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "sizeY: " << GREEN << params.sizeY << RESET);
+   
+   params.sizeZ = cfg.getParameter<real_t>("sizeZ", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "sizeZ: " << GREEN << params.sizeZ << RESET);
+   
+   params.oopp = cfg.getParameter<bool>("oopp", false );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "oopp: " << GREEN << params.oopp << RESET);
+   
+   params.numBlocksX = cfg.getParameter<uint_t>("numBlocksX", 1 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "numBlocksX: " << GREEN << params.numBlocksX << RESET);
+   
+   params.numBlocksY = cfg.getParameter<uint_t>("numBlocksY", 1 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "numBlocksY: " << GREEN << params.numBlocksY << RESET);
+   
+   params.numBlocksZ = cfg.getParameter<uint_t>("numBlocksZ", 1 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "numBlocksZ: " << GREEN << params.numBlocksZ << RESET);
+   
+   params.min_OOP = cfg.getParameter<real_t>("min_OOP", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "min_OOP: " << GREEN << params.min_OOP << RESET);
+   
+   params.max_OOP = cfg.getParameter<real_t>("max_OOP", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "max_OOP: " << GREEN << params.max_OOP << RESET);
+   
+   params.numCNTs = cfg.getParameter<int>("numCNTs", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "numCNTs: " << GREEN << params.numCNTs << RESET);
+   
+   params.numSegs = cfg.getParameter<int>("numSegs", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "numSegs: " << GREEN << params.numSegs << RESET);
+   
+   params.spacing = cfg.getParameter<real_t>("spacing", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "spacing: " << GREEN << params.spacing << RESET);
+   
+   params.localDamping = cfg.getParameter<real_t>("localDamping", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "localDamping: " << GREEN << params.localDamping << RESET);
+   
+   params.viscousDamping = cfg.getParameter<real_t>("viscousDamping", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "viscousDamping: " << GREEN << params.viscousDamping << RESET);
+   
+   params.seed = cfg.getParameter<uint_t>("seed", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "seed: " << GREEN << params.seed << RESET);
+   
+   params.vdW = cfg.getParameter<int>("vdW", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "vdW: " << GREEN << params.vdW << RESET);
+   
+   params.simulationSteps = cfg.getParameter<int>("simulationSteps", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "simulationSteps: " << GREEN << params.simulationSteps << RESET);
+   
+   params.saveVTKEveryNthStep = cfg.getParameter<int>("saveVTKEveryNthStep", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "saveVTKEveryNthStep: " << GREEN << params.saveVTKEveryNthStep << RESET);
+   
+   params.saveEnergyEveryNthStep = cfg.getParameter<int>("saveEnergyEveryNthStep", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "saveEnergyEveryNthStep: " << GREEN << params.saveEnergyEveryNthStep << RESET);
+   
+   params.saveConfEveryNthStep = cfg.getParameter<int>("saveConfEveryNthStep", 0 );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "saveConfEveryNthStep: " << GREEN << params.saveConfEveryNthStep << RESET);
+   
+   params.vtkFolder = cfg.getParameter<std::string>("vtkFolder", "." );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "vtkFolder: " << GREEN << params.vtkFolder << RESET);
+   
+   params.energyFolder = cfg.getParameter<std::string>("energyFolder", "." );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "energyFolder: " << GREEN << params.energyFolder << RESET);
+   
+   params.confFolder = cfg.getParameter<std::string>("confFolder", "." );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "confFolder: " << GREEN << params.confFolder << RESET);
+   
+   params.useMPIIO = cfg.getParameter<bool>("useMPIIO", false );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "useMPIIO: " << GREEN << params.useMPIIO << RESET);
+   
+   params.sqlFile = cfg.getParameter<std::string>("sqlFile", "cnt.sqlite" );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "sqlFile: " << GREEN << params.sqlFile << RESET);
+   
+   params.initialConfigurationFile = cfg.getParameter<std::string>("initialConfigurationFile", "" );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "initialConfigurationFile: " << GREEN << params.initialConfigurationFile << RESET);
+   
+}
+
+} //namespace mesa_pd
+} //namespace walberla
+

apps/benchmarks/CNT/FilmSpecimen.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+//======================================================================================================================
+//
+//  THIS FILE IS GENERATED - PLEASE CHANGE THE TEMPLATE !!!
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "core/config/Config.h"
+
+#include <string>
+
+namespace walberla {
+namespace mesa_pd {
+
+struct FilmSpecimen
+{
+   real_t sizeX = 0; //Specimen length (x direction)
+   real_t sizeY = 0; //Specimen length (y direction)
+   real_t sizeZ = 0; //Specimen length (z direction)
+   bool oopp = false; //Specimen out-of-plane periodicity (0 - film, 1 - cnt material)
+   uint_t numBlocksX = 1; //Number of blocks in x direction
+   uint_t numBlocksY = 1; //Number of blocks in y direction
+   uint_t numBlocksZ = 1; //Number of blocks in z direction
+   real_t min_OOP = 0; //Out-of-plane angle minimum
+   real_t max_OOP = 0; //Out-of-plane angle maximum
+   int numCNTs = 0; //Number of CNTs
+   int numSegs = 0; //Number of segments in a CNT
+   real_t spacing = 0; //Segment half-spacing
+   real_t localDamping = 0; //Local damping coefficient
+   real_t viscousDamping = 0; //Viscous damping coefficient
+   uint_t seed = 0; //random generator seed
+   int vdW = 0; //type of vdW interaction model
+   int simulationSteps = 0; //Relaxation duration
+   int saveVTKEveryNthStep = 0; //timesteps between saving VTK outputs
+   int saveEnergyEveryNthStep = 0; //timesteps between saving energies
+   int saveConfEveryNthStep = 0; //timesteps between saving confs
+   std::string vtkFolder = "."; //Folder for VTK files
+   std::string energyFolder = "."; //Folder for energy files
+   std::string confFolder = "."; //Folder for conf files
+   bool useMPIIO = false; //Write a single file instead of one file per process
+   std::string sqlFile = "cnt.sqlite"; //database file
+   std::string initialConfigurationFile = ""; //restart from checkpoint
+};
+
+void loadFromConfig(FilmSpecimen& params,
+                    const Config::BlockHandle& cfg);
+
+} //namespace mesa_pd
+} //namespace walberla
\ No newline at end of file

apps/benchmarks/CNT/FilmSpecimen.templ.cpp

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+//======================================================================================================================
+//
+//  THIS FILE IS GENERATED - PLEASE CHANGE THE TEMPLATE !!!
+//
+//======================================================================================================================
+
+#include "FilmSpecimen.h"
+#include "TerminalColors.h"
+
+#include <core/logging/Logging.h>
+
+namespace walberla {
+namespace mesa_pd {
+
+void loadFromConfig(FilmSpecimen& params, const Config::BlockHandle& cfg)
+{
+   WALBERLA_LOG_INFO_ON_ROOT(GREEN << "Loading CNT film setup: " << RESET)
+   {%- for param in parameters %}
+   params.{{param.name}} = cfg.getParameter<{{param.type}}>("{{param.name}}", {{param.defValue}} );
+   WALBERLA_LOG_INFO_ON_ROOT(YELLOW << "{{param.name}}: " << GREEN << params.{{param.name}} << RESET);
+   {% endfor %}
+}
+
+} //namespace mesa_pd
+} //namespace walberla
+
+

apps/benchmarks/CNT/FilmSpecimen.templ.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+//======================================================================================================================
+//
+//  THIS FILE IS GENERATED - PLEASE CHANGE THE TEMPLATE !!!
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "core/config/Config.h"
+
+#include <string>
+
+namespace walberla {
+namespace mesa_pd {
+
+struct FilmSpecimen
+{
+   {%- for param in parameters %}
+   {{param.type}} {{param.name}} = {{param.defValue}}; //{{param.comment}}
+   {%- endfor %}
+};
+
+void loadFromConfig(FilmSpecimen& params,
+                    const Config::BlockHandle& cfg);
+
+} //namespace mesa_pd
+} //namespace walberla

apps/benchmarks/CNT/FilmSpecimenGenerator.py

+# -*- coding: utf-8 -*-
+
+from jinja2 import Environment, FileSystemLoader
+import os
+
+
+class Parameter:
+    def __init__(self, name, type, defValue="", comment=""):
+        """Property of a data structure
+
+        Parameters
+        ----------
+        name : str
+           name of the property
+        type : str
+           type of the property
+        defValue : str
+           default value the property should be initialized with
+        comment : str
+           comment appended to the properties in the header file
+        """
+
+        self.name = name
+        self.type = type
+        self.defValue = defValue
+        self.comment = comment
+
+    def __str__(self):
+        return "name: {}, type: {}, defValue: {} \\\\{}".format(self.name, self.type, self.defValue, self.comment)
+
+
+class Config:
+    def __init__(self):
+        self.parameters = []
+
+    def parameter_exists(self, name):
+        for v in self.parameters:
+            if v.name == name:
+                return True
+        return False
+
+    def add_parameter(self, name, type, defValue, comment):
+        if self.parameter_exists(name):
+            print("parameters already added: " + name)
+        else:
+            self.parameters.append(Parameter(name, type, defValue, comment))
+
+    def generate_file(self, template):
+        context = dict()
+        context["parameters"] = self.parameters
+
+        path = ""
+        filename = template.replace(".templ", "")
+        dirname = os.path.dirname(__file__)
+        env = Environment(loader=FileSystemLoader(dirname))
+        print("generating: " + path + filename)
+        with open(path + filename, "wb") as fout:
+            content = env.get_template(template).render(context)
+            fout.write(content.encode('utf8'))
+
+    def generate(self):
+        print("=" * 90)
+        print("Config File:")
+        print("")
+        print("{0: <30}{1: <30}{2: <30}".format("Name", "Type", "Def. Value"))
+        print("=" * 90)
+        for param in self.parameters:
+            print("{0: <30.29}{1: <30.29}{2: <30.29}".format(param.name, param.type, param.defValue))
+        print("=" * 90)
+
+        self.generate_file("FilmSpecimen.templ.h")
+        self.generate_file("FilmSpecimen.templ.cpp")
+
+

apps/benchmarks/CNT/GenerateFilmSpecimen.py

+#! /usr/bin/env python3
+# -*- coding: utf-8 -*-
+
+from FilmSpecimenGenerator import Config
+
+cfg = Config()
+cfg.add_parameter("sizeX",                    "real_t",     "0",     "Specimen length (x direction)")
+cfg.add_parameter("sizeY",                    "real_t",     "0",     "Specimen length (y direction)")
+cfg.add_parameter("sizeZ",                    "real_t",     "0",     "Specimen length (z direction)")
+cfg.add_parameter("oopp",                     "bool",       "false", "Specimen out-of-plane periodicity (0 - film, 1 - cnt material)")
+cfg.add_parameter("numBlocksX",               "uint_t",     "1",     "Number of blocks in x direction")
+cfg.add_parameter("numBlocksY",               "uint_t",     "1",     "Number of blocks in y direction")
+cfg.add_parameter("numBlocksZ",               "uint_t",     "1",     "Number of blocks in z direction")
+cfg.add_parameter("min_OOP",                  "real_t",     "0",     "Out-of-plane angle minimum")
+cfg.add_parameter("max_OOP",                  "real_t",     "0",     "Out-of-plane angle maximum")
+cfg.add_parameter("numCNTs",                  "int",        "0",     "Number of CNTs")
+cfg.add_parameter("numSegs",                  "int",        "0",     "Number of segments in a CNT")
+cfg.add_parameter("spacing",                  "real_t",     "0",     "Segment half-spacing")
+cfg.add_parameter("localDamping",             "real_t",     "0",     "Local damping coefficient")
+cfg.add_parameter("viscousDamping",           "real_t",     "0",     "Viscous damping coefficient")
+cfg.add_parameter("seed",                     "uint_t",     "0",     "random generator seed")
+cfg.add_parameter("vdW",                      "int",        "0",     "type of vdW interaction model")
+cfg.add_parameter("simulationSteps",          "int",        "0",     "Relaxation duration")
+cfg.add_parameter("saveVTKEveryNthStep",      "int",        "0",     "timesteps between saving VTK outputs")
+cfg.add_parameter("saveEnergyEveryNthStep",   "int",        "0",     "timesteps between saving energies")
+cfg.add_parameter("saveConfEveryNthStep",     "int",        "0",     "timesteps between saving confs")
+cfg.add_parameter("vtkFolder",                "std::string", '"."',   "Folder for VTK files")
+cfg.add_parameter("energyFolder",             "std::string", '"."',   "Folder for energy files")
+cfg.add_parameter("confFolder",               "std::string", '"."',   "Folder for conf files")
+cfg.add_parameter("useMPIIO",                 "bool",       "false", "Write a single file instead of one file per process")
+cfg.add_parameter("sqlFile",                  "std::string", '"cnt.sqlite"',   "database file")
+cfg.add_parameter("initialConfigurationFile", "std::string", '""',   "restart from checkpoint")
+
+cfg.generate()
+
+

apps/benchmarks/CNT/InitializeCNTs.cpp

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#include "InitializeCNTs.h"
+#include "TerminalColors.h"
+
+#include "core/math/Constants.h"
+#include "core/math/Random.h"
+#include "core/mpi/MPIManager.h"
+#include "mesa_pd/kernel/cnt/Parameters.h"
+
+#include <cmath>
+
+namespace walberla{
+namespace mesa_pd {
+
+int64_t generateCNTs(const FilmSpecimen& spec,
+                     const shared_ptr<data::ParticleStorage>& ps,
+                     const domain::IDomain& domain)
+{
+   auto myRank = mpi::MPIManager::instance()->rank();
+
+   auto CNT_length = 2_r * spec.spacing * real_c(spec.numSegs);
+   // Fixed random seed is necessary for coordinated generation on all MPI proc.
+   auto rand0_1 = math::RealRandom<real_t>(static_cast<std::mt19937::result_type>(spec.seed));
+   // Create an assembly of CNTs
+   int64_t numParticles = 0;
+   for (int tube = 0; tube < spec.numCNTs; tube++)
+   {
+      // This definition of theta provides uniform distribution of random orientations on a sphere when spec.OutOfPlane = 1.
+      real_t theta = 0.5_r * math::pi * spec.min_OOP + (spec.max_OOP - spec.min_OOP) * std::acos(1_r - rand0_1());
+      real_t phi = 2_r * math::pi * rand0_1();
+      real_t pos_x = spec.sizeX * rand0_1();
+      real_t pos_y = spec.sizeY * rand0_1();
+      real_t pos_z = spec.sizeZ * rand0_1();
+
+      for (int segment = 0; segment < spec.numSegs; segment++)
+      {
+         // Generate segment position
+         Vec3 pos;
+         auto r = -0.5_r * CNT_length + 2_r * spec.spacing * real_c(segment);
+         pos[0] = pos_x + r * std::sin(theta) * std::cos(phi);
+         pos[1] = pos_y + r * std::sin(theta) * std::sin(phi);
+         pos[2] = pos_z + r * std::cos(theta);
+
+         // Account for periodicity (will not work if CNTs are longer than domain size)
+         if ((pos[0] > spec.sizeX)) pos[0] -= spec.sizeX;
+         if ((pos[0] < 0)) pos[0] += spec.sizeX;
+         if ((pos[1] > spec.sizeY)) pos[1] -= spec.sizeY;
+         if ((pos[1] < 0)) pos[1] += spec.sizeY;
+         if (spec.oopp)
+         {
+            if ((pos[2] > spec.sizeZ)) pos[2] -= spec.sizeZ;
+            if ((pos[2] < 0)) pos[2] += spec.sizeZ;
+         }
+
+         //Only create segment if it is located on the process.
+         if (domain.isContainedInProcessSubdomain(uint_c(myRank), pos))
+         {
+            data::Particle &&sp = *ps->create();
+            sp.setPosition(pos);
+            sp.setOwner(myRank);
+            sp.setInteractionRadius(kernel::cnt::outer_radius);
+            sp.setSegmentID(segment);
+            sp.setClusterID(tube);
+            sp.getRotationRef().rotate( Vec3(0_r, 1_r, 0_r), -0.5_r * math::pi + theta);
+            sp.getRotationRef().rotate( Vec3(0_r, 0_r, 1_r), phi);
+            numParticles++;
+         }
+      }
+   }
+
+   return numParticles;
+}
+
+int64_t loadCNTs(const std::string& filename,
+                 const shared_ptr<data::ParticleStorage>& ps)
+{
+   WALBERLA_LOG_INFO_ON_ROOT(GREEN << "Loading configuration (binary format): " << filename);
+   const auto numProcesses = mpi::MPIManager::instance()->numProcesses();
+   const auto rank = mpi::MPIManager::instance()->rank();
+   //---------Generation of WaLBerla model---------------------------------------------------------------------------------------------------------------------------------
+   // Note that walberla primitives are created below only on MPI process that is responsible
+   // for this branch of blockforest
+
+   int64_t numParticles = 0;
+
+   for (auto i = 0; i < numProcesses; ++i)
+   {
+      WALBERLA_MPI_BARRIER();
+      if (i != rank) continue; //bad practice but with the current file format we have to do this
+      std::ifstream binfile;
+      binfile.open(filename.c_str(), std::ios::in | std::ios::binary);
+      int size;
+      binfile.read((char *) &size, sizeof(int));
+      std::cout << RED << "size read form binary file is" << size << RESET << std::endl;
+      for (int id = 0; id < size; id++)
+      {
+         int ID;
+         int sID;
+         int cID;
+         int gID;
+         double x;
+         double y;
+         double z;
+         double theta;
+         double phi;
+         double vx;
+         double vy;
+         double vz;
+         double wx;
+         double wy;
+         double wz;
+         binfile.read((char *) &ID, sizeof(int));
+         binfile.read((char *) &sID, sizeof(int));
+         binfile.read((char *) &cID, sizeof(int));
+         binfile.read((char *) &gID, sizeof(int));
+         binfile.read((char *) &x, sizeof(double));
+         binfile.read((char *) &y, sizeof(double));
+         binfile.read((char *) &z, sizeof(double));
+         binfile.read((char *) &theta, sizeof(double));
+         binfile.read((char *) &phi, sizeof(double));
+         binfile.read((char *) &vx, sizeof(double));
+         binfile.read((char *) &vy, sizeof(double));
+         binfile.read((char *) &vz, sizeof(double));
+         binfile.read((char *) &wx, sizeof(double));
+         binfile.read((char *) &wy, sizeof(double));
+         binfile.read((char *) &wz, sizeof(double));
+         Vec3 pos;
+         pos[0] = real_c(x);
+         pos[1] = real_c(y);
+         pos[2] = real_c(z);
+
+         data::Particle &&sp = *ps->create();
+         sp.setPosition(pos);
+         sp.setOwner(rank);
+         sp.setInteractionRadius(kernel::cnt::outer_radius);
+         sp.setSegmentID(sID);
+         sp.setClusterID(cID);
+         sp.getRotationRef().rotate( Vec3(0_r, 1_r, 0_r), -0.5_r * math::pi + real_c(theta));
+         sp.getRotationRef().rotate( Vec3(0_r, 0_r, 1_r), real_c(phi));
+         sp.setLinearVelocity( Vec3(real_c(vx), real_c(vy), real_c(vz)) );
+         sp.setAngularVelocity( Vec3(real_c(wx), real_c(wy), real_c(wz)) );
+         numParticles++;
+      }
+   }
+
+   mpi::reduceInplace(numParticles, mpi::SUM);
+   WALBERLA_LOG_INFO_ON_ROOT("#particles created: " << numParticles);
+
+   return numParticles;
+}
+
+} //namespace mesa_pd
+} //namespace walberla

apps/benchmarks/CNT/InitializeCNTs.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "FilmSpecimen.h"
+
+#include "mesa_pd/data/ParticleStorage.h"
+#include "mesa_pd/domain/IDomain.h"
+
+namespace walberla{
+namespace mesa_pd {
+
+int64_t generateCNTs(const FilmSpecimen& spec,
+                     const shared_ptr<data::ParticleStorage>& ps,
+                     const domain::IDomain& domain);
+
+int64_t loadCNTs(const std::string& filename,
+                 const shared_ptr<data::ParticleStorage>& ps);
+
+} //namespace mesa_pd
+} //namespace walberla

apps/benchmarks/CNT/SQLProperties.cpp

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#include "SQLProperties.h"
+
+#include <blockforest/BlockForest.h>
+#include <core/DataTypes.h>
+#include <core/waLBerlaBuildInfo.h>
+#include <core/mpi/MPIManager.h>
+
+namespace walberla {
+namespace mesa_pd {
+
+void addBuildInfoToSQL(std::map<std::string, int64_t> &       /*integerProperties*/,
+                       std::map<std::string, double> &        /*realProperties*/,
+                       std::map<std::string, std::string> &stringProperties)
+{
+   stringProperties["walberla_sha"] = walberla::core::buildinfo::gitSHA1();
+   stringProperties["build_type"] = walberla::core::buildinfo::buildType();
+   stringProperties["compiler_flags"] = walberla::core::buildinfo::compilerFlags();
+   stringProperties["build_machine"] = walberla::core::buildinfo::buildMachine();
+}
+
+void addDomainPropertiesToSQL(const ::walberla::blockforest::BlockForest &forest,
+                              std::map<std::string, int64_t> &integerProperties,
+                              std::map<std::string, double> &realProperties,
+                              std::map<std::string, std::string> &    /*stringProperties*/ )
+{
+   integerProperties["numMPIProcesses"] = ::walberla::mpi::MPIManager::instance()->numProcesses();
+
+   realProperties["domainXMin"] = forest.getDomain().xMin();
+   realProperties["domainXMax"] = forest.getDomain().xMax();
+   realProperties["domainYMin"] = forest.getDomain().yMin();
+   realProperties["domainYMax"] = forest.getDomain().yMax();
+   realProperties["domainZMin"] = forest.getDomain().zMin();
+   realProperties["domainZMax"] = forest.getDomain().zMax();
+
+   integerProperties["xBlocks"] = ::walberla::int_c(forest.getXSize());
+   integerProperties["yBlocks"] = ::walberla::int_c(forest.getYSize());
+   integerProperties["zBlocks"] = ::walberla::int_c(forest.getZSize());
+
+   integerProperties["xPeriodic"] = (forest.isXPeriodic() ? 1 : 0);
+   integerProperties["yPeriodic"] = (forest.isYPeriodic() ? 1 : 0);
+   integerProperties["zPeriodic"] = (forest.isZPeriodic() ? 1 : 0);
+}
+
+std::string envToString(const char *env)
+{
+   return env != nullptr ? std::string(env) : "";
+}
+
+void addSlurmPropertiesToSQL(std::map<std::string, int64_t> &        /*integerProperties*/,
+                             std::map<std::string, double> &         /*realProperties*/,
+                             std::map<std::string, std::string> &stringProperties)
+{
+   stringProperties["SLURM_CLUSTER_NAME"] = envToString(std::getenv("SLURM_CLUSTER_NAME"));
+   stringProperties["SLURM_CPUS_ON_NODE"] = envToString(std::getenv("SLURM_CPUS_ON_NODE"));
+   stringProperties["SLURM_CPUS_PER_TASK"] = envToString(std::getenv("SLURM_CPUS_PER_TASK"));
+   stringProperties["SLURM_JOB_ACCOUNT"] = envToString(std::getenv("SLURM_JOB_ACCOUNT"));
+   stringProperties["SLURM_JOB_ID"] = envToString(std::getenv("SLURM_JOB_ID"));
+   stringProperties["SLURM_JOB_CPUS_PER_NODE"] = envToString(std::getenv("SLURM_JOB_CPUS_PER_NODE"));
+   stringProperties["SLURM_JOB_NAME"] = envToString(std::getenv("SLURM_JOB_NAME"));
+   stringProperties["SLURM_JOB_NUM_NODES"] = envToString(std::getenv("SLURM_JOB_NUM_NODES"));
+   stringProperties["SLURM_NTASKS"] = envToString(std::getenv("SLURM_NTASKS"));
+   stringProperties["SLURM_NTASKS_PER_CORE"] = envToString(std::getenv("SLURM_NTASKS_PER_CORE"));
+   stringProperties["SLURM_NTASKS_PER_NODE"] = envToString(std::getenv("SLURM_NTASKS_PER_NODE"));
+   stringProperties["SLURM_NTASKS_PER_SOCKET"] = envToString(std::getenv("SLURM_NTASKS_PER_SOCKET"));
+   stringProperties["SLURM_CPU_BIND_TYPE"] = envToString(std::getenv("SLURM_CPU_BIND_TYPE"));
+}
+
+} //namespace mesa_pd
+} //namespace walberla
+

apps/benchmarks/CNT/SQLProperties.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#pragma once
+
+#include <blockforest/BlockForest.h>
+#include <core/DataTypes.h>
+
+#include <map>
+#include <string>
+
+namespace walberla {
+namespace mesa_pd {
+
+void addBuildInfoToSQL(std::map<std::string, int64_t> &integerProperties,
+                       std::map<std::string, double> &realProperties,
+                       std::map<std::string, std::string> &stringProperties);
+
+void addDomainPropertiesToSQL(const walberla::blockforest::BlockForest &forest,
+                              std::map<std::string, int64_t> &integerProperties,
+                              std::map<std::string, double> &realProperties,
+                              std::map<std::string, std::string> & /*stringProperties*/ );
+
+void addSlurmPropertiesToSQL(std::map<std::string, int64_t> &integerProperties,
+                             std::map<std::string, double> &realProperties,
+                             std::map<std::string, std::string> &stringProperties);
+
+} //namespace mesa_pd
+} //namespace walberla

apps/benchmarks/CNT/SphericalSegmentAccessor.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "mesa_pd/data/ParticleAccessor.h"
+#include "mesa_pd/kernel/cnt/Parameters.h"
+
+namespace walberla {
+namespace mesa_pd {
+
+class SphericalSegmentAccessor : public data::ParticleAccessor
+{
+public:
+   SphericalSegmentAccessor(std::shared_ptr<data::ParticleStorage>& ps)
+         : ParticleAccessor(ps)
+   {}
+
+   constexpr auto getMass(const size_t /*p_idx*/) const {return kernel::cnt::mass_T;}
+   walberla::real_t& getMassRef(const size_t /*p_idx*/) {WALBERLA_ABORT("not available");}
+   void setMass(const size_t /*p_idx*/, walberla::real_t const & /*v*/) {WALBERLA_ABORT("not available");}
+
+   constexpr auto getInvMass(const size_t /*p_idx*/) const {return 1_r / kernel::cnt::mass_T;}
+   walberla::real_t& getInvMassRef(const size_t /*p_idx*/) {WALBERLA_ABORT("not available");}
+   void setInvMass(const size_t /*p_idx*/, walberla::real_t const & /*v*/) {WALBERLA_ABORT("not available");}
+
+   constexpr auto& getInertiaBF(const size_t /*p_idx*/) const {return Ia;}
+   Mat3& getInertiaBFRef(const size_t /*p_idx*/) {WALBERLA_ABORT("not available");}
+   void setInertiaBF(const size_t /*p_idx*/, Mat3 const & /*v*/) {WALBERLA_ABORT("not available");}
+
+   constexpr auto& getInvInertiaBF(const size_t /*p_idx*/) const {return invI;}
+   Mat3& getInvInertiaBFRef(const size_t /*p_idx*/) {WALBERLA_ABORT("not available");}
+   void setInvInertiaBF(const size_t /*p_idx*/, Mat3 const & /*v*/) {WALBERLA_ABORT("not available");}
+private:
+   //  - sphere   :  I = (2/5)*mass*radius^2
+   static constexpr auto Ia = 0.4_r * kernel::cnt::mass_T * kernel::cnt::inner_radius * kernel::cnt::inner_radius;
+   static constexpr auto invI = Mat3(1_r/Ia, 0_r, 0_r, 0_r, 1_r/Ia, 0_r, 0_r, 0_r, 1_r/Ia);
+};
+
+} //namespace mesa_pd
+} //namespace walberla

src/pe/collision/GJKEPAHelper.cpp → apps/benchmarks/CNT/Statistics.cpp

@@ -13,65 +13,58 @@
 //  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 GJKHelper.cpp
+//! \file
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
 //! \author Sebastian Eibl <sebastian.eibl@fau.de>
 //
 //======================================================================================================================
 
-//*************************************************************************************************
-// Includes
-//*************************************************************************************************
-
-#include "GJKEPAHelper.h"
-
-#include "pe/rigidbody/GeomPrimitive.h"
-
-extern "C" {
-   #include "pe/extern/libccd/ccd/ccd.h"
-   #include "pe/extern/libccd/ccd/quat.h"
-}
+#include "Statistics.h"
+#include "TerminalColors.h"
 
 #include "core/logging/Logging.h"
 
-namespace walberla {
-namespace pe {
+#include <fstream>
 
-Vec3       convertVec3(const ccd_vec3_t& vec) { return Vec3(vec.v[0], vec.v[1], vec.v[2]); }
-ccd_vec3_t convertVec3(const Vec3& vec)       { ccd_vec3_t ret; ret.v[0] = vec[0]; ret.v[1] = vec[1]; ret.v[2] = vec[2]; return ret; }
+namespace walberla {
+namespace mesa_pd {
 
-void support(const void *obj, const ccd_vec3_t *dir, ccd_vec3_t *vec)
+void saveStatistics(const std::string& filename,
+                    const std::vector<Statistics>& stats)
 {
-    ConstGeomID bd = reinterpret_cast<ConstGeomID> (obj);
-    Vec3 d = convertVec3(*dir);
-    Vec3 sup = bd->support( d );
-    *vec = convertVec3(sup);
+   WALBERLA_LOG_INFO_ON_ROOT(CYAN << "saving statistics: " << filename << RESET);
+   std::ofstream fout(filename);
+   fout.precision(10);
+   fout << "# kineticEnergy vdwEnergy stretchEnergy tensileEnergy shearEnergy bendingEnergy twistingEnergy numAlignedPairs" << std::endl;
+   for (const auto& stat : stats)
+   {
+      fout << stat.kineticEnergy << " "
+           << stat.vdwEnergy << " "
+           << stat.stretchEnergy << " "
+           << stat.tensileEnergy << " "
+           << stat.shearEnergy << " "
+           << stat.bendingEnergy << " "
+           << stat.twistingEnergy << " "
+           << stat.numAlignedPairs << std::endl;
+   }
+   fout.close();
 }
 
-bool collideGJK( ConstGeomID bd1,
-                 ConstGeomID bd2,
-                 Vec3& contactPoint,
-                 Vec3& contactNormal,
-                 real_t& penetrationDepth,
-                 const unsigned long numIterations,
-                 const real_t epaTol )
-{
-    ccd_t ccd;
-    CCD_INIT(&ccd); // initialize ccd_t struct
-
-    // set up ccd_t struct
-    ccd.support1       = support;       // support function for first object
-    ccd.support2       = support;       // support function for second object
-    ccd.max_iterations = numIterations; // maximal number of iterations
-    ccd.epa_tolerance  = epaTol;
+} //namespace mesa_pd
 
-    ccd_vec3_t dir, pos;
-    int intersect = ccdGJKPenetration(reinterpret_cast<const void*> (bd1), reinterpret_cast<const void*> (bd2), &ccd, &penetrationDepth, &dir, &pos);
-    contactPoint  = convertVec3(pos);
-    contactNormal = -convertVec3(dir);
-    penetrationDepth *= -1;
-
-    return (intersect == 0);
+namespace mpi {
+void reduceInplace( mesa_pd::Statistics& stats, mpi::Operation operation, int recvRank, MPI_Comm comm )
+{
+   reduceInplace(stats.kineticEnergy, operation, recvRank, comm);
+   reduceInplace(stats.vdwEnergy, operation, recvRank, comm);
+   reduceInplace(stats.stretchEnergy, operation, recvRank, comm);
+   reduceInplace(stats.tensileEnergy, operation, recvRank, comm);
+   reduceInplace(stats.shearEnergy, operation, recvRank, comm);
+   reduceInplace(stats.bendingEnergy, operation, recvRank, comm);
+   reduceInplace(stats.twistingEnergy, operation, recvRank, comm);
+   reduceInplace(stats.numAlignedPairs, operation, recvRank, comm);
 }
+} //namespace mpi
 
-} // namespace pe
-} // namespace walberla
+} //namespace walberla
\ No newline at end of file

apps/benchmarks/CNT/Statistics.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "core/mpi/Reduce.h"
+
+#include <string>
+#include <vector>
+
+namespace walberla {
+namespace mesa_pd {
+
+struct Statistics
+{
+   real_t kineticEnergy   = 0_r;
+   real_t vdwEnergy       = 0_r;
+   real_t stretchEnergy   = 0_r;
+   real_t tensileEnergy   = 0_r;
+   real_t shearEnergy     = 0_r;
+   real_t bendingEnergy   = 0_r;
+   real_t twistingEnergy  = 0_r;
+   int64_t numAlignedPairs = 0;
+};
+
+void saveStatistics(const std::string& filename,
+                    const std::vector<Statistics>& stats);
+
+} //namespace mesa_pd
+
+namespace mpi {
+void reduceInplace( mesa_pd::Statistics& stats, mpi::Operation operation, int recvRank = 0, MPI_Comm comm = MPI_COMM_WORLD );
+} //namespace mpi
+
+} //namespace walberla

apps/benchmarks/CNT/TerminalColors.h

+//======================================================================================================================
+//
+//  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
+//! \author Igor Ostanin <i.ostanin@skoltech.ru>
+//! \author Grigorii Drozdov <drozd013@umn.edu>
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//
+//======================================================================================================================
+
+#pragma once
+
+#include <string>
+
+namespace walberla {
+
+#ifdef PLAIN_OUTPUT
+//----Terminal color manipulators-----
+const std::string RED("");      // errors and problems
+const std::string GREEN("");    // loading operations
+const std::string YELLOW("");   // information (non-Walberla)
+const std::string PURPLE("");   // stages of program implementation
+const std::string CYAN("");     // save/write operations
+const std::string RESET("");
+//====================================
+#else
+//----Terminal color manipulators-----
+const std::string RED("\033[1;31m");      // errors and problems
+const std::string GREEN("\033[1;32m");    // loading operations
+const std::string YELLOW("\033[1;33m");   // information (non-Walberla)
+const std::string PURPLE("\033[1;35m");   // stages of program implementation
+const std::string CYAN("\033[1;36m");     // save/write operations
+const std::string RESET("\033[0m");
+//====================================
+#endif
+
+} //namespace walberla
\ No newline at end of file

apps/benchmarks/CNT/specimen_setup.dat

+//--------------CNT material generator input file (do not change anything except numerical values!)-
+// This  file  sets  up  generation of a rectangular specimen of a CNT material
+// At the generation stage, triple periodic boundary conditions are imposed.
+// Out-of-plane (xy) angle span tunes transverse anizotropy of the film (1 means isotropy)
+// Note that the numerical values should be single space separated
+
+CNT
+{
+   sizeX 20000; // A - Specimen length (x direction).
+   sizeY 20000; // A - Specimen width (y direction).
+   sizeZ 200; // A - Specimen height (z direction).
+   oopp 0; // - out of plane periodicity
+   numBlocksX  2; // - Number of blocks in x direction (blocks are distributed over available MPI processes,
+   numBlocksY  2; // - Number of blocks in y direction  nx x ny x nz blocks total
+   numBlocksZ  1; // - Number of blocks in z direction (nx, ny, nz can not be less than 2 )
+   min_OOP 0.995; // (in PI/2 rads) - minimum theta angle (set 0 for isotropy)
+   max_OOP 1.0; // (in PI/2 rads) - maximum theta angle (set 1 for isotropy)
+   numCNTs 1000; // - Number of CNTs.
+   numSegs 800; // - Number of segments in a CNT.
+   spacing 6.78; // - Segment spacing
+   localDamping 0.2; // - "local" damping
+   viscousDamping 0.00; // - viscous damping
+   seed 1201; // - random seed - needs to be changed for a different realization
+   vdW 1; // - vdW model (0 - isotr.vdW (comp speed-1.0), 1 - aniz. vdW (comp speed 0.54), 2 - numer.integrated vdW (comp-speed 0.14))
+   simulationSteps 1000000; // - timesteps relaxation duration
+   saveEveryNthStep 10000; // - timesteps between checkpoints
+   saveConfEveryNthStep 100000; // - timesteps between saving confs
+   file serie_3_spec_1; // - filename to save the resulting conf (saved to ../confs/)
+}
+//==================================================================================================

apps/benchmarks/ComplexGeometry/CMakeLists.txt

+if ( WALBERLA_BUILD_WITH_OPENMESH )
+   
+	waLBerla_link_geometry_to_builddir( "*.obj" )             
+	waLBerla_link_files_to_builddir( "*.conf" )  
+
+   waLBerla_add_executable( NAME ComplexGeometry FILES ComplexGeometry.cpp DEPENDS walberla::boundary walberla::core walberla::lbm walberla::mesh walberla::vtk )
+
+	##############
+	# Some tests #
+	##############
+
+   waLBerla_execute_test( NO_MODULE_LABEL NAME ComplexGeometry LABELS longrun COMMAND $<TARGET_FILE:ComplexGeometry> test.conf DEPENDS_ON_TARGETS ComplexGeometry )
+endif()
+