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/functors/lattice/superLatticeLocalF2D.hh | 994 +++++++++++++++++++++++++++
1 file changed, 994 insertions(+)
create mode 100644 src/functors/lattice/superLatticeLocalF2D.hh
(limited to 'src/functors/lattice/superLatticeLocalF2D.hh')
diff --git a/src/functors/lattice/superLatticeLocalF2D.hh b/src/functors/lattice/superLatticeLocalF2D.hh
new file mode 100644
index 0000000..15473cd
--- /dev/null
+++ b/src/functors/lattice/superLatticeLocalF2D.hh
@@ -0,0 +1,994 @@
+/* This file is part of the OpenLB library
+ *
+ * Copyright (C) 2014 Albert Mink, Mathias J. Krause
+ * 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.
+ */
+
+#ifndef SUPER_LATTICE_LOCAL_F_2D_HH
+#define SUPER_LATTICE_LOCAL_F_2D_HH
+
+#include // for generic i/o
+#include // for lpnorm
+#include
+
+#include "superLatticeLocalF2D.h"
+#include "blockLatticeLocalF2D.h"
+#include "dynamics/lbHelpers.h" // for computation of lattice rho and velocity
+#include "geometry/superGeometry2D.h"
+#include "indicator/superIndicatorF2D.h"
+
+namespace olb {
+
+template
+SuperLatticeDissipation2D::SuperLatticeDissipation2D(
+ SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticeF2D(sLattice, 1), _converter(converter)
+{
+ this->getName() = "dissipation";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticeDissipation2D(this->_sLattice.getBlockLattice(iC),this->_converter));
+ }
+}
+
+template
+bool SuperLatticeDissipation2D::operator()(T output[],
+ const int input[])
+{
+
+ // int globIC = input[0];
+ // int locix= input[1];
+ // int lociy= input[2];
+ //// int lociz= input[3];
+ // if ( this->sLattice.getLoadBalancer().rank(globIC) == singleton::mpi().getRank() ) {
+ // // local coords are given, fetch local cell and compute value(s)
+ // T rho, uTemp[DESCRIPTOR::d], pi[util::TensorVal::n];
+ // int overlap = this->sLattice.getOverlap();
+ // this->sLattice.getExtendedBlockLattice(this->sLattice.getLoadBalancer().loc(globIC)).get(locix+overlap, lociy+overlap/*, lociz+overlap*/).computeAllMomenta(rho, uTemp, pi);
+
+ // T PiNeqNormSqr = pi[0]*pi[0] + 2.*pi[1]*pi[1] + pi[2]*pi[2];
+ // if (util::TensorVal::n == 6)
+ // PiNeqNormSqr += pi[2]*pi[2] + pi[3]*pi[3] + 2.*pi[4]*pi[4] +pi[5]*pi[5];
+
+ // T nuLattice = converter.getLatticeNu();
+ // T omega = converter.getOmega();
+ // T finalResult = PiNeqNormSqr*nuLattice*pow(omega*descriptors::invCs2(),2)/rho/2.;
+
+ // return std::vector(1,finalResult);
+ // } else {
+ // return std::vector(); // empty vector
+ // }
+ return false;
+}
+
+template
+SuperLatticePhysDissipation2D::SuperLatticePhysDissipation2D(
+ SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 1)
+{
+ this->getName() = "physDissipation";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticePhysDissipation2D(this->_sLattice.getBlockLattice(iC),this->_converter));
+ }
+}
+
+template
+bool SuperLatticePhysDissipation2D::operator()(
+ T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticeDensity2D::SuperLatticeDensity2D(
+ SuperLattice2D& sLattice)
+ : SuperLatticeF2D(sLattice, 1)
+{
+ this->getName() = "density";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back( new BlockLatticeDensity2D(this->_sLattice.getBlockLattice(iC)) );
+ }
+}
+
+template
+bool SuperLatticeDensity2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticeVelocity2D::SuperLatticeVelocity2D(
+ SuperLattice2D& sLattice)
+ : SuperLatticeF2D(sLattice, 2)
+{
+ this->getName() = "velocity";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticeVelocity2D(this->_sLattice.getBlockLattice(iC)));
+ }
+}
+
+template
+bool SuperLatticeVelocity2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysStrainRate2D::SuperLatticePhysStrainRate2D(
+ SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 4)
+{
+ this->getName() = "physStrainRate";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticePhysStrainRate2D(this->_sLattice.getBlockLattice(iC),this->_converter));
+ }
+}
+
+template
+bool SuperLatticePhysStrainRate2D::operator()(
+ T output[], const int input[])
+{
+
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysWallShearStress2D::SuperLatticePhysWallShearStress2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ const int material, const UnitConverter& converter,
+ IndicatorF2D& indicator)
+ : SuperLatticePhysF2D(sLattice, converter, 1),
+ _superGeometry(superGeometry), _material(material)
+{
+ this->getName() = "physWallShearStress";
+ const int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(
+ new BlockLatticePhysWallShearStress2D(
+ this->_sLattice.getExtendedBlockLattice(iC),
+ _superGeometry.getExtendedBlockGeometry(iC),
+ this->_sLattice.getOverlap(),
+ _material,
+ this->_converter,
+ indicator)
+ );
+ }
+}
+
+template
+bool SuperLatticePhysWallShearStress2D::operator() (T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticeGeometry2D::SuperLatticeGeometry2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ const int material)
+ : SuperLatticeF2D(sLattice, 1), _superGeometry(superGeometry),
+ _material(material)
+{
+ this->getName() = "geometry";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back( new BlockLatticeGeometry2D(
+ this->_sLattice.getBlockLattice(iC),
+ this->_superGeometry.getBlockGeometry(iC),
+ _material) );
+ }
+}
+
+template
+bool SuperLatticeGeometry2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticeRank2D::SuperLatticeRank2D(
+ SuperLattice2D& sLattice)
+ : SuperLatticeF2D(sLattice, 1)
+{
+ this->getName() = "rank";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back( new BlockLatticeRank2D(this->_sLattice.getBlockLattice(iC)) );
+ }
+}
+
+template
+bool SuperLatticeRank2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ this->getBlockF( this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ return true;
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticeCuboid2D::SuperLatticeCuboid2D(
+ SuperLattice2D& sLattice)
+ : SuperLatticeF2D(sLattice, 1)
+{
+ this->getName() = "cuboid";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back( new BlockLatticeCuboid2D(this->_sLattice.getBlockLattice(iC),
+ this->_sLattice.getLoadBalancer().glob(iC)) );
+ }
+}
+
+template
+bool SuperLatticeCuboid2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ this->getBlockF( this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ return true;
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysPressure2D::SuperLatticePhysPressure2D(
+ SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 1)
+{
+ this->getName() = "physPressure";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticePhysPressure2D(this->_sLattice.getBlockLattice(iC), this->_converter));
+ }
+}
+
+template
+bool SuperLatticePhysPressure2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF( this->_sLattice.getLoadBalancer().loc(input[0]) )(output, &input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysVelocity2D::SuperLatticePhysVelocity2D(
+ SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 2)
+{
+ this->getName() = "physVelocity";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticePhysVelocity2D(this->_sLattice.getBlockLattice(iC),
+ this->_converter));
+ }
+}
+
+template
+bool SuperLatticePhysVelocity2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+
+
+}
+
+
+template
+SuperLatticeField2D::SuperLatticeField2D(
+ SuperLattice2D& sLattice)
+ : SuperLatticeF2D(sLattice, DESCRIPTOR::template size())
+{
+ this->getName() = "ExtField";
+ for (int iC = 0; iC < this->_sLattice.getLoadBalancer().size(); iC++ ) {
+ this->_blockF.emplace_back(
+ new BlockLatticeField2D(this->_sLattice.getBlockLattice(iC)));
+ }
+}
+
+template
+bool SuperLatticeField2D::operator()(
+ T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().isLocal(input[0])) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]))(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysExternalPorosity2D::SuperLatticePhysExternalPorosity2D
+(SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 1)
+{
+ this->getName() = "ExtPorosityField";
+ const int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(
+ new BlockLatticePhysExternalPorosity2D(
+ this->_sLattice.getExtendedBlockLattice(iC),
+ this->_sLattice.getOverlap(),
+ this->_converter)
+ );
+ }
+}
+
+template
+bool SuperLatticePhysExternalPorosity2D::operator()(
+ T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ const int loc = this->_sLattice.getLoadBalancer().loc(input[0]);
+ return this->getBlockF(loc)(output, &input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticeVolumeFractionApproximation2D::SuperLatticeVolumeFractionApproximation2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ IndicatorF2D& indicator, int refinementLevel,
+ const UnitConverter& converter, bool insideOut)
+ : SuperLatticeF2D(sLattice, 1)
+{
+ this->getName() = "volumeFractionApproximation";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticeVolumeFractionApproximation2D(this->_sLattice.getBlockLattice(iC),
+ superGeometry.getBlockGeometry(iC),
+ indicator, refinementLevel,
+ converter, insideOut));
+ }
+}
+
+template
+bool SuperLatticeVolumeFractionApproximation2D::operator()(
+ T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticeVolumeFractionPolygonApproximation2D::SuperLatticeVolumeFractionPolygonApproximation2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ IndicatorF2D& indicator, const UnitConverter& converter, bool insideOut)
+ : SuperLatticeF2D(sLattice, 1)
+{
+ this->getName() = "volumeFractionPolygonApproximation";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticeVolumeFractionPolygonApproximation2D(this->_sLattice.getBlockLattice(iC),
+ superGeometry.getBlockGeometry(iC),
+ indicator, converter, insideOut));
+ }
+}
+
+template
+bool SuperLatticeVolumeFractionPolygonApproximation2D::operator()(
+ T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+
+template
+SuperLatticePhysExternalVelocity2D::SuperLatticePhysExternalVelocity2D(
+ SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 2)
+{
+ this->getName() = "ExtVelocityField";
+ const int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(
+ new BlockLatticePhysExternalVelocity2D(
+ this->_sLattice.getExtendedBlockLattice(iC),
+ this->_sLattice.getOverlap(),
+ this->_converter)
+ );
+ }
+}
+
+template
+bool SuperLatticePhysExternalVelocity2D::operator()(
+ T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ const int loc = this->_sLattice.getLoadBalancer().loc(input[0]);
+ return this->getBlockF(loc)(output, &input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysExternalParticleVelocity2D::SuperLatticePhysExternalParticleVelocity2D(
+ SuperLattice2D& sLattice, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 2)
+{
+ this->getName() = "ExtPartVelField";
+
+ for (int iC = 0; iC < sLattice.getLoadBalancer().size(); ++iC) {
+ this->_blockF.emplace_back(
+ new BlockLatticePhysExternalParticleVelocity2D(
+ sLattice.getExtendedBlockLattice(iC),
+ converter)
+ );
+ }
+}
+
+template
+bool SuperLatticePhysExternalParticleVelocity2D::operator()(
+ T output[], const int input[])
+{
+ auto& load = this->_sLattice.getLoadBalancer();
+ const int& globIC = input[0];
+
+ if (load.rank(globIC) == singleton::mpi().getRank()) {
+ const int overlap = this->_sLattice.getOverlap();
+
+ int inputLocal[2] = { };
+ inputLocal[0] = input[1] + overlap;
+ inputLocal[1] = input[2] + overlap;
+
+ return this->getBlockF(load.loc(globIC))(output, inputLocal);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysBoundaryForce2D::SuperLatticePhysBoundaryForce2D(
+ SuperLattice2D& sLattice,
+ FunctorPtr>&& indicatorF,
+ const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 2),
+ _indicatorF(std::move(indicatorF))
+{
+ this->getName() = "physBoundaryForce";
+ for (int iC = 0; iC < this->_sLattice.getLoadBalancer().size(); ++iC) {
+ this->_blockF.emplace_back(
+ new BlockLatticePhysBoundaryForce2D(
+ this->_sLattice.getBlockLattice(iC),
+ _indicatorF->getBlockIndicatorF(iC),
+ this->_converter));
+ }
+}
+
+template
+SuperLatticePhysBoundaryForce2D::SuperLatticePhysBoundaryForce2D(
+ SuperLattice2D& sLattice,
+ SuperGeometry2D& superGeometry, const int material,
+ const UnitConverter& converter)
+ : SuperLatticePhysBoundaryForce2D(sLattice,
+ superGeometry.getMaterialIndicator(material),
+ converter)
+{ }
+
+template
+bool SuperLatticePhysBoundaryForce2D::operator() (T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]))(output, &input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePSMPhysForce2D::SuperLatticePSMPhysForce2D(
+ SuperLattice2D& sLattice,
+ const UnitConverter& converter,
+ int mode_)
+ : SuperLatticePhysF2D(sLattice, converter, 2)
+{
+ this->getName() = "PSMPhysForce";
+ for (int iC = 0; iC < this->_sLattice.getLoadBalancer().size(); ++iC) {
+ this->_blockF.emplace_back(
+ new BlockLatticePSMPhysForce2D(
+ this->_sLattice.getBlockLattice(iC),
+ this->_converter,
+ mode_));
+ }
+}
+
+template
+bool SuperLatticePSMPhysForce2D::operator() (T output[], const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]))(output, &input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+/*****************************NEW*****************************/
+
+template
+SuperLatticePhysCorrBoundaryForce2D::SuperLatticePhysCorrBoundaryForce2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ const int material, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 2),
+ _superGeometry(superGeometry), _material(material)
+{
+ this->getName() = "physCorrBoundaryForce";
+}
+
+template
+bool SuperLatticePhysCorrBoundaryForce2D::operator()(
+ T output[], const int input[])
+{
+ // int globIC = input[0];
+ // int locix= input[1];
+ // int lociy= input[2];
+ // int lociz= input[3];
+
+ // std::vector force(3, T());
+ // if ( this->sLattice.getLoadBalancer().rank(globIC) == singleton::mpi().getRank() ) {
+ // int globX = (int)this->sLattice.getCuboidGeometry().get(globIC).get_globPosX() + locix;
+ // int globY = (int)this->sLattice.getCuboidGeometry().get(globIC).get_globPosY() + lociy;
+ // int globZ = (int)this->sLattice.getCuboidGeometry().get(globIC).get_globPosZ() + lociz;
+ // if(superGeometry.getMaterial(globX, globY, globZ) == material) {
+ // for (int iPop = 1; iPop < DESCRIPTOR::q ; ++iPop) {
+ // // Get direction
+ // const int* c = DESCRIPTOR::c[iPop];
+ // // Get next cell located in the current direction
+ // // Check if the next cell is a fluid node
+ // if (superGeometry.getMaterial(globX + c[0], globY + c[1], globZ + c[2]) == 1) {
+ // int overlap = this->sLattice.getOverlap();
+ // // Get f_q of next fluid cell where l = opposite(q)
+ // T f = this->sLattice.getExtendedBlockLattice(this->sLattice.getLoadBalancer().loc(globIC)).get(locix+overlap + c[0], lociy+overlap + c[1], lociz+overlap + c[2])[iPop];
+ // // Get f_l of the boundary cell
+ // // Add f_q and f_opp
+ // f += this->sLattice.getExtendedBlockLattice(this->sLattice.getLoadBalancer().loc(globIC)).get(locix+overlap, lociy+overlap, lociz+overlap)[util::opposite(iPop)];
+ // // Update force
+ // force[0] -= c[0]*(f-2.*descriptors::t(iPop));
+ // force[1] -= c[1]*(f-2.*descriptors::t(iPop));
+ // force[2] -= c[2]*(f-2.*descriptors::t(iPop));
+ // }
+ // }
+ // force[0]=this->converter.physForce(force[0]);
+ // force[1]=this->converter.physForce(force[1]);
+ // force[2]=this->converter.physForce(force[2]);
+ // return force;
+ // }
+ // else {
+ // return force;
+ // }
+ // }
+ // else {
+ // return force;
+ // }
+ return false;
+}
+
+template
+SuperLatticePorosity2D::SuperLatticePorosity2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ const int material, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 1),
+ _superGeometry(superGeometry), _material(material)
+{
+ this->getName() = "porosity";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticePorosity2D(
+ this->_sLattice.getBlockLattice(iC),
+ this->_superGeometry.getBlockGeometry(iC),
+ _material,
+ converter));
+ }
+}
+
+template
+bool SuperLatticePorosity2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysPermeability2D::SuperLatticePhysPermeability2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ const int material, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 1),
+ _superGeometry(superGeometry), _material(material)
+{
+ this->getName() = "permeability";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back( new BlockLatticePhysPermeability2D(
+ this->_sLattice.getBlockLattice(iC),
+ this->_superGeometry.getBlockGeometry(iC),
+ _material,
+ this->getConverter() ) );
+ }
+}
+
+template
+bool SuperLatticePhysPermeability2D::operator()(
+ T output[], const int input[])
+{
+ // int globIC = input[0];
+ // int locix= input[1];
+ // int lociy= input[2];
+ //// int lociz= input[3];
+
+ // T* value = new T[1];
+ // int overlap = this->sLattice.getOverlap();
+ // this->sLattice.getExtendedBlockLattice(this->sLattice.getLoadBalancer().loc(globIC)).get(locix+overlap, lociy+overlap/*, lociz+overlap*/).computeField(0,1,value);
+ // std::vector result(1,this->converter.physPermeability(value[0]));
+ // delete value;
+ // if (!(result[0]<42)&&!(result[0]>42)&&!(result[0]==42)) result[0] = 999999;
+ // if (isinf(result[0])) result[0] = 1e6;
+ // return result;
+ return false;
+}
+
+template
+SuperLatticePhysDarcyForce2D::SuperLatticePhysDarcyForce2D(
+ SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ const int material, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice, converter, 2),
+ _superGeometry(superGeometry), _material(material)
+{
+ this->getName() = "alphaU";
+}
+
+template
+bool SuperLatticePhysDarcyForce2D::operator()(
+ T output[], const int input[])
+{
+ // SuperLatticePhysPermeability2D permeability(this->sLattice,this->superGeometry,this->material,this->converter);
+ // SuperLatticeVelocity2D velocity(this->sLattice);
+
+ // T nu = this->converter.getCharNu();
+ // T K = permeability(input)[0];
+ // std::vector u = velocity(input);
+
+ // std::vector result(2,T());
+ // result[0] = -nu/K*u[0];
+ // result[1] = -nu/K*u[1];
+ //// result[2] = -nu/K*u[2];
+
+ // return result;
+ return false;
+}
+
+template
+SuperEuklidNorm2D::SuperEuklidNorm2D(
+ SuperLatticeF2D& f)
+ : SuperLatticeF2D(f.getSuperLattice(), 1), _f(f)
+{
+ this->getName() = "l2(" + _f.getName() + ")";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockEuklidNorm2D(f.getBlockF(iC)));
+ }
+}
+
+template
+bool SuperEuklidNorm2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysTemperature2D::SuperLatticePhysTemperature2D(
+ SuperLattice2D& sLattice, ThermalUnitConverter const& converter)
+ : SuperLatticeThermalPhysF2D(sLattice, converter, 1)
+{
+ this->getName() = "physTemperature";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticePhysTemperature2D(this->_sLattice.getBlockLattice(iC), this->_converter));
+ }
+}
+
+template
+bool SuperLatticePhysTemperature2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF( this->_sLattice.getLoadBalancer().loc(input[0]) )(output, &input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePhysHeatFlux2D::SuperLatticePhysHeatFlux2D(
+ SuperLattice2D& sLattice, const ThermalUnitConverter& converter)
+ : SuperLatticeThermalPhysF2D(sLattice, converter, 2)
+{
+ this->getName() = "physHeatFlux";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back(new BlockLatticePhysHeatFlux2D(this->_sLattice.getBlockLattice(iC),
+ this->_converter));
+ }
+}
+
+template
+bool SuperLatticePhysHeatFlux2D::operator()(T output[],
+ const int input[])
+{
+ if (this->_sLattice.getLoadBalancer().rank(input[0]) == singleton::mpi().getRank()) {
+ return this->getBlockF(this->_sLattice.getLoadBalancer().loc(input[0]) )(output,&input[1]);
+ }
+ else {
+ return false;
+ }
+}
+
+template
+SuperLatticePorousMomentumLossForce2D::SuperLatticePorousMomentumLossForce2D
+(SuperLattice2D& sLattice, SuperGeometry2D& superGeometry,
+ std::vector* >& indicator, const UnitConverter& converter)
+ : SuperLatticePhysF2D(sLattice,converter,4*indicator.size())
+{
+ this->getName() = "physPorousMomentumLossForce";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back( new BlockLatticePorousMomentumLossForce2D(this->_sLattice.getBlockLattice(iC), superGeometry.getBlockGeometry(iC), indicator, converter));
+ }
+}
+
+template
+bool SuperLatticePorousMomentumLossForce2D::operator() (T output[],
+ const int input[])
+{
+ for (int i=0; igetTargetDim(); i++) {
+ output[i] = 0.;
+ }
+ for (int iC = 0; iC < this->_sLattice.getLoadBalancer().size(); ++iC) {
+ int globiC = this->_sLattice.getLoadBalancer().glob(iC);
+ if ( this->_sLattice.getLoadBalancer().rank(globiC) == singleton::mpi().getRank() ) {
+ this->getBlockF(iC)(output,&input[1]);
+ }
+ }
+
+#ifdef PARALLEL_MODE_MPI
+ for (int i = 0; i < this->getTargetDim(); ++i) {
+ singleton::mpi().reduceAndBcast(output[i], MPI_SUM);
+ }
+#endif
+ return true;
+
+}
+
+template
+SuperLatticeIndicatorSmoothIndicatorIntersection2D::SuperLatticeIndicatorSmoothIndicatorIntersection2D (
+ SuperLattice2D& sLattice,
+ SuperGeometry2D& superGeometry,
+ IndicatorF2D& normalInd, SmoothIndicatorF2D& smoothInd )
+ : SuperLatticeF2D(sLattice, 1)
+{
+ this->getName() = "Indicator-SmoothIndicator Intersection";
+ int maxC = this->_sLattice.getLoadBalancer().size();
+ this->_blockF.reserve(maxC);
+ for (int iC = 0; iC < maxC; iC++) {
+ this->_blockF.emplace_back( new BlockLatticeIndicatorSmoothIndicatorIntersection2D