OpenLB 1.7
Loading...
Searching...
No Matches
setRtlbmDirectedTemperatureBoundary3D.hh
Go to the documentation of this file.
1/* This file is part of the OpenLB library
2 *
3 * Copyright (C) 2020 Alexander Schulz
4 * E-mail contact: info@openlb.net
5 * The most recent release of OpenLB can be downloaded at
6 * <http://www.openlb.net/>
7 *
8 * This program is free software; you can redistribute it and/or
9 * modify it under the terms of the GNU General Public License
10 * as published by the Free Software Foundation; either version 2
11 * of the License, or (at your option) any later version.
12 *
13 * This program is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 * GNU General Public License for more details.
17 *
18 * You should have received a copy of the GNU General Public
19 * License along with this program; if not, write to the Free
20 * Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
21 * Boston, MA 02110-1301, USA.
22*/
23
24//This file contains the RtlbmDirectedTemperature Boundary
25//This is a new version of the Boundary, which only contains free floating functions
26#ifndef SET_RTLBMDIRECTED_TEMPERATURE_BOUNDARY_3D_HH
27#define SET_RTLBMDIRECTED_TEMPERATURE_BOUNDARY_3D_HH
28
30
31namespace olb {
32
34template<typename T, typename DESCRIPTOR>
35void setRtlbmDirectedTemperatureBoundary(SuperLattice<T, DESCRIPTOR>& sLattice, T omega, SuperGeometry<T,3>& superGeometry, int material)
36{
37 setRtlbmDirectedTemperatureBoundary<T,DESCRIPTOR>(sLattice, omega, superGeometry.getMaterialIndicator(material));
38}
39
41template<typename T, typename DESCRIPTOR>
43{
44 int _overlap = 1;
45 OstreamManager clout(std::cout, "setRtlbmDirectedTemperatureBoundary");
46 bool includeOuterCells = false;
47 if (indicator->getSuperGeometry().getOverlap() == 1) {
48 includeOuterCells = true;
49 clout << "WARNING: overlap == 1, boundary conditions set on overlap despite unknown neighbor materials" << std::endl;
50 }
51 for (int iCloc = 0; iCloc < sLattice.getLoadBalancer().size(); iCloc++) {
52 setRtlbmDirectedTemperatureBoundary<T,DESCRIPTOR>(sLattice.getBlock(iCloc),
53 indicator->getBlockIndicatorF(iCloc), omega, includeOuterCells);
54 }
56 addPoints2CommBC(sLattice, std::forward<decltype(indicator)>(indicator), _overlap);
57}
58
59
61template<typename T, typename DESCRIPTOR>
63 bool includeOuterCells)
64{
65 using namespace boundaryhelper;
66 auto& blockGeometryStructure = indicator.getBlockGeometry();
67 const int margin = includeOuterCells ? 0 : 1;
68 std::vector<int> discreteNormal(4,0);
69 blockGeometryStructure.forSpatialLocations([&](auto iX, auto iY, auto iZ) {
70 if (blockGeometryStructure.getNeighborhoodRadius({iX, iY, iZ}) >= margin
71 && indicator(iX, iY, iZ)) {
72 Dynamics<T,DESCRIPTOR>* dynamics = nullptr;
73 PostProcessorGenerator3D<T,DESCRIPTOR>* postProcessor = nullptr;
74 discreteNormal = blockGeometryStructure.getStatistics().getType(iX, iY, iZ);
75 if (discreteNormal[0] == 0) { // flat //set momenta, dynamics and postProcessors for indicated cells
76 if (discreteNormal[1] != 0 && discreteNormal[1] == -1) {
77 dynamics = _block.template getDynamics<RtlbmDirectedBoundaryDynamics<T,DESCRIPTOR,
79 >>();
80 }
81 else if (discreteNormal[1] != 0 && discreteNormal[1] == 1) {
82 dynamics = _block.template getDynamics<RtlbmDirectedBoundaryDynamics<T,DESCRIPTOR,
84 >>();
85 }
86 else if (discreteNormal[2] != 0 && discreteNormal[2] == -1) {
87 dynamics = _block.template getDynamics<RtlbmDirectedBoundaryDynamics<T,DESCRIPTOR,
89 >>();
90 }
91 else if (discreteNormal[2] != 0 && discreteNormal[2] == 1) {
92 dynamics = _block.template getDynamics<RtlbmDirectedBoundaryDynamics<T,DESCRIPTOR,
94 >>();
95 }
96 else if (discreteNormal[3] != 0 && discreteNormal[3] == -1) {
97 dynamics = _block.template getDynamics<RtlbmDirectedBoundaryDynamics<T,DESCRIPTOR,
99 >>();
100 }
101 else if (discreteNormal[3] != 0 && discreteNormal[3] == 1) {
102 dynamics = _block.template getDynamics<RtlbmDirectedBoundaryDynamics<T,DESCRIPTOR,
104 >>();
105 }
106 }
107 else if (discreteNormal[0] == 1) { // corner //set momenta, dynamics and postProcessors on indicated boundery corner cells
108 dynamics = _block.getDynamics(NormalDynamicsForPlainMomenta<T,DESCRIPTOR,
110 >::construct(Vector<int,3>(discreteNormal.data() + 1)));
111 }
112 else if (discreteNormal[0] == 3) { // edge //set momenta, dynamics and postProcessors on indicated boundary edge cells
113 dynamics = _block.getDynamics(NormalSpecialDynamicsForPlainMomenta<T,DESCRIPTOR,
115 >::construct(Vector<int,3>(discreteNormal.data() + 1)));
116 }
117 dynamics->getParameters(_block).template set<descriptors::OMEGA>(omega);
118 setBoundary(_block, iX,iY,iZ, dynamics, postProcessor);
119 }
120 });
121}
122
123
124}//namespace olb
125#endif
126
Base block indicator functor.
BlockGeometry< T, 3 > & getBlockGeometry()
Get underlying block geometry structure.
Platform-abstracted block lattice for external access and inter-block interaction.
virtual Dynamics< T, DESCRIPTOR > * getDynamics(DynamicsPromise< T, DESCRIPTOR > &&)=0
Return pointer to dynamics yielded by promise.
Smart pointer for managing the various ways of passing functors around.
Definition functorPtr.h:60
class for marking output with some text
Representation of a statistic for a parallel 2D geometry.
std::unique_ptr< SuperIndicatorF< T, D > > getMaterialIndicator(std::vector< int > &&materials)
Returns a material indicator using the given vector of materials.
Base indicator functor (discrete)
Super class maintaining block lattices for a cuboid decomposition.
BlockLattice< T, DESCRIPTOR > & getBlock(int locC)
Return BlockLattice with local index locC.
LoadBalancer< T > & getLoadBalancer()
Read and write access to the load balancer.
Plain old scalar vector.
Definition vector.h:47
Tuple< FixedDensity, FixedVelocityMomentumGeneric, ZeroStress, DefineSeparately > EquilibriumBoundaryTuple
Velocity and density are stored in external fields.
Definition aliases.h:100
Top level namespace for all of OpenLB.
void addPoints2CommBC(SuperLattice< T, DESCRIPTOR > &sLattice, FunctorPtr< SuperIndicatorF2D< T > > &&indicator, int _overlap)
Adds needed Cells to the Communicator _commBC in SuperLattice.
void setRtlbmDirectedTemperatureBoundary(SuperLattice< T, DESCRIPTOR > &sLattice, T omega, SuperGeometry< T, 3 > &superGeometry, int material)
Initialising the setRtlbmDirectedTemperatureBoundary function on the superLattice domain.
void setBoundary(BlockLattice< T, DESCRIPTOR > &block, int iX, int iY, Dynamics< T, DESCRIPTOR > *dynamics, PostProcessorGenerator2D< T, DESCRIPTOR > *postProcessor)
Interface for per-cell dynamics.
Definition interface.h:54
virtual AbstractParameters< T, DESCRIPTOR > & getParameters(BlockLattice< T, DESCRIPTOR > &block)=0
Parameters access for legacy post processors.