[ADD] kernel for the VBondModel, viscous damping

src/mesa_pd/kernel/VBondModel/ViscousDamping.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/common/ParticleFunctions.h>
+#include <mesa_pd/data/DataTypes.h>
+#include <mesa_pd/data/IAccessor.h>
+
+#include <core/math/Constants.h>
+#include <core/logging/Logging.h>
+
+#include <vector>
+
+namespace walberla {
+namespace mesa_pd {
+namespace kernel {
+namespace VBondModel {
+
+class ViscousDamping
+{
+public:
+   ViscousDamping(real_t forceDampingFactor, real_t torqueDampingFactor)
+   : forceDampingFactor_(forceDampingFactor)
+   , torqueDampingFactor_(torqueDampingFactor)
+   {}
+
+   template<typename Accessor>
+   void operator()(const size_t p_idx1,
+                   const size_t p_idx2,
+                   Accessor &ac) const;
+
+   auto getForceDampingFactor() const {return forceDampingFactor_;}
+   auto getTorqueDampingFactor() const {return torqueDampingFactor_;}
+private:
+   const real_t forceDampingFactor_;
+   const real_t torqueDampingFactor_;
+};
+
+template<typename Accessor>
+inline void ViscousDamping::operator()(const size_t p_idx1,
+                                       const size_t p_idx2,
+                                       Accessor &ac) const
+{
+   Vec3 velDampingForce = (ac.getLinearVelocity(p_idx1) - ac.getLinearVelocity(p_idx2)) * forceDampingFactor_;
+   addForceAtomic(p_idx1, ac, -velDampingForce);
+   addForceAtomic(p_idx2, ac,  velDampingForce);
+
+   Vec3 angDampingTorque = (ac.getAngularVelocity(p_idx1) - ac.getAngularVelocity(p_idx2)) * torqueDampingFactor_;
+   addTorqueAtomic(p_idx1, ac, -angDampingTorque);
+   addTorqueAtomic(p_idx2, ac,  angDampingTorque);
+
+}
+
+} //namespace VBondModel
+} //namespace kernel
+} //namespace mesa_pd
+} //namespace walberla
\ No newline at end of file

tests/mesa_pd/CMakeLists.txt

@@ -164,6 +164,9 @@ waLBerla_execute_test( NAME   MESA_PD_Kernel_SyncNextNeighborsBlockForest PROCES
 waLBerla_compile_test( NAME   MESA_PD_Kernel_TemperatureIntegration FILES kernel/TemperatureIntegration.cpp DEPENDS core )
 waLBerla_execute_test( NAME   MESA_PD_Kernel_TemperatureIntegration )
 
+waLBerla_compile_test( NAME   MESA_PD_Kernel_VBondModel_ViscousDamping FILES kernel/VBondModel/ViscousDamping.test.cpp DEPENDS core )
+waLBerla_execute_test( NAME   MESA_PD_Kernel_VBondModel_ViscousDamping )
+
 waLBerla_compile_test( NAME   MESA_PD_Kernel_VBondModel_WallContact FILES kernel/VBondModel/WallContact.test.cpp DEPENDS core )
 waLBerla_execute_test( NAME   MESA_PD_Kernel_VBondModel_WallContact )
 

tests/mesa_pd/kernel/VBondModel/ViscousDamping.test.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 <mesa_pd/data/ParticleAccessor.h>
+#include <mesa_pd/data/ParticleStorage.h>
+
+#include <mesa_pd/kernel/VBondModel/ViscousDamping.h>
+
+#include <core/Environment.h>
+#include <core/logging/Logging.h>
+
+#include <iostream>
+
+namespace walberla {
+namespace mesa_pd {
+
+int main( int argc, char ** argv )
+{
+
+   Environment env(argc, argv);
+   WALBERLA_UNUSED(env);
+   mpi::MPIManager::instance()->useWorldComm();
+
+   if (std::is_same<real_t, float>::value)
+   {
+      WALBERLA_LOG_WARNING("waLBerla build in sp mode: skipping test due to low precision");
+      return EXIT_SUCCESS;
+   }
+
+   //init data structures
+   auto ps = std::make_shared<data::ParticleStorage>(2);
+   data::Particle &&p1 = *ps->create();
+   data::Particle &&p2 = *ps->create();
+
+   data::ParticleAccessor ac(ps);
+
+   //init kernels
+   kernel::VBondModel::ViscousDamping damping(0.2_r,0.5_r);
+   WALBERLA_CHECK_FLOAT_EQUAL(damping.getForceDampingFactor(),
+                              0.2_r);
+   WALBERLA_CHECK_FLOAT_EQUAL(damping.getTorqueDampingFactor(),
+                              0.5_r);
+
+   p1.setLinearVelocity(Vec3(1_r, 2_r, 3_r));
+   p2.setLinearVelocity(Vec3(3_r, 2_r, 1_r));
+   p1.setAngularVelocity(Vec3(2_r, 3_r, 4_r));
+   p2.setAngularVelocity(Vec3(3_r, 2_r, 1_r));
+   damping(0, 1, ac);
+   WALBERLA_CHECK_FLOAT_EQUAL(p1.getForce(),
+                              Vec3(0.4_r, 0_r, -0.4_r));
+   WALBERLA_CHECK_FLOAT_EQUAL(p2.getForce(),
+                              Vec3(-0.4_r, 0_r, 0.4_r));
+   WALBERLA_CHECK_FLOAT_EQUAL(p1.getTorque(),
+                              Vec3(0.5_r, -0.5_r, -1.5_r));
+   WALBERLA_CHECK_FLOAT_EQUAL(p2.getTorque(),
+                              Vec3(-0.5_r, 0.5_r, 1.5_r));
+
+   return EXIT_SUCCESS;
+}
+
+} //namespace mesa_pd
+} //namespace walberla
+
+int main( int argc, char ** argv )
+{
+   return walberla::mesa_pd::main(argc, argv);
+}