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/boundary/advectionDiffusionBoundaries.hh | 313 +++++++++++++++++++++++++++ 1 file changed, 313 insertions(+) create mode 100644 src/boundary/advectionDiffusionBoundaries.hh (limited to 'src/boundary/advectionDiffusionBoundaries.hh') diff --git a/src/boundary/advectionDiffusionBoundaries.hh b/src/boundary/advectionDiffusionBoundaries.hh new file mode 100644 index 0000000..7c7f391 --- /dev/null +++ b/src/boundary/advectionDiffusionBoundaries.hh @@ -0,0 +1,313 @@ +/* This file is part of the OpenLB library + * + * Copyright (C) 2008 Orestis Malaspinas, Andrea Parmigiani + * 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 ADVECTION_DIFFUSION_BOUNDARIES_HH +#define ADVECTION_DIFFUSION_BOUNDARIES_HH + +#include "advectionDiffusionBoundaries.h" +#include "dynamics/latticeDescriptors.h" +#include "core/util.h" +#include "dynamics/lbHelpers.h" + +namespace olb { + + + +//================================================================================================== +//==================== For regularized Advection Diffusion Boundary Condition ====================== +//============================================================================================ + + +// For flat Walls + +template +AdvectionDiffusionBoundariesDynamics:: +AdvectionDiffusionBoundariesDynamics( T omega_, Momenta& momenta_) + : BasicDynamics(momenta_), boundaryDynamics(omega_, momenta_) +{ +} + +template +T AdvectionDiffusionBoundariesDynamics:: +computeEquilibrium(int iPop, T rho, const T u[DESCRIPTOR::d], T uSqr) const +{ + return lbHelpers::equilibriumFirstOrder( iPop, rho, u ); +} + + +template +void AdvectionDiffusionBoundariesDynamics:: +collide(Cell& cell,LatticeStatistics& statistics) +{ + typedef DESCRIPTOR L; + typedef lbHelpers lbH; + + T dirichletTemperature = this->_momenta.computeRho(cell); + T* u = cell.template getFieldPointer(); + + std::vector unknownIndexes = util::subIndexOutgoing(); + std::vector knownIndexes = util::remainingIndexes(unknownIndexes); + + int missingNormal = 0; + + if ((L::d == 3 && L::q == 7)||(L::d == 2 && L::q == 5)) { + T sum = T(); + for (unsigned i = 0; i < knownIndexes.size(); ++i) { + sum += cell[knownIndexes[i]]; + } + + T difference = dirichletTemperature - (T) 1 - sum; // on cell there are non-shiftet values -> temperature has to be changed + + // here I know all missing and non missing f_i + for (unsigned i = 0; i < unknownIndexes.size(); ++i) { + int numOfNonNullComp = 0; + for (int iDim = 0; iDim < L::d; ++iDim) { + numOfNonNullComp += abs(descriptors::c(unknownIndexes[i],iDim)); + } + if (numOfNonNullComp == 1) { + missingNormal = unknownIndexes[i]; + // here missing diagonal directions are erased + // just the normal direction stays (D3Q7) + unknownIndexes.erase(unknownIndexes.begin() + i); + break; + + } + } + cell[missingNormal] = difference; // on cell there are non-shiftet values -> temperature has to be changed + boundaryDynamics.collide(cell, statistics); // only for D3Q7 + } else { + // part for q=19 copied from AdvectionDiffusionEdgesDynamics.collide() + // but here just all missing directions, even at border of inlet area + // has to be checked! + for (unsigned iteratePop = 0; iteratePop < unknownIndexes.size(); + ++iteratePop) { + cell[unknownIndexes[iteratePop]] = + lbH::equilibriumFirstOrder(unknownIndexes[iteratePop], dirichletTemperature, u) + - (cell[util::opposite(unknownIndexes[iteratePop])] + - lbH::equilibriumFirstOrder( + util::opposite(unknownIndexes[iteratePop]), + dirichletTemperature, u)); + } + } +} + +template +T AdvectionDiffusionBoundariesDynamics:: +getOmega() const +{ + return boundaryDynamics.getOmega(); +} + +template +void AdvectionDiffusionBoundariesDynamics:: +setOmega(T omega_) +{ + boundaryDynamics.setOmega(omega_); +} + +//================================================================= +// For 2D Corners with regularized Dynamic ============================================== +//================================================================= +template +AdvectionDiffusionCornerDynamics2D::AdvectionDiffusionCornerDynamics2D( + T omega_, Momenta& momenta_) + : BasicDynamics(momenta_), + boundaryDynamics(omega_, momenta_) +{ +} + +template +T AdvectionDiffusionCornerDynamics2D::computeEquilibrium(int iPop, T rho, const T u[DESCRIPTOR::d], T uSqr) const +{ + return lbHelpers::equilibriumFirstOrder( iPop, rho, u ); +} + + +template +void AdvectionDiffusionCornerDynamics2D::collide(Cell& cell,LatticeStatistics& statistics) +{ + typedef DESCRIPTOR L; + typedef lbHelpers lbH; + + T temperature = this->_momenta.computeRho(cell); + T* u = cell.template getFieldPointer(); + // I need to get Missing information on the corners !!!! + std::vector unknownIndexes = util::subIndexOutgoing2DonCorners(); + // here I know all missing and non missing f_i + + + // The collision procedure for D2Q5 and D3Q7 lattice is the same ... + // Given the rule f_i_neq = -f_opposite(i)_neq + // I have the right number of equations for the number of unknowns using these lattices + + for (unsigned iPop = 0; iPop < unknownIndexes.size(); ++iPop) { + cell[unknownIndexes[iPop]] = lbH::equilibriumFirstOrder(unknownIndexes[iPop], temperature, u) + -(cell[util::opposite(unknownIndexes[iPop])] + - lbH::equilibriumFirstOrder(util::opposite(unknownIndexes[iPop]), temperature, u) ) ; + } + + // Once all the f_i are known, I can call the collision for the Regularized Model. + boundaryDynamics.collide(cell, statistics); + +} + +template +T AdvectionDiffusionCornerDynamics2D::getOmega() const +{ + return boundaryDynamics.getOmega(); +} + +template +void AdvectionDiffusionCornerDynamics2D::setOmega(T omega_) +{ + boundaryDynamics.setOmega(omega_); +} + + + +//================================================================= +// For 3D Corners with regularized Dynamic ============================================== +//================================================================= +template +AdvectionDiffusionCornerDynamics3D::AdvectionDiffusionCornerDynamics3D( + T omega_, Momenta& momenta_) + : BasicDynamics(momenta_), + boundaryDynamics(omega_, momenta_) +{ +} + +template +T AdvectionDiffusionCornerDynamics3D::computeEquilibrium(int iPop, T rho, const T u[DESCRIPTOR::d], T uSqr) const +{ + return lbHelpers::equilibriumFirstOrder( iPop, rho, u ); +} + + +template +void AdvectionDiffusionCornerDynamics3D::collide(Cell& cell,LatticeStatistics& statistics) +{ + typedef DESCRIPTOR L; + typedef lbHelpers lbH; + + T temperature = this->_momenta.computeRho(cell); + T* u = cell.template getFieldPointer(); + // I need to get Missing information on the corners !!!! + std::vector unknownIndexes = util::subIndexOutgoing3DonCorners(); + // here I know all missing and non missing f_i + + + // The collision procedure for D2Q5 and D3Q7 lattice is the same ... + // Given the rule f_i_neq = -f_opposite(i)_neq + // I have the right number of equations for the number of unknowns using these lattices + + for (unsigned iPop = 0; iPop < unknownIndexes.size(); ++iPop) { + cell[unknownIndexes[iPop]] = lbH::equilibriumFirstOrder(unknownIndexes[iPop], temperature, u) + -(cell[util::opposite(unknownIndexes[iPop])] + - lbH::equilibriumFirstOrder(util::opposite(unknownIndexes[iPop]), temperature, u) ) ; + } + + // Once all the f_i are known, I can call the collision for the Regularized Model. + boundaryDynamics.collide(cell, statistics); + +} + +template +T AdvectionDiffusionCornerDynamics3D::getOmega() const +{ + return boundaryDynamics.getOmega(); +} + +template +void AdvectionDiffusionCornerDynamics3D::setOmega(T omega_) +{ + boundaryDynamics.setOmega(omega_); +} + +//================================================================= +// For 3D Edges with regularized Dynamic ============================================== +//================================================================= +template +AdvectionDiffusionEdgesDynamics::AdvectionDiffusionEdgesDynamics( + T omega_, Momenta& momenta_) + : BasicDynamics(momenta_), + boundaryDynamics(omega_, momenta_) +{ +} + +template +T AdvectionDiffusionEdgesDynamics::computeEquilibrium(int iPop, T rho, const T u[DESCRIPTOR::d], T uSqr) const +{ + return lbHelpers::equilibriumFirstOrder( iPop, rho, u ); +} + + +template +void AdvectionDiffusionEdgesDynamics::collide(Cell& cell,LatticeStatistics& statistics) +{ + typedef DESCRIPTOR L; + typedef lbHelpers lbH; + + T temperature = this->_momenta.computeRho(cell); + T* u = cell.template getFieldPointer(); + // I need to get Missing information on the corners !!!! + std::vector unknownIndexes = util::subIndexOutgoing3DonEdges(); + // here I know all missing and non missing f_i + + + // The collision procedure for D2Q5 and D3Q7 lattice is the same ... + // Given the rule f_i_neq = -f_opposite(i)_neq + // I have the right number of equations for the number of unknowns using these lattices + + for (unsigned iPop = 0; iPop < unknownIndexes.size(); ++iPop) { + cell[unknownIndexes[iPop]] = lbH::equilibriumFirstOrder(unknownIndexes[iPop], temperature, u) + -(cell[util::opposite(unknownIndexes[iPop])] + - lbH::equilibriumFirstOrder(util::opposite(unknownIndexes[iPop]), temperature, u) ) ; + } + + // Once all the f_i are known, I can call the collision for the Regularized Model. + boundaryDynamics.collide(cell, statistics); + +} + +template +T AdvectionDiffusionEdgesDynamics::getOmega() const +{ + return boundaryDynamics.getOmega(); +} + +template +void AdvectionDiffusionEdgesDynamics::setOmega(T omega_) +{ + boundaryDynamics.setOmega(omega_); +} + + + + +} // namespace olb + + + + +#endif -- cgit v1.2.3