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 --- ...shanChenForcedSingleComponentPostProcessor2D.hh | 171 +++++++++++++++++++++ 1 file changed, 171 insertions(+) create mode 100644 src/dynamics/shanChenForcedSingleComponentPostProcessor2D.hh (limited to 'src/dynamics/shanChenForcedSingleComponentPostProcessor2D.hh') diff --git a/src/dynamics/shanChenForcedSingleComponentPostProcessor2D.hh b/src/dynamics/shanChenForcedSingleComponentPostProcessor2D.hh new file mode 100644 index 0000000..1c0cfb7 --- /dev/null +++ b/src/dynamics/shanChenForcedSingleComponentPostProcessor2D.hh @@ -0,0 +1,171 @@ +/* This file is part of the OpenLB library + * + * Copyright (C) 2008 Orestis Malaspinas, Andrea Parmigiani, 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. +*/ + +#ifndef SHAN_CHEN_FORCED_SINGLE_COMPONENT_POST_PROCESSOR_2D_HH +#define SHAN_CHEN_FORCED_SINGLE_COMPONENT_POST_PROCESSOR_2D_HH + +#include "shanChenForcedSingleComponentPostProcessor2D.h" +#include "interactionPotential.h" +#include "core/blockLattice2D.h" +#include "core/util.h" +#include "core/finiteDifference2D.h" + +namespace olb { + +//////// ShanChenForcedSingleComponentPostProcessor2D /////////////////////////////////// + + +template +ShanChenForcedSingleComponentPostProcessor2D :: +ShanChenForcedSingleComponentPostProcessor2D(int x0_, int x1_, int y0_, int y1_, T G_, + std::vector rho0_, + AnalyticalF1D& iP_, + std::vector partners_) + : x0(x0_), x1(x1_), y0(y0_), y1(y1_), G(G_), rho0(rho0_), interactionPotential(iP_), partners(partners_) +{ } + +template +ShanChenForcedSingleComponentPostProcessor2D :: +ShanChenForcedSingleComponentPostProcessor2D(T G_, + std::vector rho0_, + AnalyticalF1D& iP_, + std::vector partners_) + : x0(0), x1(0), y0(0), y1(0), G(G_), rho0(rho0_), interactionPotential(iP_), partners(partners_) +{ } + +template +void ShanChenForcedSingleComponentPostProcessor2D:: +processSubDomain( BlockLattice2D& blockLattice, + int x0_, int x1_, int y0_, int y1_ ) +{ + typedef DESCRIPTOR L; + + int newX0, newX1, newY0, newY1; + if ( util::intersect ( x0, x1, y0, y1, + x0_, x1_, y0_, y1_, + newX0, newX1, newY0, newY1 ) ) { + int nx = newX1-newX0+3; // include a one-cell boundary + int ny = newY1-newY0+3; // include a one-cell boundary + int offsetX = newX0-1; + int offsetY = newY0-1; + + BlockData2D rhoField1(nx, ny); + + // Compute density and velocity on every site of first lattice, and store result + // in external scalars; envelope cells are included, because they are needed + // to compute the interaction potential in what follows. + for (int iX=newX0-1; iX<=newX1+1; ++iX) { + for (int iY=newY0-1; iY<=newY1+1; ++iY) { + Cell& cell = blockLattice.get(iX,iY); + rhoField1.get(iX-offsetX, iY-offsetY) = cell.computeRho()*rho0[0]; + } + } + + for (int iX=newX0; iX<=newX1; ++iX) { + for (int iY=newY0; iY<=newY1; ++iY) { + Cell& blockCell = blockLattice.get(iX,iY); + + T* j = blockCell.template getFieldPointer(); + lbHelpers::computeJ(blockCell,j); + + T blockOmega = blockCell.getDynamics()->getOmega(); + + // Computation of the common velocity, shared among the two populations + T rhoTot = rhoField1.get(iX-offsetX, iY-offsetY)*blockOmega; + + T uTot[DESCRIPTOR::d]; + T *blockU = blockCell.template getFieldPointer(); // contains precomputed value rho*u + for (int iD = 0; iD < DESCRIPTOR::d; ++iD) { + uTot[iD] = (blockU[iD]*rho0[0]*blockOmega) / rhoTot; + } + + // Computation of the interaction potential + T rhoBlockContribution[L::d] = {T(), T()}; + T psi; + interactionPotential(&psi, &rhoField1.get(iX-offsetX, iY-offsetY)); + for (int iPop = 0; iPop < L::q; ++iPop) { + int nextX = iX + descriptors::c(iPop,0); + int nextY = iY + descriptors::c(iPop,1); + T blockRho; + interactionPotential(&blockRho, &rhoField1.get(nextX-offsetX, nextY-offsetY));//rho0[0]; + for (int iD = 0; iD < L::d; ++iD) { + rhoBlockContribution[iD] += psi * blockRho * descriptors::c(iPop,iD)* descriptors::t(iPop); + } + } + + // Computation and storage of the final velocity, consisting + // of u and the momentum difference due to interaction + // potential plus external force + T *blockForce = blockCell.template getFieldPointer(); + T *externalBlockForce = blockCell.template getFieldPointer(); + + for (int iD = 0; iD < L::d; ++iD) { + blockU[iD] = uTot[iD]; + blockForce[iD] = externalBlockForce[iD] - G*rhoBlockContribution[iD]/rhoField1.get(iX-offsetX, iY-offsetY); + } + } + } + } +} + +template +void ShanChenForcedSingleComponentPostProcessor2D:: +process(BlockLattice2D& blockLattice) +{ + processSubDomain(blockLattice, x0, x1, y0, y1); +} + + +/// LatticeCouplingGenerator for NS coupling + +template +ShanChenForcedSingleComponentGenerator2D::ShanChenForcedSingleComponentGenerator2D ( + int x0_, int x1_, int y0_, int y1_, T G_, std::vector rho0_, AnalyticalF1D& iP_ ) + : LatticeCouplingGenerator2D(x0_, x1_, y0_, y1_), G(G_), rho0(rho0_), interactionPotential(iP_) +{ } + +template +ShanChenForcedSingleComponentGenerator2D::ShanChenForcedSingleComponentGenerator2D ( + T G_, std::vector rho0_, AnalyticalF1D& iP_ ) + : LatticeCouplingGenerator2D(0, 0, 0, 0), G(G_), rho0(rho0_), interactionPotential(iP_) +{ } + +template +PostProcessor2D* ShanChenForcedSingleComponentGenerator2D::generate ( + std::vector partners) const +{ + return new ShanChenForcedSingleComponentPostProcessor2D( + this->x0,this->x1,this->y0,this->y1,G, rho0, interactionPotential, partners); +} + +template +LatticeCouplingGenerator2D* ShanChenForcedSingleComponentGenerator2D::clone() const +{ + return new ShanChenForcedSingleComponentGenerator2D(*this); +} + + + +} // namespace olb + +#endif -- cgit v1.2.3