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/core/blockLatticeStructure2D.hh | 330 ++++++++++++++++++++++++++++++++++++
1 file changed, 330 insertions(+)
create mode 100644 src/core/blockLatticeStructure2D.hh
(limited to 'src/core/blockLatticeStructure2D.hh')
diff --git a/src/core/blockLatticeStructure2D.hh b/src/core/blockLatticeStructure2D.hh
new file mode 100644
index 0000000..e2577a4
--- /dev/null
+++ b/src/core/blockLatticeStructure2D.hh
@@ -0,0 +1,330 @@
+/* This file is part of the OpenLB library
+ *
+ * Copyright (C) 2006-2018 Jonas Latt, Adrian Kummerlaender
+ * 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
+ * Dynamics for a generic 2D block structure -- header file.
+ */
+#ifndef BLOCK_LATTICE_STRUCTURE_2D_HH
+#define BLOCK_LATTICE_STRUCTURE_2D_HH
+
+#include "blockLatticeStructure2D.h"
+#include "functors/lattice/indicator/blockIndicatorBaseF2D.hh"
+#include "functors/lattice/indicator/blockIndicatorF2D.hh"
+
+
+namespace olb {
+
+
+template
+void BlockLatticeStructure2D::defineRho(
+ BlockIndicatorF2D& indicator, AnalyticalF2D& rho)
+{
+ T physR[2] = { };
+ T rhoTmp = T();
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ if (indicator(iX, iY)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY);
+ rho(&rhoTmp, physR);
+ get(iX, iY).defineRho(rhoTmp);
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure2D::defineRho(
+ BlockGeometryStructure2D& blockGeometry, int material, AnalyticalF2D& rho)
+{
+ BlockIndicatorMaterial2D indicator(blockGeometry, material);
+ defineRho(indicator, rho);
+}
+
+template
+void BlockLatticeStructure2D::defineU(
+ BlockIndicatorF2D& indicator, AnalyticalF2D& u)
+{
+ T physR[2] = { };
+ T uTmp[2] = { };
+ const auto& geometry = indicator.getBlockGeometryStructure();
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ if (indicator(iX, iY)) {
+ geometry.getPhysR(physR, iX, iY);
+ u(uTmp, physR);
+ get(iX, iY).defineU(uTmp);
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure2D::defineU(
+ BlockGeometryStructure2D& blockGeometry, int material, AnalyticalF2D& u)
+{
+ BlockIndicatorMaterial2D indicator(blockGeometry, material);
+ defineU(indicator, u);
+}
+
+template
+void BlockLatticeStructure2D::defineRhoU(
+ BlockIndicatorF2D& indicator,
+ AnalyticalF2D& rho, AnalyticalF2D& u)
+{
+ T physR[2] = { };
+ T uTmp[2] = { };
+ T rhoTmp = T();
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ if (indicator(iX, iY)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY);
+ rho(&rhoTmp, physR);
+ u(uTmp, physR);
+ get(iX, iY).defineRhoU(rhoTmp, uTmp);
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure2D::defineRhoU(
+ BlockGeometryStructure2D& blockGeometry, int material,
+ AnalyticalF2D& rho, AnalyticalF2D& u)
+{
+ BlockIndicatorMaterial2D indicator(blockGeometry, material);
+ defineRhoU(indicator, rho, u);
+}
+
+template
+void BlockLatticeStructure2D::definePopulations(
+ BlockIndicatorF2D& indicator, AnalyticalF2D& Pop)
+{
+ T physR[2] = { };
+ T PopTmp[DESCRIPTOR::q];
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ if (indicator(iX, iY)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY);
+ Pop(PopTmp, physR);
+ get(iX, iY).definePopulations(PopTmp);
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure2D::definePopulations(
+ BlockGeometryStructure2D& blockGeometry, int material, AnalyticalF2D& Pop)
+{
+ BlockIndicatorMaterial2D indicator(blockGeometry, material);
+ definePopulations(indicator, Pop);
+}
+
+template
+template
+void BlockLatticeStructure2D::defineField(
+ BlockIndicatorF2D& indicator, AnalyticalF2D& field)
+{
+ T* fieldTmp = new T[DESCRIPTOR::template size()];
+ T physR[2] = { };
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ if (indicator(iX, iY)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY);
+ field(fieldTmp, physR);
+ get(iX, iY).template defineField(fieldTmp);
+ }
+ }
+ }
+ delete[] fieldTmp;
+}
+
+template
+template
+void BlockLatticeStructure2D::defineField(
+ BlockGeometryStructure2D& blockGeometry, int material, AnalyticalF2D& field)
+{
+ BlockIndicatorMaterial2D indicator(blockGeometry, material);
+ defineField(indicator, field);
+}
+
+template
+template
+void BlockLatticeStructure2D::defineField(
+ BlockGeometryStructure2D& blockGeometry, IndicatorF2D& indicatorF,
+ AnalyticalF2D& field)
+{
+ BlockIndicatorFfromIndicatorF2D indicator(indicatorF, blockGeometry);
+ defineField(indicator, field);
+}
+
+template
+template
+void BlockLatticeStructure2D::addField(
+ BlockGeometryStructure2D& blockGeometry, IndicatorF2D& indicator,
+ AnalyticalF2D& field)
+{
+ T* fieldTmp = new T[DESCRIPTOR::template size()];
+ T physR[2] = { };
+ bool inside;
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ blockGeometry.getPhysR(physR, iX, iY);
+ indicator(&inside, &(physR[0]));
+ if (inside) {
+ field(fieldTmp,physR);
+ get(iX, iY).template addField(fieldTmp);
+ }
+ }
+ }
+ delete[] fieldTmp;
+}
+
+template
+template
+void BlockLatticeStructure2D::addField(
+ BlockGeometryStructure2D& blockGeometry, IndicatorF2D& indicator,
+ AnalyticalF2D& field, AnalyticalF2D& porous)
+{
+ T* fieldTmp = new T[DESCRIPTOR::template size()];
+ bool inside;
+ T physR[2] = { };
+ T porousA[1] = { };
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ blockGeometry.getPhysR(physR, iX, iY);
+ indicator(&inside, physR);
+ if (inside) {
+ porous(porousA, physR);
+ field(fieldTmp,physR);
+ for (int i = 0; i < DESCRIPTOR::template size(); ++i) {
+ fieldTmp[i] *= porousA[0];
+ }
+ get(iX, iY).template addField(fieldTmp);
+ }
+ }
+ }
+ delete[] fieldTmp;
+}
+
+#ifndef OLB_PRECOMPILED
+template
+void BlockLatticeStructure2D::setExternalParticleField(BlockGeometryStructure2D& blockGeometry, AnalyticalF2D& velocity, SmoothIndicatorF2D& sIndicator)
+{
+
+ int start[2] = {0};
+ int end[2] = {0};
+ // check for intersection of cuboid and indicator
+ Cuboid2D tmpCuboid(blockGeometry.getOrigin()[0], blockGeometry.getOrigin()[1], blockGeometry.getDeltaR(), blockGeometry.getNx(), blockGeometry.getNy());
+ T posXmin = sIndicator.getPos()[0] - sIndicator.getCircumRadius();
+ T posXmax = sIndicator.getPos()[0] + sIndicator.getCircumRadius();
+ T posYmin = sIndicator.getPos()[1] - sIndicator.getCircumRadius();
+ T posYmax = sIndicator.getPos()[1] + sIndicator.getCircumRadius();
+ if(tmpCuboid.checkInters(posXmin, posXmax, posYmin, posYmax, start[0], end[0], start[1], end[1])) {
+
+ for (int k=0; k<2; k++) {
+ start[k] -= 1;
+ if(start[k] < 0) start[k] = 0;
+ end[k] += 2;
+ if(end[k] > blockGeometry.getExtend()[k]) end[k] = blockGeometry.getExtend()[k];
+ }
+
+ T foo[3] = { }; /// Contains foo[0]=vel0; foo[1]=vel1; foo[2]=porosity
+ T physR[2]= { };
+ T porosity[1] = { };
+ for (int iX = start[0]; iX < end[0]; ++iX) {
+ for (int iY = start[1]; iY < end[1]; ++iY) {
+ blockGeometry.getPhysR(physR, iX, iY);
+ sIndicator(porosity, physR);
+ if (!util::nearZero(porosity[0])) {
+ // TODO: Check / adapt to use descriptor fields
+ velocity(foo,physR);
+ foo[0] *= porosity[0];
+ foo[1] *= porosity[0];
+ foo[2] = porosity[0];
+ get(iX, iY).template addField(foo);
+ get(iX, iY).template addField(&foo[2]);
+ porosity[0] = 1.-porosity[0];
+ *(get(iX, iY).template getFieldPointer()) *= porosity[0];
+ }
+ }
+ }
+ }
+}
+#endif
+
+template
+template
+void BlockLatticeStructure2D::multiplyField(
+ BlockGeometryStructure2D& blockGeometry, IndicatorF2D& indicator,
+ AnalyticalF2D& field)
+{
+ T* fieldTmp = new T [DESCRIPTOR::template size()];
+ bool inside;
+ T physR[3] = { };
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ blockGeometry.getPhysR(physR, iX, iY);
+ indicator(&inside, &(physR[0]));
+ if (inside) {
+ field(fieldTmp, physR);
+ get(iX, iY).template multiplyField(fieldTmp);
+ }
+ }
+ }
+ delete[] fieldTmp;
+}
+
+template
+void BlockLatticeStructure2D::iniEquilibrium(
+ BlockIndicatorF2D& indicator,
+ AnalyticalF2D& rho, AnalyticalF2D& u)
+{
+ T physR[2] = { };
+ T uTmp[2] = { };
+ T rhoTmp = T();
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ if (indicator(iX, iY)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY);
+ u(uTmp, physR);
+ rho(&rhoTmp, physR);
+ get(iX, iY).iniEquilibrium(rhoTmp, uTmp);
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure2D::iniEquilibrium(
+ BlockGeometryStructure2D& blockGeometry, int material,
+ AnalyticalF2D& rho, AnalyticalF2D& u)
+{
+ BlockIndicatorMaterial2D indicator(blockGeometry, material);
+ iniEquilibrium(indicator, rho, u);
+}
+
+
+} // namespace olb
+
+#endif
--
cgit v1.2.3