From 94d3e79a8617f88dc0219cfdeedfa3147833719d Mon Sep 17 00:00:00 2001
From: Adrian Kummerlaender
Date: Mon, 24 Jun 2019 14:43:36 +0200
Subject: Initialize at openlb-1-3
---
src/dynamics/mrtDynamics.hh | 203 ++++++++++++++++++++++++++++++++++++++++++++
1 file changed, 203 insertions(+)
create mode 100644 src/dynamics/mrtDynamics.hh
(limited to 'src/dynamics/mrtDynamics.hh')
diff --git a/src/dynamics/mrtDynamics.hh b/src/dynamics/mrtDynamics.hh
new file mode 100644
index 0000000..b62ed0f
--- /dev/null
+++ b/src/dynamics/mrtDynamics.hh
@@ -0,0 +1,203 @@
+/* This file is part of the OpenLB library
+ *
+ * Copyright (C) 2006, 2007 Jonas Latt
+ * E-mail contact: info@openlb.net
+ * The most recent release of OpenLB can be downloaded at
+ *
+ *
+ * This program 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 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program 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 this program; if not, write to the Free
+ * Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
+ * Boston, MA 02110-1301, USA.
+*/
+
+/** \file
+ * A collection of dynamics classes (e.g. BGK) with which a Cell object
+ * can be instantiated -- generic implementation.
+ */
+#ifndef MRT_DYNAMICS_HH
+#define MRT_DYNAMICS_HH
+
+#include
+#include
+#include "mrtHelpers.h"
+
+namespace olb {
+
+//==============================================================================//
+/////////////////////////// Class MRTdynamics ///////////////////////////////
+//==============================================================================//
+/** \param omega_ relaxation parameter, related to the dynamic viscosity
+ * \param momenta_ a Momenta object to know how to compute velocity momenta
+ * \param lambda_ will be used as an
+ */
+
+// Original implementation based on:
+// D'Humieres et al., "Multiple-relaxation-time lattice Boltzmann models in three dimensions",
+// Phil: Trans. R. soc. Lond. A (2002) 360, 437-451
+// and
+// Yu et al,, "LES of turbulent square jet flow using an MRT lattice Boltzmann model",
+// Computers & Fluids 35 (2006), 957-965
+template
+MRTdynamics::MRTdynamics (
+ T omega_, Momenta& momenta_ )
+ : BasicDynamics(momenta_), omega(omega_), lambda(omega_)
+{
+ T rt[DESCRIPTOR::q]; // relaxation times vector.
+ for (int iPop = 0; iPop < DESCRIPTOR::q; ++iPop) {
+ rt[iPop] = descriptors::s(iPop);
+ }
+ for (int iPop = 0; iPop < descriptors::shearIndexes(); ++iPop) {
+ rt[descriptors::shearViscIndexes(iPop)] = omega;
+ }
+ for (int iPop = 0; iPop < DESCRIPTOR::q; ++iPop) {
+ for (int jPop = 0; jPop < DESCRIPTOR::q; ++jPop) {
+ invM_S[iPop][jPop] = T();
+ for (int kPop = 0; kPop < DESCRIPTOR::q; ++kPop) {
+ if (kPop == jPop) {
+ invM_S[iPop][jPop] += descriptors::invM(iPop,kPop) *
+ rt[kPop];
+ }
+ }
+ }
+ }
+
+}
+
+template
+T MRTdynamics::computeEquilibrium(int iPop, T rho, const T u[DESCRIPTOR::d], T uSqr) const
+{
+ return lbHelpers::equilibrium(iPop, rho, u, uSqr);
+}
+
+template
+void MRTdynamics::computeAllEquilibrium(T momentaEq[DESCRIPTOR::q],
+ T rho, const T u[DESCRIPTOR::d],
+ const T uSqr)
+{
+ mrtHelpers::computeEquilibrium(momentaEq, rho, u, uSqr);
+}
+
+template
+void MRTdynamics::collide (
+ Cell& cell,
+ LatticeStatistics& statistics )
+{
+ typedef DESCRIPTOR L;
+ typedef mrtHelpers mrtH;
+
+ T rho, u[L::d];
+ this->_momenta.computeRhoU(cell, rho, u);
+
+ T uSqr = mrtH::mrtCollision(cell,rho,u,invM_S);
+
+ statistics.incrementStats(rho, uSqr);
+}
+
+template
+T MRTdynamics::getOmega() const
+{
+ return omega;
+}
+
+template
+void MRTdynamics::setOmega(T omega_)
+{
+ omega = omega_;
+}
+
+template
+T MRTdynamics::getLambda() const
+{
+ return lambda;
+}
+
+template
+void MRTdynamics::setLambda(T lambda_)
+{
+ lambda = lambda_;
+}
+
+
+
+
+template
+MRTdynamics2::MRTdynamics2 (
+ T omega_, Momenta& momenta_ )
+ : MRTdynamics(omega_, momenta_)
+{
+ T rt[DESCRIPTOR::q]; // relaxation times vector.
+ for (int iPop = 0; iPop < DESCRIPTOR::q; ++iPop) {
+ rt[iPop] = descriptors::s_2(iPop);
+ }
+ for (int iPop = 0; iPop < descriptors::shearIndexes(); ++iPop) {
+ rt[descriptors::shearViscIndexes(iPop)] = omega;
+ }
+ for (int iPop = 0; iPop < DESCRIPTOR::q; ++iPop) {
+ for (int jPop = 0; jPop < DESCRIPTOR::q; ++jPop) {
+ invM_S_2[iPop][jPop] = T();
+ for (int kPop = 0; kPop < DESCRIPTOR::q; ++kPop) {
+ if (kPop == jPop) {
+ invM_S_2[iPop][jPop] += descriptors::invM(iPop,kPop) *
+ rt[kPop];
+ }
+ }
+ }
+ }
+}
+
+
+// Stabalized MRT scheme with uniform relaxation times
+template
+void MRTdynamics2::collide (
+ Cell& cell,
+ LatticeStatistics& statistics )
+{
+ typedef mrtHelpers mrtH;
+
+ T rho, u[DESCRIPTOR::d];
+ this->_momenta.computeRhoU(cell, rho, u);
+
+ T uSqr = mrtH::mrtCollision(cell,rho,u,invM_S_2);
+
+ statistics.incrementStats(rho, uSqr);
+}
+
+
+
+template
+ForcedMRTdynamics::ForcedMRTdynamics (
+ T omega_, Momenta& momenta_ )
+ : MRTdynamics(omega_, momenta_)
+{
+}
+
+template
+void ForcedMRTdynamics::collide (
+ Cell& cell,
+ LatticeStatistics& statistics )
+{
+ T rho, u[DESCRIPTOR::d];
+ this->_momenta.computeRhoU(cell, rho, u);
+
+ T uSqr = mrtHelpers::mrtCollision(cell, rho, u, this->invM_S);
+ mrtHelpers::addExternalForce(cell, rho, u, this->invM_S);
+
+ statistics.incrementStats(rho, uSqr);
+}
+
+
+} // end namespace
+
+#endif
+
--
cgit v1.2.3