Commit a5400d89 authored by Sebastian Eibl's avatar Sebastian Eibl
Browse files

reworked benchmark structure

parent a40b3c5c
//======================================================================================================================
//
// 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 Accessor.h
//! \author Sebastian Eibl <sebastian.eibl@fau.de>
//
//======================================================================================================================
#include <mesa_pd/data/ParticleAccessor.h>
#include <mesa_pd/data/ParticleStorage.h>
#include <mesa_pd/data/ShapeStorage.h>
namespace walberla {
namespace mesa_pd {
class ParticleAccessorWithShape : public data::ParticleAccessor
{
public:
ParticleAccessorWithShape(std::shared_ptr<data::ParticleStorage>& ps, std::shared_ptr<data::ShapeStorage>& ss)
: ParticleAccessor(ps)
, ss_(ss)
{}
const walberla::real_t& getInvMass(const size_t p_idx) const {return ss_->shapes[ps_->getShapeIDRef(p_idx)]->getInvMass();}
walberla::real_t& getInvMassRef(const size_t p_idx) {return ss_->shapes[ps_->getShapeIDRef(p_idx)]->getInvMass();}
void setInvMass(const size_t p_idx, const walberla::real_t& v) { ss_->shapes[ps_->getShapeIDRef(p_idx)]->getInvMass() = v;}
const auto& getInvInertiaBF(const size_t p_idx) const {return ss_->shapes[ps_->getShapeIDRef(p_idx)]->getInvInertiaBF();}
auto& getInvInertiaBFRef(const size_t p_idx) {return ss_->shapes[ps_->getShapeIDRef(p_idx)]->getInvInertiaBF();}
void setInvInertiaBF(const size_t p_idx, const Mat3& v) { ss_->shapes[ps_->getShapeIDRef(p_idx)]->getInvInertiaBF() = v;}
data::BaseShape* getShape(const size_t p_idx) const {return ss_->shapes[ps_->getShapeIDRef(p_idx)].get();}
private:
std::shared_ptr<data::ShapeStorage> ss_;
};
} // namespace mesa_pd
} // namespace walberla
......@@ -2,13 +2,13 @@ waLBerla_link_files_to_builddir( *.cfg )
waLBerla_link_files_to_builddir( *.py )
waLBerla_add_executable ( NAME PE_GranularGas
FILES PE_GranularGas.cpp
FILES PE_GranularGas.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp
DEPENDS blockforest core pe postprocessing )
waLBerla_add_executable ( NAME MESA_PD_GranularGas
FILES MESA_PD_GranularGas.cpp
FILES MESA_PD_GranularGas.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
DEPENDS blockforest core pe mesa_pd postprocessing vtk )
waLBerla_add_executable ( NAME MESA_PD_KernelBenchmark
FILES MESA_PD_KernelBenchmark.cpp
FILES MESA_PD_KernelBenchmark.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
DEPENDS blockforest core pe mesa_pd postprocessing vtk )
# -*- coding: utf-8 -*-
from jinja2 import Environment, FileSystemLoader
import os
class Parameter:
def __init__(self, name, type, defValue=""):
"""Propery of a data strcuture
Parameters
----------
name : str
name of the property
type : str
type of the property
defValue : str
default value the property should be initialized with
"""
self.name = name
self.type = type
self.defValue = defValue
def __str__(self):
return "name: {}, type: {}, defValue: {}".format(self.name, self.type, self.defValue)
class Config:
def __init__(self):
self.parameters = []
def parameterExists(self, name):
for v in self.parameters:
if v.name==name:
return True
return False
def addParameter(self, name, type, defValue):
if self.parameterExists( name ):
print("parameters already added: " + name)
else:
self.parameters.append( Parameter(name, type, defValue) )
def generateFile(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)
fout = open(path + filename, "wb")
content = env.get_template(template).render(context)
fout.write(content.encode('utf8'))
fout.close()
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.generateFile("Parameters.templ.h")
self.generateFile("Parameters.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 Contact.cpp
//! \author Sebastian Eibl <sebastian.eibl@fau.de>
//
//======================================================================================================================
#include <mesa_pd/data/DataTypes.h>
namespace walberla {
namespace mesa_pd {
struct Contact
{
Contact(const size_t idx1,
const size_t idx2,
const Vec3 contactPoint,
const Vec3 contactNormal,
const real_t penetrationDepth)
: idx1_(idx1)
, idx2_(idx2)
, contactPoint_(contactPoint)
, contactNormal_(contactNormal)
, penetrationDepth_(penetrationDepth) {}
size_t idx1_;
size_t idx2_;
Vec3 contactPoint_;
Vec3 contactNormal_;
real_t penetrationDepth_;
};
} // namespace mesa_pd
} // namespace walberla
//======================================================================================================================
//
// 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 CreateParticles.h
//! \author Sebastian Eibl <sebastian.eibl@fau.de>
//
//======================================================================================================================
#include "CreateParticles.h"
namespace walberla {
namespace mesa_pd {
data::ParticleStorage::iterator createPlane( data::ParticleStorage& ps,
data::ShapeStorage& ss,
const Vec3& pos,
const Vec3& normal )
{
auto p0 = ps.create(true);
p0->getPositionRef() = pos;
p0->getShapeIDRef() = ss.create<data::HalfSpace>( normal );
p0->getOwnerRef() = walberla::mpi::MPIManager::instance()->rank();
p0->getTypeRef() = 0;
data::particle_flags::set(p0->getFlagsRef(), data::particle_flags::INFINITE);
data::particle_flags::set(p0->getFlagsRef(), data::particle_flags::FIXED);
data::particle_flags::set(p0->getFlagsRef(), data::particle_flags::NON_COMMUNICATING);
return p0;
}
data::ParticleStorage::iterator createSphere( data::ParticleStorage& ps,
const Vec3& pos,
const real_t& radius,
const uint64_t shapeID)
{
auto p = ps.create();
p->getPositionRef() = pos;
p->getInteractionRadiusRef() = radius;
p->getShapeIDRef() = shapeID;
p->getOwnerRef() = walberla::MPIManager::instance()->rank();
p->getTypeRef() = 0;
return p;
}
} // namespace mesa_pd
} // namespace walberla
//======================================================================================================================
//
// 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 CreateParticles.cpp
//! \author Sebastian Eibl <sebastian.eibl@fau.de>
//
//======================================================================================================================
#include <mesa_pd/data/DataTypes.h>
#include <mesa_pd/data/ParticleStorage.h>
#include <mesa_pd/data/ShapeStorage.h>
#include <core/mpi/MPIManager.h>
namespace walberla {
namespace mesa_pd {
data::ParticleStorage::iterator createPlane( data::ParticleStorage& ps,
data::ShapeStorage& ss,
const Vec3& pos,
const Vec3& normal );
data::ParticleStorage::iterator createSphere( data::ParticleStorage& ps,
const Vec3& pos,
const real_t& radius,
const uint64_t shapeID);
} // namespace mesa_pd
} // namespace walberla
GranularGas
{
simulationCorner < 0, 0, 0 >;
simulationDomain < 6, 6, 6 >;
simulationDomain < 40, 40, 40 >;
blocks < 2,2,2 >;
isPeriodic < 1, 1, 1 >;
initialRefinementLevel 0;
sorting none;
radius 0.6;
spacing 1.0;
......
......@@ -18,6 +18,16 @@
//
//======================================================================================================================
#include "Accessor.h"
#include "check.h"
#include "Contact.h"
#include "CreateParticles.h"
#include "NodeTimings.h"
#include "Parameters.h"
#include "SelectProperty.h"
#include "sortParticleStorage.h"
#include "SQLProperties.h"
#include <mesa_pd/vtk/ParticleVtkOutput.h>
#include <mesa_pd/collision_detection/AnalyticContactDetection.h>
......@@ -59,54 +69,6 @@
namespace walberla {
namespace mesa_pd {
class SelectRank
{
public:
using return_type = int;
int operator()(const data::Particle& /*p*/) const { return rank_; }
int operator()(const data::Particle&& /*p*/) const { return rank_; }
private:
int rank_ = walberla::mpi::MPIManager::instance()->rank();
};
class ParticleAccessorWithShape : public data::ParticleAccessor
{
public:
ParticleAccessorWithShape(std::shared_ptr<data::ParticleStorage>& ps, std::shared_ptr<data::ShapeStorage>& ss)
: ParticleAccessor(ps)
, ss_(ss)
{}
const auto& getInvMass(const size_t p_idx) const {return ss_->shapes[ps_->getShapeID(p_idx)]->getInvMass();}
const auto& getInvInertiaBF(const size_t p_idx) const {return ss_->shapes[ps_->getShapeID(p_idx)]->getInvInertiaBF();}
data::BaseShape* getShape(const size_t p_idx) const {return ss_->shapes[ps_->getShapeID(p_idx)].get();}
private:
std::shared_ptr<data::ShapeStorage> ss_;
};
void createPlane( data::ParticleStorage& ps,
data::ShapeStorage& ss,
const Vec3& pos,
const Vec3& normal )
{
auto p0 = ps.create(true);
p0->getPositionRef() = pos;
p0->getShapeIDRef() = ss.create<data::HalfSpace>( normal );
p0->getOwnerRef() = walberla::mpi::MPIManager::instance()->rank();
p0->getTypeRef() = 0;
data::particle_flags::set(p0->getFlagsRef(), data::particle_flags::INFINITE);
data::particle_flags::set(p0->getFlagsRef(), data::particle_flags::FIXED);
data::particle_flags::set(p0->getFlagsRef(), data::particle_flags::NON_COMMUNICATING);
}
std::string envToString(const char* env)
{
return env != nullptr ? std::string(env) : "";
}
int main( int argc, char ** argv )
{
using namespace walberla::timing;
......@@ -128,40 +90,8 @@ int main( int argc, char ** argv )
if (cfg == nullptr) WALBERLA_ABORT("No config specified!");
const Config::BlockHandle mainConf = cfg->getBlock( "GranularGas" );
const std::string host = mainConf.getParameter<std::string>("host", "none" );
WALBERLA_LOG_INFO_ON_ROOT("host: " << host);
const int jobid = mainConf.getParameter<int>("jobid", 0 );
WALBERLA_LOG_INFO_ON_ROOT("jobid: " << jobid);
const real_t spacing = mainConf.getParameter<real_t>("spacing", real_t(1.0) );
WALBERLA_LOG_INFO_ON_ROOT("spacing: " << spacing);
const real_t radius = mainConf.getParameter<real_t>("radius", real_t(0.5) );
WALBERLA_LOG_INFO_ON_ROOT("radius: " << radius);
bool bBarrier = mainConf.getParameter<bool>("bBarrier", false );
WALBERLA_LOG_INFO_ON_ROOT("bBarrier: " << bBarrier);
int64_t numOuterIterations = mainConf.getParameter<int64_t>("numOuterIterations", 10 );
WALBERLA_LOG_INFO_ON_ROOT("numOuterIterations: " << numOuterIterations);
int64_t initialRefinementLevel = mainConf.getParameter<int64_t>("initialRefinementLevel", 0 );
WALBERLA_LOG_INFO_ON_ROOT("initialRefinementLevel: " << initialRefinementLevel);
int64_t simulationSteps = mainConf.getParameter<int64_t>("simulationSteps", 10 );
WALBERLA_LOG_INFO_ON_ROOT("simulationSteps: " << simulationSteps);
real_t dt = mainConf.getParameter<real_t>("dt", real_c(0.01) );
WALBERLA_LOG_INFO_ON_ROOT("dt: " << dt);
const int visSpacing = mainConf.getParameter<int>("visSpacing", 1000 );
WALBERLA_LOG_INFO_ON_ROOT("visSpacing: " << visSpacing);
const std::string path = mainConf.getParameter<std::string>("path", "vtk_out" );
WALBERLA_LOG_INFO_ON_ROOT("path: " << path);
const std::string sqlFile = mainConf.getParameter<std::string>("sqlFile", "benchmark.sqlite" );
WALBERLA_LOG_INFO_ON_ROOT("sqlFile: " << sqlFile);
Parameters params;
loadFromConfig(params, mainConf);
WALBERLA_LOG_INFO_ON_ROOT("*** BLOCKFOREST ***");
// create forest
......@@ -186,29 +116,25 @@ int main( int argc, char ** argv )
auto ps = std::make_shared<data::ParticleStorage>(100);
auto ss = std::make_shared<data::ShapeStorage>();
ParticleAccessorWithShape accessor(ps, ss);
data::LinkedCells lc(localDomain.getExtended(spacing), spacing );
data::LinkedCells lc(localDomain.getExtended(params.spacing), params.spacing );
auto smallSphere = ss->create<data::Sphere>( radius );
auto smallSphere = ss->create<data::Sphere>( params.radius );
ss->shapes[smallSphere]->updateMassAndInertia(real_t(2707));
for (auto& iBlk : *forest)
{
for (auto pt : grid_generator::SCGrid(iBlk.getAABB(), Vector3<real_t>(spacing, spacing, spacing) * real_c(0.5), spacing))
for (auto pt : grid_generator::SCGrid(iBlk.getAABB(),
Vector3<real_t>(params.spacing) * real_c(0.5),
params.spacing))
{
WALBERLA_CHECK(iBlk.getAABB().contains(pt));
auto p = ps->create();
p->getPositionRef() = pt;
p->getInteractionRadiusRef() = radius;
p->getShapeIDRef() = smallSphere;
p->getOwnerRef() = mpiManager->rank();
p->getTypeRef() = 0;
createSphere(*ps, pt, params.radius, smallSphere);
}
}
int64_t numParticles = int64_c(ps->size());
walberla::mpi::reduceInplace(numParticles, walberla::mpi::SUM);
WALBERLA_LOG_INFO_ON_ROOT("#particles created: " << numParticles);
auto shift = (spacing - radius - radius) * real_t(0.5);
auto shift = (params.spacing - params.radius - params.radius) * real_t(0.5);
auto confiningDomain = simulationDomain.getExtended(shift);
if (!forest->isPeriodic(0))
......@@ -241,7 +167,7 @@ int main( int argc, char ** argv )
WALBERLA_LOG_INFO_ON_ROOT("*** SIMULATION - START ***");
// Init kernels
kernel::ExplicitEulerWithShape explicitEulerWithShape( dt );
kernel::ExplicitEulerWithShape explicitEulerWithShape( params.dt );
kernel::InsertParticleIntoLinkedCells ipilc;
kernel::SpringDashpot dem(1);
dem.setStiffness(0, 0, real_t(0));
......@@ -256,8 +182,10 @@ int main( int argc, char ** argv )
// initial sync
SNN(*ps, domain);
sortParticleStorage(*ps, params.sorting, lc.domain_, uint_c(lc.numCellsPerDim_[0]));
// vtkWriter->write();
for (int64_t outerIteration = 0; outerIteration < numOuterIterations; ++outerIteration)
for (int64_t outerIteration = 0; outerIteration < params.numOuterIterations; ++outerIteration)
{
WALBERLA_LOG_INFO_ON_ROOT("*** RUNNING OUTER ITERATION " << outerIteration << " ***");
......@@ -274,9 +202,9 @@ int main( int argc, char ** argv )
int64_t contactsChecked = 0;
int64_t contactsDetected = 0;
int64_t contactsTreated = 0;
if (bBarrier) WALBERLA_MPI_BARRIER();
if (params.bBarrier) WALBERLA_MPI_BARRIER();
timer.start();
for (int64_t i=0; i < simulationSteps; ++i)
for (int64_t i=0; i < params.simulationSteps; ++i)
{
// if (i % visSpacing == 0)
// {
......@@ -286,7 +214,7 @@ int main( int argc, char ** argv )
tp["GenerateLinkedCells"].start();
lc.clear();
ps->forEachParticle(true, kernel::SelectAll(), accessor, ipilc, accessor, lc);
if (bBarrier) WALBERLA_MPI_BARRIER();
if (params.bBarrier) WALBERLA_MPI_BARRIER();
tp["GenerateLinkedCells"].end();
tp["DEM"].start();
......@@ -310,23 +238,23 @@ int main( int argc, char ** argv )
}
},
accessor );
if (bBarrier) WALBERLA_MPI_BARRIER();
if (params.bBarrier) WALBERLA_MPI_BARRIER();
tp["DEM"].end();
tp["ReduceForce"].start();
RP.operator()<ForceTorqueNotification>(*ps);
if (bBarrier) WALBERLA_MPI_BARRIER();
if (params.bBarrier) WALBERLA_MPI_BARRIER();
tp["ReduceForce"].end();
tp["Euler"].start();
//ps->forEachParticle(false, [&](const size_t idx){WALBERLA_CHECK_EQUAL(ps->getForce(idx), Vec3(0,0,0), *(*ps)[idx] << "\n" << idx);});
ps->forEachParticle(true, kernel::SelectLocal(), accessor, explicitEulerWithShape, accessor);
if (bBarrier) WALBERLA_MPI_BARRIER();
if (params.bBarrier) WALBERLA_MPI_BARRIER();
tp["Euler"].end();
tp["SNN"].start();
SNN(*ps, domain);
if (bBarrier) WALBERLA_MPI_BARRIER();
if (params.bBarrier) WALBERLA_MPI_BARRIER();
tp["SNN"].end();
}
timer.end();
......@@ -362,7 +290,7 @@ int main( int argc, char ** argv )
{
WALBERLA_LOG_INFO_ON_ROOT(*timer_reduced);
WALBERLA_LOG_INFO_ON_ROOT("runtime: " << timer_reduced->max());
PUpS = double_c(numParticles) * double_c(simulationSteps) / double_c(timer_reduced->max());
PUpS = double_c(numParticles) * double_c(params.simulationSteps) / double_c(timer_reduced->max());
WALBERLA_LOG_INFO_ON_ROOT("PUpS: " << PUpS);
}
......@@ -370,19 +298,10 @@ int main( int argc, char ** argv )
WALBERLA_LOG_INFO_ON_ROOT(*tp_reduced);
WALBERLA_LOG_INFO_ON_ROOT("*** SIMULATION - END ***");
WALBERLA_LOG_INFO_ON_ROOT("*** CHECKING RESULT - START ***");
auto pIt = ps->begin();
for (auto& iBlk : *forest)
if (params.checkSimulation)
{
for (auto it = grid_generator::SCIterator(iBlk.getAABB(), Vector3<real_t>(spacing, spacing, spacing) * real_c(0.5), spacing);
it != grid_generator::SCIterator();
++it, ++pIt)
{
WALBERLA_CHECK_UNEQUAL(pIt, ps->end());
WALBERLA_CHECK_FLOAT_EQUAL((*pIt).getPositionRef(), *it);
}
check(*ps, *forest, params.spacing);
}
WALBERLA_LOG_INFO_ON_ROOT("*** CHECKING RESULT - END ***");
WALBERLA_LOG_INFO_ON_ROOT("*** SQL OUTPUT - START ***");
numParticles = 0;
......@@ -416,6 +335,8 @@ int main( int argc, char ** argv )
walberla::mpi::reduceInplace(local_aabbs, walberla::mpi::SUM);
walberla::mpi::reduceInplace(neighbor_subdomains, walberla::mpi::SUM);
walberla::mpi::reduceInplace(neighbor_processes, walberla::mpi::SUM);
uint_t runId = uint_c(-1);
WALBERLA_ROOT_SECTION()
{
std::map< std::string, walberla::int64_t > integerProperties;
......@@ -424,26 +345,14 @@ int main( int argc, char ** argv )
stringProperties["walberla_git"] = WALBERLA_GIT_SHA1;
stringProperties["tag"] = "mesa_pd";
stringProperties["host"] = host;
integerProperties["jobid"] = jobid;
integerProperties["mpi_num_processes"] = mpiManager->numProcesses();
integerProperties["omp_max_threads"] = omp_get_max_threads();
integerProperties["numOuterIterations"] = numOuterIterations;
integerProperties["simulationSteps"] = simulationSteps;
integerProperties["bBarrier"] = int64_c(bBarrier);
realProperties["PUpS"] = double_c(PUpS);
integerProperties["outerIteration"] = int64_c(outerIteration);
integerProperties["num_particles"] = numParticles;
integerProperties["num_ghost_particles"] = numGhostParticles;
integerProperties["contacts_checked"] = contactsChecked;
integerProperties["contacts_detected"] = contactsDetected;
integerProperties["contacts_treated"] = contactsTreated;
integerProperties["blocks_x"] = int64_c(forest->getXSize());
integerProperties["blocks_y"] = int64_c(forest->getYSize());
integerProperties["blocks_z"] = int64_c(forest->getZSize());
integerProperties["initialRefinementLevel"] = int64_c(initialRefinementLevel);
realProperties["domain_x"] = double_c(forest->getDomain().xSize());
realProperties["domain_y"] = double_c(forest->getDomain().ySize());
realProperties["domain_z"] = double_c(forest->getDomain().zSize());
integerProperties["local_aabbs"] = int64_c(local_aabbs);
integerProperties["neighbor_subdomains"] = int64_c(neighbor_subdomains);
integerProperties["neighbor_processes"] = int64_c(neighbor_processes);
......@@ -457,26 +366,18 @@ int main( int argc, char ** argv )
integerProperties["RPReceives"] = RPReceives;
realProperties["linkedCellsVolume"] = linkedCellsVolume;
integerProperties["numLinkedCells"] = int64_c(numLinkedCells);
realProperties["timer_min"] = timer_reduced->min();