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/blockLatticeStructure3D.hh | 305 ++++++++++++++++++++++++++++++++++++
1 file changed, 305 insertions(+)
create mode 100644 src/core/blockLatticeStructure3D.hh
(limited to 'src/core/blockLatticeStructure3D.hh')
diff --git a/src/core/blockLatticeStructure3D.hh b/src/core/blockLatticeStructure3D.hh
new file mode 100644
index 0000000..a0af9ac
--- /dev/null
+++ b/src/core/blockLatticeStructure3D.hh
@@ -0,0 +1,305 @@
+/* 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 3D block lattice view -- generic implementation.
+ */
+#ifndef BLOCK_LATTICE_STRUCTURE_3D_HH
+#define BLOCK_LATTICE_STRUCTURE_3D_HH
+
+#include "blockLatticeStructure3D.h"
+#include "functors/lattice/indicator/blockIndicatorBaseF3D.h"
+#include "functors/lattice/indicator/blockIndicatorF3D.h"
+
+namespace olb {
+
+template
+void BlockLatticeStructure3D::defineRho(
+ BlockIndicatorF3D& indicator, AnalyticalF3D& rho)
+{
+ T physR[3] = { };
+ T rhoTmp = T();
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ for (int iZ = 0; iZ < getNz(); ++iZ) {
+ if (indicator(iX, iY, iZ)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY, iZ);
+ rho(&rhoTmp,physR);
+ get(iX, iY, iZ).defineRho(rhoTmp);
+ }
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure3D::defineRho(
+ BlockGeometryStructure3D& blockGeometry, int material, AnalyticalF3D& rho)
+{
+ BlockIndicatorMaterial3D indicator(blockGeometry, material);
+ defineRho(indicator, rho);
+}
+
+template
+void BlockLatticeStructure3D::defineU(
+ BlockIndicatorF3D& indicator, AnalyticalF3D& u)
+{
+ T physR[3] = { };
+ T uTmp[3] = { };
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ for (int iZ = 0; iZ < getNz(); ++iZ) {
+ if (indicator(iX, iY, iZ)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY, iZ);
+ u(uTmp,physR);
+ get(iX, iY, iZ).defineU(uTmp);
+ }
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure3D::defineU(
+ BlockGeometryStructure3D& blockGeometry, int material, AnalyticalF3D& u)
+{
+ BlockIndicatorMaterial3D indicator(blockGeometry, material);
+ defineU(indicator, u);
+}
+
+template
+void BlockLatticeStructure3D::defineRhoU(
+ BlockIndicatorF3D& indicator,
+ AnalyticalF3D& rho, AnalyticalF3D& u)
+{
+ T physR[3] = { };
+ T uTmp[3] = { };
+ T rhoTmp = T();
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ for (int iZ = 0; iZ < getNz(); ++iZ) {
+ if (indicator(iX, iY, iZ)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY, iZ);
+ rho(&rhoTmp,physR);
+ u(uTmp,physR);
+ get(iX, iY, iZ).defineRhoU(rhoTmp,uTmp);
+ }
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure3D::defineRhoU(
+ BlockGeometryStructure3D& blockGeometry, int material,
+ AnalyticalF3D& rho, AnalyticalF3D& u)
+{
+ BlockIndicatorMaterial3D indicator(blockGeometry, material);
+ defineRhoU(indicator, rho, u);
+}
+
+template
+void BlockLatticeStructure3D::definePopulations(
+ BlockIndicatorF3D& indicator, AnalyticalF3D& Pop)
+{
+ T physR[3] = { };
+ T PopTmp[DESCRIPTOR::q];
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ for (int iZ = 0; iZ < getNz(); ++iZ) {
+ if (indicator(iX, iY, iZ)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY, iZ);
+ Pop(PopTmp,physR);
+ get(iX, iY, iZ).definePopulations(PopTmp);
+ }
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure3D::definePopulations(
+ BlockGeometryStructure3D& blockGeometry, int material, AnalyticalF3D& Pop)
+{
+ BlockIndicatorMaterial3D indicator(blockGeometry, material);
+ definePopulations(indicator, Pop);
+}
+
+template
+void BlockLatticeStructure3D::definePopulations(
+ BlockIndicatorF3D& indicator, BlockF3D& Pop)
+{
+ int latticeR[3];
+ T PopTmp[DESCRIPTOR::q];
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ for (int iZ = 0; iZ < getNz(); ++iZ) {
+ if (indicator(iX, iY, iZ)) {
+ latticeR[0] = iX;
+ latticeR[1] = iY;
+ latticeR[2] = iZ;
+ Pop(PopTmp,latticeR);
+ get(iX, iY, iZ).definePopulations(PopTmp);
+ }
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure3D::definePopulations(
+ BlockGeometryStructure3D& blockGeometry, int material, BlockF3D& Pop)
+{
+ BlockIndicatorMaterial3D indicator(blockGeometry, material);
+ definePopulations(indicator, Pop);
+}
+
+template
+template
+void BlockLatticeStructure3D::defineField(
+ BlockIndicatorF3D& indicator, AnalyticalF3D& field)
+{
+ T* fieldTmp = new T[DESCRIPTOR::template size()];
+ T physR[3] = { };
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ for (int iZ = 0; iZ < getNz(); ++iZ) {
+ if (indicator(iX, iY, iZ)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY, iZ);
+ field(fieldTmp, physR);
+ get(iX, iY, iZ).template defineField(fieldTmp);
+ }
+ }
+ }
+ }
+ delete[] fieldTmp;
+}
+
+template
+template
+void BlockLatticeStructure3D::defineField(
+ BlockGeometryStructure3D& blockGeometry, int material,
+ AnalyticalF3D& field)
+{
+ BlockIndicatorMaterial3D indicator(blockGeometry, material);
+ defineField(indicator, field);
+}
+
+template
+template
+void BlockLatticeStructure3D::defineField(
+ BlockGeometryStructure3D& blockGeometry, IndicatorF3D& indicatorF,
+ AnalyticalF3D& field)
+{
+ BlockIndicatorFfromIndicatorF3D indicator(indicatorF, blockGeometry);
+ defineField(indicator, field);
+}
+
+template
+void BlockLatticeStructure3D::iniEquilibrium(
+ BlockIndicatorF3D& indicator,
+ AnalyticalF3D& rho, AnalyticalF3D& u)
+{
+ T physR[3] = { };
+ T uTmp[3] = { };
+ T rhoTmp = T();
+ for (int iX = 0; iX < getNx(); ++iX) {
+ for (int iY = 0; iY < getNy(); ++iY) {
+ for (int iZ = 0; iZ < getNz(); ++iZ) {
+ if (indicator(iX, iY, iZ)) {
+ indicator.getBlockGeometryStructure().getPhysR(physR, iX, iY, iZ);
+ u(uTmp, physR);
+ rho(&rhoTmp, physR);
+ get(iX, iY, iZ).iniEquilibrium(rhoTmp, uTmp);
+ }
+ }
+ }
+ }
+}
+
+template
+void BlockLatticeStructure3D::iniEquilibrium(
+ BlockGeometryStructure3D& blockGeometry, int material,
+ AnalyticalF3D& rho, AnalyticalF3D& u)
+{
+ BlockIndicatorMaterial3D indicator(blockGeometry, material);
+ iniEquilibrium(indicator, rho, u);
+}
+
+#ifndef OLB_PRECOMPILED
+template
+void BlockLatticeStructure3D::setExternalParticleField(
+ BlockGeometryStructure3D& blockGeometry,
+ AnalyticalF3D& velocity, SmoothIndicatorF3D& sIndicator)
+{
+
+ int start[3] = {0};
+ int end[3] = {0};
+ // check for intersection of cuboid and indicator
+ Cuboid3D tmpCuboid(blockGeometry.getOrigin()[0], blockGeometry.getOrigin()[1], blockGeometry.getOrigin()[2], blockGeometry.getDeltaR(), blockGeometry.getNx(), blockGeometry.getNy(), blockGeometry.getNz());
+ 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();
+ T posZmin = sIndicator.getPos()[2] - sIndicator.getCircumRadius();
+ T posZmax = sIndicator.getPos()[2] + sIndicator.getCircumRadius();
+ if(tmpCuboid.checkInters(posXmin, posXmax, posYmin, posYmax, posZmin, posZmax,
+ start[0], end[0], start[1], end[1], start[2], end[2]))
+ {
+ for (int k=0; k<3; 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[4] = { }; /// Contains foo[0]=vel0; foo[1]=vel1; foo[2]=vel2; foo[3]=porosity
+ T physR[3] = { };
+ T porosity[1] = { };
+
+ for (int iX = start[0]; iX < end[0]; ++iX) {
+ for (int iY = start[1]; iY < end[1]; ++iY) {
+ for (int iZ = start[2]; iZ < end[2]; ++iZ) {
+ blockGeometry.getPhysR(physR, iX, iY, iZ);
+ sIndicator(porosity, physR);
+ if (!util::nearZero(porosity[0])) {
+ velocity(foo, physR);
+ foo[0] *= porosity[0];
+ foo[1] *= porosity[0];
+ foo[2] *= porosity[0];
+ foo[3] = porosity[0];
+ get(iX, iY, iZ).template addField(foo);
+ get(iX, iY, iZ).template addField(&foo[3]);
+ porosity[0] = 1. - porosity[0];
+ *(get(iX, iY, iZ).template getFieldPointer()) *= porosity[0];
+ }
+ }
+ }
+ }
+ }
+}
+#endif
+
+} // namespace olb
+
+#endif
--
cgit v1.2.3