OpenLB 1.7
Loading...
Searching...
No Matches
navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh
Go to the documentation of this file.
1/* This file is part of the OpenLB library
2 *
3 * Copyright (C) 2008 Orestis Malaspinas, Andrea Parmigiani, Jonas Latt
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#ifndef NAVIER_STOKES_INTO_ADVECTION_DIFFUSION_COUPLING_POST_PROCESSOR_3D_HH
25#define NAVIER_STOKES_INTO_ADVECTION_DIFFUSION_COUPLING_POST_PROCESSOR_3D_HH
26
27#include "latticeDescriptors.h"
28#include "navierStokesAdvectionDiffusionCouplingPostProcessor3D.h"
29#include "core/util.h"
30#include "utilities/finiteDifference3D.h"
31#include "advectionDiffusionForces.hh"
32#include "advectionDiffusionForces.h"
33
34#include "descriptorAlias.h"
35
36namespace olb {
37
38//=====================================================================================
39//============== TotalEnthalpyPhaseChangeCouplingPostProcessor3D ===============
40//=====================================================================================
41
42template<typename T, typename DESCRIPTOR, typename DYNAMICS>
43TotalEnthalpyPhaseChangeCouplingPostProcessor3D<T,DESCRIPTOR,DYNAMICS>::
44TotalEnthalpyPhaseChangeCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
45 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_,
46 std::vector<BlockStructureD<3>* > partners_)
47 : x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_),
48 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_),
49 dir(dir_), partners(partners_)
50{
51 this->getName() = "TotalEnthalpyPhaseChangeCouplingPostProcessor3D";
52 // we normalize the direction of force vector
53 T normDir = T();
54 for (unsigned iD = 0; iD < dir.size(); ++iD) {
55 normDir += dir[iD]*dir[iD];
56 }
57 normDir = util::sqrt(normDir);
58 for (unsigned iD = 0; iD < dir.size(); ++iD) {
59 dir[iD] /= normDir;
60 }
61
62 for (unsigned iD = 0; iD < dir.size(); ++iD) {
63 forcePrefactor[iD] = gravity * dir[iD];
64 }
65
66 tPartner = static_cast<BlockLattice<T,descriptors::D3Q7<descriptors::VELOCITY,descriptors::TEMPERATURE>> *>(partners[0]);
67}
68
69template<typename T, typename DESCRIPTOR, typename DYNAMICS>
70void TotalEnthalpyPhaseChangeCouplingPostProcessor3D<T,DESCRIPTOR,DYNAMICS>::
71processSubDomain(BlockLattice<T,DESCRIPTOR>& blockLattice,
72 int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
73{
74
75 int newX0, newX1, newY0, newY1, newZ0, newZ1;
76 if ( util::intersect (
77 x0, x1, y0, y1, z0, z1,
78 x0_, x1_, y0_, y1_, z0_, z1_,
79 newX0, newX1, newY0, newY1, newZ0, newZ1 ) ) {
80 auto* dynamics = static_cast<DYNAMICS*>(tPartner->template getDynamics<DYNAMICS>());
81 auto& parameters = static_cast<ParametersOfDynamicsD<DYNAMICS>&>(
82 tPartner->template getData<OperatorParameters<DYNAMICS>>());
83
84 for (int iX=newX0; iX<=newX1; ++iX) {
85 for (int iY=newY0; iY<=newY1; ++iY) {
86 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
87 auto cell = blockLattice.get(iX,iY,iZ);
88 auto partnerCell = tPartner->get(iX,iY,iZ);
89
90 T enthalpy = partnerCell.computeRho();
91
92 cell.template setField<descriptors::POROSITY>(
93 dynamics->template computeLiquidFraction<T>(parameters, enthalpy));
94 auto temperature = partnerCell.template getFieldPointer<descriptors::TEMPERATURE>();
95 temperature[0] = dynamics->template computeTemperature<T>(parameters, enthalpy);
96
97 // computation of the bousinessq force
98 auto force = blockLattice.get(iX,iY,iZ).template getFieldPointer<descriptors::FORCE>();
99 T temperatureDifference = temperature[0] - T0;
100 for (unsigned iD = 0; iD < L::d; ++iD) {
101 force[iD] = forcePrefactor[iD] * temperatureDifference;
102 }
103 // Velocity coupling
104 FieldD<T,DESCRIPTOR,descriptors::VELOCITY> u;
105 blockLattice.get(iX,iY,iZ).computeU(u.data());
106 tPartner->get(iX,iY,iZ).template setField<descriptors::VELOCITY>(u);
107 }
108 }
109 }
110 }
111}
112
113template<typename T, typename DESCRIPTOR, typename DYNAMICS>
114void TotalEnthalpyPhaseChangeCouplingPostProcessor3D<T,DESCRIPTOR,DYNAMICS>::
115process(BlockLattice<T,DESCRIPTOR>& blockLattice)
116{
117 processSubDomain(blockLattice, x0, x1, y0, y1, z0, z1);
118}
119
121
122template<typename T, typename DESCRIPTOR, typename DYNAMICS>
123TotalEnthalpyPhaseChangeCouplingGenerator3D<T,DESCRIPTOR,DYNAMICS>::
124TotalEnthalpyPhaseChangeCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
125 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_)
126 : LatticeCouplingGenerator3D<T,DESCRIPTOR>(x0_, x1_, y0_, y1_, z0_, z1_),
127 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_), dir(dir_)
128{ }
129
130template<typename T, typename DESCRIPTOR, typename DYNAMICS>
131PostProcessor3D<T,DESCRIPTOR>* TotalEnthalpyPhaseChangeCouplingGenerator3D<T,DESCRIPTOR,DYNAMICS>::generate (
132 std::vector<BlockStructureD<3>* > partners) const
133{
134 return new TotalEnthalpyPhaseChangeCouplingPostProcessor3D<T,DESCRIPTOR,DYNAMICS>(
135 this->x0,this->x1,this->y0,this->y1,this->z0,this->z1, gravity, T0, deltaTemp, dir,partners);
136}
137
138template<typename T, typename DESCRIPTOR, typename DYNAMICS>
143
144
145//=====================================================================================
146//============== PhaseFieldCouplingPostProcessor3D ===============
147//=====================================================================================
148
149template<typename T, typename DESCRIPTOR>
150PhaseFieldCouplingPostProcessor3D<T,DESCRIPTOR>::
151PhaseFieldCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
152 T rho_L, T rho_H, T mu_L, T mu_H, T surface_tension, T interface_thickness,
153 std::vector<BlockStructureD<3>* > partners_)
154 : x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_),
155 _rho_L(rho_L), _rho_H(rho_H), _delta_rho(rho_H - rho_L), _mu_L(mu_L), _mu_H(mu_H), _surface_tension(surface_tension), _interface_thickness(interface_thickness),
156 _beta(12.0 * surface_tension / interface_thickness), _kappa(1.5 * surface_tension * interface_thickness)
157{
158 this->getName() = "PhaseFieldCouplingPostProcessor3D";
159 tPartner = static_cast<BlockLattice<T,descriptors::D3Q7<descriptors::VELOCITY,descriptors::INTERPHASE_NORMAL>> *>(partners_[0]);
160}
161
162template<typename T, typename DESCRIPTOR>
163void PhaseFieldCouplingPostProcessor3D<T,DESCRIPTOR>::
164processSubDomain(BlockLattice<T,DESCRIPTOR>& blockLattice,
165 int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
166{
167
168 int newX0, newX1, newY0, newY1, newZ0, newZ1;
169 if ( util::intersect (
170 x0, x1, y0, y1, z0, z1,
171 x0_, x1_, y0_, y1_, z0_, z1_,
172 newX0, newX1, newY0, newY1, newZ0, newZ1 ) ) {
173
174 // generate phi cache
175 auto& phi_cache = blockLattice.template getField<PHI_CACHE>()[0];
176 for (int iX=newX0-1; iX<=newX1+1; ++iX) {
177 for (int iY=newY0-1; iY<=newY1+1; ++iY) {
178 for (int iZ=newZ0-1; iZ<=newZ1+1; ++iZ) {
179 phi_cache[blockLattice.getCellId(iX,iY,iZ)] = util::max(util::min(tPartner->get(iX,iY,iZ).computeRho(), 1.0), 0.0);
180 }
181 }
182 }
183
184 for (int iX=newX0; iX<=newX1; ++iX) {
185 for (int iY=newY0; iY<=newY1; ++iY) {
186 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
187 T phi = phi_cache[blockLattice.getCellId(iX, iY, iZ)];
188
189 // compute rho from phi
190 T rho = _rho_L + phi * _delta_rho;
191
192 // compute dynamic viscosity
193 T viscosity = _mu_L + phi * (_mu_H - _mu_L);
194
195 // get relaxation time
196 auto tau = blockLattice.get(iX,iY,iZ).template getFieldPointer<descriptors::TAU_EFF>();
197
198 // compute grad phi and laplace phi
199 Vector<T,3> grad_phi;
200 T laplace_phi = 0.0;
201 for (int iPop = 1; iPop < L::q; ++iPop) {
202 int nextX = iX + descriptors::c<L>(iPop,0);
203 int nextY = iY + descriptors::c<L>(iPop,1);
204 int nextZ = iZ + descriptors::c<L>(iPop,2);
205 T neighbor_phi = phi_cache[blockLattice.getCellId(nextX, nextY, nextZ)];
206
207 laplace_phi += (neighbor_phi - phi) * descriptors::t<T,L>(iPop);
208 neighbor_phi *= descriptors::t<T,L>(iPop);
209 grad_phi += neighbor_phi * descriptors::c<L>(iPop);
210 }
211 grad_phi *= descriptors::invCs2<T,L>();
212 laplace_phi *= 2.0 * descriptors::invCs2<T,L>();
213
214 // compute grad rho
215 Vector<T,3> grad_rho = _delta_rho * grad_phi;
216
217 // compute interphase normal, save to external field
218 T norm_grad_phi = util::max(norm(grad_phi), std::numeric_limits<T>::epsilon());
219 tPartner->get(iX,iY,iZ).template setField<descriptors::INTERPHASE_NORMAL>(grad_phi / norm_grad_phi);
220
221 // compute _forces (F_s, F_b, F_p, F_nu)
222 // F_s (surface tension)
223 T chemical_potential = (4.0 * _beta * (phi - 0.0) * (phi - 0.5) * (phi - 1.0)) - _kappa * laplace_phi;
224 Vector<T,3> surface_tension_force = chemical_potential * grad_phi;
225
226 // F_b (body force, e.g. bouyancy)
227 Vector<T,3> body_force;
228
229 // F_p (pressure)
230 T pressure = blockLattice.get(iX,iY,iZ).computeRho();
231 Vector<T,3> pressure_force = -pressure / descriptors::invCs2<T,L>() * grad_rho;
232
233 // F_nu (viscous)
234 Vector<T,3> viscous_force;
235 T rho_tmp, u_tmp[3];
236 blockLattice.get(iX,iY,iZ).computeRhoU( rho_tmp, u_tmp );
237 T p_tmp = rho_tmp / descriptors::invCs2<T,DESCRIPTOR>();
238 T uSqr_tmp = util::normSqr<T,DESCRIPTOR::d>(u_tmp);
239 for (int iPop = 0; iPop < L::q; ++iPop) {
240 T fEq = blockLattice.get(iX,iY,iZ).getDynamics()->computeEquilibrium( iPop, p_tmp, u_tmp );
241 T fNeq = blockLattice.get(iX,iY,iZ)[iPop] - fEq;
242 for (int iD = 0; iD < L::d; ++iD) {
243 for (int jD = 0; jD < L::d; ++jD) {
244 viscous_force[iD] += descriptors::c<L>(iPop,iD) * descriptors::c<L>(iPop,jD) * fNeq * grad_rho[jD];
245 }
246 }
247 }
248 viscous_force *= -viscosity / rho / tau[0] * descriptors::invCs2<T,L>();
249
250 // save force/rho to external field
251 blockLattice.get(iX,iY,iZ).template setField<descriptors::FORCE>(
252 (surface_tension_force + body_force + pressure_force + viscous_force) / rho
253 );
254
255 // compute u, save to external field
256 auto u = tPartner->get(iX,iY,iZ).template getField<descriptors::VELOCITY>();
257 blockLattice.get(iX,iY,iZ).computeU(u.data());
258 tPartner->get(iX,iY,iZ).template setField<descriptors::VELOCITY>(u);
259
260 // compute relaxation time, save to external field
261 tau[0] = viscosity / rho * descriptors::invCs2<T,L>() + 0.5;
262 }
263 }
264 }
265 }
266}
267
268template<typename T, typename DESCRIPTOR>
269void PhaseFieldCouplingPostProcessor3D<T,DESCRIPTOR>::
270process(BlockLattice<T,DESCRIPTOR>& blockLattice)
271{
272 processSubDomain(blockLattice, x0, x1, y0, y1, z0, z1);
273}
274
276
277template<typename T, typename DESCRIPTOR>
278PhaseFieldCouplingGenerator3D<T,DESCRIPTOR>::
279PhaseFieldCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
280 T rho_L, T rho_H, T mu_L, T mu_H, T surface_tension, T interface_thickness)
281 : LatticeCouplingGenerator3D<T,DESCRIPTOR>(x0_, x1_, y0_, y1_, z0_, z1_),
282 _rho_L(rho_L), _rho_H(rho_H), _mu_L(mu_L), _mu_H(mu_H), _surface_tension(surface_tension), _interface_thickness(interface_thickness)
283{ }
284
285template<typename T, typename DESCRIPTOR>
286PostProcessor3D<T,DESCRIPTOR>* PhaseFieldCouplingGenerator3D<T,DESCRIPTOR>::generate (
287 std::vector<BlockStructureD<3>* > partners) const
288{
289 return new PhaseFieldCouplingPostProcessor3D<T,DESCRIPTOR>(
290 this->x0,this->x1,this->y0,this->y1,this->z0,this->z1, _rho_L, _rho_H, _mu_L, _mu_H, _surface_tension, _interface_thickness, partners);
291}
292
293template<typename T, typename DESCRIPTOR>
298
299
300//=====================================================================================
301//============== SmagorinskyBoussinesqCouplingPostProcessor3D ===============
302//=====================================================================================
303
304template<typename T, typename DESCRIPTOR>
305SmagorinskyBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>::
306SmagorinskyBoussinesqCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
307 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_, T PrTurb_, T smagoPrefactor_,
308 std::vector<BlockStructureD<3>* > partners_)
309 : x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_),
310 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_),
311 dir(dir_), PrTurb(PrTurb_), smagoPrefactor(smagoPrefactor_), partners(partners_)
312{
313 this->getName() = "SmagorinskyBoussinesqCouplingPostProcessor3D";
314 // we normalize the direction of force vector
315 T normDir = T();
316 for (unsigned iD = 0; iD < dir.size(); ++iD) {
317 normDir += dir[iD]*dir[iD];
318 }
319 normDir = util::sqrt(normDir);
320 for (unsigned iD = 0; iD < dir.size(); ++iD) {
321 dir[iD] /= normDir;
322 }
323
324 for (unsigned iD = 0; iD < dir.size(); ++iD) {
325 forcePrefactor[iD] = gravity * dir[iD];
326 }
327
328 tauTurbADPrefactor = descriptors::invCs2<T,descriptors::D3Q7<descriptors::VELOCITY,descriptors::TAU_EFF>>() / descriptors::invCs2<T,DESCRIPTOR>() / PrTurb;
329 tPartner = static_cast<BlockLattice<T,descriptors::D3Q7<descriptors::VELOCITY,descriptors::TAU_EFF>> *>(partners[0]);
330}
331
332
333template<typename T, typename DESCRIPTOR>
334void SmagorinskyBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>::
335processSubDomain(BlockLattice<T,DESCRIPTOR>& blockLattice,
336 int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
337{
338
339 int newX0, newX1, newY0, newY1, newZ0, newZ1;
340 if ( util::intersect (
341 x0, x1, y0, y1, z0, z1,
342 x0_, x1_, y0_, y1_, z0_, z1_,
343 newX0, newX1, newY0, newY1, newZ0, newZ1 ) ) {
344
345 for (int iX=newX0; iX<=newX1; ++iX) {
346 for (int iY=newY0; iY<=newY1; ++iY) {
347 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
348
349 // computation of the bousinessq force
350 auto force = blockLattice.get(iX,iY,iZ).template getFieldPointer<descriptors::FORCE>();
351 T temperatureDifference = tPartner->get(iX,iY,iZ).computeRho() - T0;
352 for (unsigned iD = 0; iD < L::d; ++iD) {
353 force[iD] = forcePrefactor[iD] * temperatureDifference;
354 }
355
356 // Velocity coupling
357 auto u = tPartner->get(iX,iY,iZ).template getField<descriptors::VELOCITY>();
358
359 // tau coupling
360 auto tauNS = blockLattice.get(iX,iY,iZ).template getFieldPointer<descriptors::TAU_EFF>();
361 auto tauAD = tPartner->get(iX,iY,iZ).template getFieldPointer<descriptors::TAU_EFF>();
362
363 T rho, pi[util::TensorVal<DESCRIPTOR >::n];
364 blockLattice.get(iX,iY,iZ).computeAllMomenta(rho, u.data(), pi);
365 tPartner->get(iX,iY,iZ).template setField<descriptors::VELOCITY>(u);
366 T PiNeqNormSqr = pi[0]*pi[0] + 2.0*pi[1]*pi[1] + pi[2]*pi[2];
367 if (util::TensorVal<DESCRIPTOR >::n == 6) {
368 PiNeqNormSqr += pi[2]*pi[2] + pi[3]*pi[3] + 2*pi[4]*pi[4] +pi[5]*pi[5];
369 }
370 T PiNeqNorm = util::sqrt(PiNeqNormSqr);
372 T tau_mol_NS = 1. / blockLattice.get(iX,iY,iZ).getDynamics()->getParameters(blockLattice).template getOrFallback<descriptors::OMEGA>(0);
373 T tau_mol_AD = 1. / tPartner->get(iX,iY,iZ).getDynamics()->getParameters(*tPartner).template getOrFallback<descriptors::OMEGA>(0);
375 T tau_turb_NS = 0.5*(util::sqrt(tau_mol_NS*tau_mol_NS + smagoPrefactor/rho*PiNeqNorm) - tau_mol_NS);
376 if (tau_turb_NS != tau_turb_NS) {
377 tau_turb_NS = 0.0;
378 }
380 tauNS[0] = tau_mol_NS+tau_turb_NS;
381
382 T tau_turb_AD = tau_turb_NS * tauTurbADPrefactor;
383 tauAD[0] = tau_mol_AD+tau_turb_AD;
384 }
385 }
386 }
387 }
388}
389
390template<typename T, typename DESCRIPTOR>
391void SmagorinskyBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>::
392process(BlockLattice<T,DESCRIPTOR>& blockLattice)
393{
394 processSubDomain(blockLattice, x0, x1, y0, y1, z0, z1);
395}
396
398
399template<typename T, typename DESCRIPTOR>
400SmagorinskyBoussinesqCouplingGenerator3D<T,DESCRIPTOR>::
401SmagorinskyBoussinesqCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
402 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_, T PrTurb_, T smagoPrefactor_)
403 : LatticeCouplingGenerator3D<T,DESCRIPTOR>(x0_, x1_, y0_, y1_, z0_, z1_),
404 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_), dir(dir_), PrTurb(PrTurb_), smagoPrefactor(smagoPrefactor_)
405{ }
406
407template<typename T, typename DESCRIPTOR>
408PostProcessor3D<T,DESCRIPTOR>* SmagorinskyBoussinesqCouplingGenerator3D<T,DESCRIPTOR>::generate (
409 std::vector<BlockStructureD<3>* > partners) const
410{
411 return new SmagorinskyBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>(
412 this->x0,this->x1,this->y0,this->y1,this->z0,this->z1, gravity, T0, deltaTemp, dir, PrTurb, smagoPrefactor, partners);
413}
414
415template<typename T, typename DESCRIPTOR>
420
421
422
423//=====================================================================================
424//============== AdvectionDiffusionParticleCouplingPostProcessor3D ===========
425//=====================================================================================
426
427template<typename T, typename DESCRIPTOR, typename ADLattice, typename FIELD_A, typename FIELD_B>
428AdvectionDiffusionParticleCouplingPostProcessor3D<T,DESCRIPTOR,ADLattice,FIELD_A,FIELD_B>::
429AdvectionDiffusionParticleCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, int iC_,
430 std::vector<BlockStructureD<3>* > partners_,
431 std::vector<std::reference_wrapper<AdvectionDiffusionForce3D<T,DESCRIPTOR,ADLattice>>> forces_)
432 : _forces(forces_),
433 x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_), iC(iC_),
434 _partnerLattice(static_cast<BlockLattice<T,ADLattice> *>(partners_[0])),
435 _cell(_partnerLattice->get(x0,y0,z0)),
436 _cellXp(_partnerLattice->get(x0+1,y0,z0)),
437 _cellXn(_partnerLattice->get(x0-1,y0,z0)),
438 _cellYp(_partnerLattice->get(x0,y0+1,z0)),
439 _cellYn(_partnerLattice->get(x0,y0-1,z0)),
440 _cellZp(_partnerLattice->get(x0,y0,z0+1)),
441 _cellZn(_partnerLattice->get(x0,y0,z0-1))
442{
443 this->getName() = "AdvectionDiffusionParticleCouplingPostProcessor3D";
444}
445
446template<typename T, typename DESCRIPTOR, typename ADLattice, typename FIELD_A, typename FIELD_B>
447void AdvectionDiffusionParticleCouplingPostProcessor3D<T,DESCRIPTOR,ADLattice,FIELD_A,FIELD_B>::
448processSubDomain(BlockLattice<T,DESCRIPTOR>& blockLattice,
449 int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
450{
451 auto vel = par ? _cell.template getField<FIELD_B>() : _cell.template getField<FIELD_A>();
452 //auto vel_new = par ? _cell.template getFieldPointer<FIELD_A>() : _cell.template getFieldPointer<FIELD_B>();
453
454 auto velXp = par ? _cellXp.template getField<FIELD_B>() : _cellXp.template getField<FIELD_A>();
455 auto velXn = par ? _cellXn.template getField<FIELD_B>() : _cellXn.template getField<FIELD_A>();
456 auto velYp = par ? _cellYp.template getField<FIELD_B>() : _cellYp.template getField<FIELD_A>();
457 auto velYn = par ? _cellYn.template getField<FIELD_B>() : _cellYn.template getField<FIELD_A>();
458 auto velZp = par ? _cellZp.template getField<FIELD_B>() : _cellZp.template getField<FIELD_A>();
459 auto velZn = par ? _cellZn.template getField<FIELD_B>() : _cellZn.template getField<FIELD_A>();
460
461 int newX0, newX1, newY0, newY1, newZ0, newZ1;
462
463 if ( util::intersect (
464 x0, x1, y0, y1, z0, z1,
465 x0_, x1_, y0_, y1_, z0_, z1_,
466 newX0, newX1, newY0, newY1, newZ0, newZ1 ) ) {
467
468 for (int iX=newX0; iX<=newX1; ++iX) {
469 for (int iY=newY0; iY<=newY1; ++iY) {
470 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
471 int latticeR[4] = {iC, iX, iY, iZ};
472 T velGrad[3] = {0.,0.,0.};
473 T forceValue[3] = {0.,0.,0.};
474 T velF[3] = {0.,0.,0.};
475
476 auto nsCell = blockLattice.get(iX,iY,iZ);
477
478 if (_forces.begin() != _forces.end()) {
479 // calculating upwind Gradient
480 // vel contains velocity information on ADlattice
481 // velGrad contains upwind vel on ADlattice
482 if (vel[0]<0.) {
483 velGrad[0] = vel[0]*(velXp[0]-vel[0]);
484 velGrad[1] = vel[0]*(velXp[1]-vel[1]);
485 velGrad[2] = vel[0]*(velXp[2]-vel[2]);
486 }
487 else {
488 velGrad[0] = vel[0]*(vel[0]-velXn[0]);
489 velGrad[1] = vel[0]*(vel[1]-velXn[1]);
490 velGrad[2] = vel[0]*(vel[2]-velXn[2]);
491 }
492 if (vel[1]<0.) {
493 velGrad[0] += vel[1]*(velYp[0]-vel[0]);
494 velGrad[1] += vel[1]*(velYp[1]-vel[1]);
495 velGrad[2] += vel[1]*(velYp[2]-vel[2]);
496 }
497 else {
498 velGrad[0] += vel[1]*(vel[0]-velYn[0]);
499 velGrad[1] += vel[1]*(vel[1]-velYn[1]);
500 velGrad[2] += vel[1]*(vel[2]-velYn[2]);
501 }
502 if (vel[2]<0.) {
503 velGrad[0] += vel[2]*(velZp[0]-vel[0]);
504 velGrad[1] += vel[2]*(velZp[1]-vel[1]);
505 velGrad[2] += vel[2]*(velZp[2]-vel[2]);
506 }
507 else {
508 velGrad[0] += vel[2]*(vel[0]-velZn[0]);
509 velGrad[1] += vel[2]*(vel[1]-velZn[1]);
510 velGrad[2] += vel[2]*(vel[2]-velZn[2]);
511 }
512
513 for (AdvectionDiffusionForce3D<T, DESCRIPTOR, ADLattice>& f : _forces) {
514 // writes force in forceValues, vel refers to ADlattice
515 auto adCell = _partnerLattice->get(x0,y0,z0);
516 f.applyForce(forceValue, &nsCell, &adCell, vel.data(), latticeR);
517 if (par) {
518 _cell.template setField<FIELD_B>(vel);
519 }
520 else {
521 _cell.template setField<FIELD_A>(vel);
522 }
523 }
524
525 // compute new particle velocity
526 Vector<T,DESCRIPTOR::d> newVel;
527 for (int i=0; i < DESCRIPTOR::d; i++) {
528 newVel[i] = vel[i] + forceValue[i] - velGrad[i];
529 }
530 if (par) {
531 _cell.template setField<FIELD_A>(newVel);
532 } else {
533 _cell.template setField<FIELD_B>(newVel);
534 }
535 }
536 else { // set particle velocity to fluid velocity
537 nsCell.computeU(velF);
538 Vector<T,DESCRIPTOR::d> newVel;
539 for (int i = 0; i < DESCRIPTOR::d; i++) {
540 newVel[i] = velF[i];
541 }
542 if (par) {
543 _cell.template setField<FIELD_A>(newVel);
544 } else {
545 _cell.template setField<FIELD_B>(newVel);
546 }
547 }
548 }
549 }
550 }
551 }
552 par = !par;
553}
554
555template<typename T, typename DESCRIPTOR, typename ADLattice, typename FIELD_A, typename FIELD_B>
556void AdvectionDiffusionParticleCouplingPostProcessor3D<T,DESCRIPTOR,ADLattice,FIELD_A,FIELD_B>::
557process(BlockLattice<T,DESCRIPTOR>& blockLattice)
558{
559 processSubDomain(blockLattice, x0, x1, y0, y1, z0, z1);
560}
561
562// LatticeCouplingGenerator for one-way advectionDiffusion coupling with Stokes drag
563
564template<typename T, typename DESCRIPTOR, typename ADLattice, typename FIELD_A, typename FIELD_B>
569
570template<typename T, typename DESCRIPTOR, typename ADLattice, typename FIELD_A, typename FIELD_B>
571PostProcessor3D<T,DESCRIPTOR>* AdvectionDiffusionParticleCouplingGenerator3D<T,DESCRIPTOR,ADLattice,FIELD_A,FIELD_B>::generate(
572 std::vector<BlockStructureD<3>* > partners) const
573{
574 return new AdvectionDiffusionParticleCouplingPostProcessor3D<T,DESCRIPTOR,ADLattice,FIELD_A,FIELD_B>(
575 this->x0,this->x1,this->y0,this->y1,this->z0,this->z1, this->iC, partners, ADforces);
576}
577
578template<typename T, typename DESCRIPTOR, typename ADLattice, typename FIELD_A, typename FIELD_B>
583
584template<typename T, typename DESCRIPTOR, typename ADLattice, typename FIELD_A, typename FIELD_B>
590
591
592//=====================================================================================
593//============== PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D ===========
594//=====================================================================================
595
596template<typename T, typename DESCRIPTOR>
597PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D<T,DESCRIPTOR>::
598PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
599 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_,
600 std::vector<BlockStructureD<3>* > partners_)
601 : x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_),
602 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_),
603 dir(dir_), partners(partners_)
604{
605 this->getName() = "PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D";
606 // we normalize the direction of force vector
607 T normDir = T();
608 for (unsigned iD = 0; iD < dir.size(); ++iD) {
609 normDir += dir[iD]*dir[iD];
610 }
611 normDir = util::sqrt(normDir);
612 for (unsigned iD = 0; iD < dir.size(); ++iD) {
613 dir[iD] /= normDir;
614 }
615}
616
617template<typename T, typename DESCRIPTOR>
618void PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D<T,DESCRIPTOR>::
619processSubDomain(BlockLattice<T,DESCRIPTOR>& blockLattice,
620 int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
621{
622 typedef DESCRIPTOR L;
623 enum {x,y,z};
624
625 BlockLattice<T,descriptors::PorousAdvectionDiffusionD3Q7Descriptor> *tPartner =
626 static_cast<BlockLattice<T,descriptors::PorousAdvectionDiffusionD3Q7Descriptor> *>(partners[0]);
627
628 int newX0, newX1, newY0, newY1, newZ0, newZ1;
629 if ( util::intersect (
630 x0, x1, y0, y1, z0, z1,
631 x0_, x1_, y0_, y1_, z0_, z1_,
632 newX0, newX1, newY0, newY1, newZ0, newZ1 ) ) {
633
634 for (int iX=newX0; iX<=newX1; ++iX) {
635 for (int iY=newY0; iY<=newY1; ++iY) {
636 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
637 // Velocity coupling
638 auto u = tPartner->get(iX,iY,iZ).template getField<descriptors::VELOCITY>();
639 blockLattice.get(iX,iY,iZ).computeU(u.data());
640 tPartner->get(iX,iY,iZ).template setField<descriptors::VELOCITY>(u);
641
642 //coupling between the temperature and navier stokes.
643
644 auto force = blockLattice.get(iX,iY,iZ).template getFieldPointer<descriptors::FORCE>();
645 // this should return the interpolated solid-fluid temperature
646 T temperature = tPartner->get(iX,iY,iZ).computeRho();
647 T rho = blockLattice.get(iX,iY,iZ).computeRho();
648 for (unsigned iD = 0; iD < L::d; ++iD) {
649 force[iD] = gravity * rho * (temperature - T0) / deltaTemp * dir[iD];
650 }
651 }
652 }
653 }
654 }
655}
656
657template<typename T, typename DESCRIPTOR>
658void PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D<T,DESCRIPTOR>::
659process(BlockLattice<T,DESCRIPTOR>& blockLattice)
660{
661 processSubDomain(blockLattice, x0, x1, y0, y1, z0, z1);
662}
663
664template<typename T, typename DESCRIPTOR>
665PorousNavierStokesAdvectionDiffusionCouplingGenerator3D<T,DESCRIPTOR>::
666PorousNavierStokesAdvectionDiffusionCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_,int z0_, int z1_,
667 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_)
668 : LatticeCouplingGenerator3D<T,DESCRIPTOR>(x0_, x1_, y0_, y1_, z0_, z1_),
669 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_), dir(dir_)
670{ }
671
672template<typename T, typename DESCRIPTOR>
673PostProcessor3D<T,DESCRIPTOR>* PorousNavierStokesAdvectionDiffusionCouplingGenerator3D<T,DESCRIPTOR>::generate (
674 std::vector<BlockStructureD<3>* > partners) const
675{
676 return new PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D<T,DESCRIPTOR>(
677 this->x0,this->x1,this->y0,this->y1,this->z0,this->z1, gravity, T0, deltaTemp, dir,partners);
678}
679
680template<typename T, typename DESCRIPTOR>
685
686
687//=====================================================================================
688//============== MixedScaleBoussinesqCouplingPostProcessor3D ===============
689//=====================================================================================
690
691template<typename T, typename DESCRIPTOR>
692MixedScaleBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>::
693MixedScaleBoussinesqCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
694 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_, T PrTurb_,
695 std::vector<BlockStructureD<3>* > partners_)
696 : x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_),
697 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_),
698 dir(dir_), PrTurb(PrTurb_), partners(partners_)
699{
700 // we normalize the direction of force vector
701 T normDir = T();
702 for (unsigned iD = 0; iD < dir.size(); ++iD) {
703 normDir += dir[iD]*dir[iD];
704 }
705 normDir = util::sqrt(normDir);
706 for (unsigned iD = 0; iD < dir.size(); ++iD) {
707 dir[iD] /= normDir;
708 }
709
710 for (unsigned iD = 0; iD < dir.size(); ++iD) {
711 forcePrefactor[iD] = gravity * dir[iD];
712 }
713
714 tauTurbADPrefactor = descriptors::invCs2<T,descriptors::D3Q7<>>() / descriptors::invCs2<T,DESCRIPTOR>() / PrTurb;
715 tPartner = static_cast<BlockLattice<T,descriptors::D3Q7<descriptors::VELOCITY,descriptors::TAU_EFF,descriptors::CUTOFF_HEAT_FLUX>> *>(partners[0]);
716}
717
718template<typename T, typename DESCRIPTOR>
719void MixedScaleBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>::
720processSubDomain(BlockLattice<T,DESCRIPTOR>& blockLattice,
721 int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
722{
723
724 const T C_nu = 0.04;
725 const T C_alpha = 0.5;
726 // const T deltaT = 1.0;
727
728 const T invCs2 = descriptors::invCs2<T,DESCRIPTOR>();
729 const T invCs2_g = descriptors::invCs2<T,descriptors::D3Q7<>>();
730
731 int newX0, newX1, newY0, newY1, newZ1, newZ0;
732 if ( util::intersect (
733 x0, x1, y0, y1, z0, z1,
734 x0_, x1_, y0_, y1_, z0_, z1_,
735 newX0, newX1, newY0, newY1, newZ1, newZ0 ) ) {
736
737 auto& heatFluxCache = blockLattice.template getField<HEAT_FLUX_CACHE>()[0];
738 for (int iX=newX0-1; iX<=newX1+1; ++iX) {
739 for (int iY=newY0-1; iY<=newY1+1; ++iY) {
740 for (int iZ=newZ0-1; iZ<=newZ1+1; ++iZ) {
741 const T temperature = tPartner->get(iX,iY,iZ).computeRho();
742 heatFluxCache[blockLattice.getCellId(iX, iY, iZ)] = temperature;
743
744 // computation of the bousinessq force
745 T temperatureDifference = temperature - T0;
746 blockLattice.get(iX,iY,iZ).template setField<descriptors::FORCE>(temperatureDifference*forcePrefactor);
747
748 FieldD<T,DESCRIPTOR,descriptors::VELOCITY> u;
749 blockLattice.get(iX,iY,iZ).computeU(u.data());
750 tPartner->get(iX,iY,iZ).template setField<descriptors::VELOCITY>(u);
751 }
752 }
753 }
754
755 for (int iX=newX0; iX<=newX1; ++iX) {
756 for (int iY=newY0; iY<=newY1; ++iY) {
757 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
758
759 auto u_ppp = tPartner->get(iX+1, iY+1, iZ+1).template getField<descriptors::VELOCITY>();
760 auto u_0pp = tPartner->get(iX, iY+1, iZ+1).template getField<descriptors::VELOCITY>();
761 auto u_npp = tPartner->get(iX-1, iY+1, iZ+1).template getField<descriptors::VELOCITY>();
762 auto u_p0p = tPartner->get(iX+1, iY, iZ+1).template getField<descriptors::VELOCITY>();
763 auto u_00p = tPartner->get(iX, iY, iZ+1).template getField<descriptors::VELOCITY>();
764 auto u_n0p = tPartner->get(iX-1, iY, iZ+1).template getField<descriptors::VELOCITY>();
765 auto u_pnp = tPartner->get(iX+1, iY-1, iZ+1).template getField<descriptors::VELOCITY>();
766 auto u_0np = tPartner->get(iX, iY-1, iZ+1).template getField<descriptors::VELOCITY>();
767 auto u_nnp = tPartner->get(iX-1, iY-1, iZ+1).template getField<descriptors::VELOCITY>();
768
769 auto u_pp0 = tPartner->get(iX+1, iY+1, iZ ).template getField<descriptors::VELOCITY>();
770 auto u_0p0 = tPartner->get(iX, iY+1, iZ ).template getField<descriptors::VELOCITY>();
771 auto u_np0 = tPartner->get(iX-1, iY+1, iZ ).template getField<descriptors::VELOCITY>();
772 auto u_p00 = tPartner->get(iX+1, iY, iZ ).template getField<descriptors::VELOCITY>();
773 auto u_000 = tPartner->get(iX, iY, iZ ).template getField<descriptors::VELOCITY>();
774 auto u_n00 = tPartner->get(iX-1, iY, iZ ).template getField<descriptors::VELOCITY>();
775 auto u_pn0 = tPartner->get(iX+1, iY-1, iZ ).template getField<descriptors::VELOCITY>();
776 auto u_0n0 = tPartner->get(iX, iY-1, iZ ).template getField<descriptors::VELOCITY>();
777 auto u_nn0 = tPartner->get(iX-1, iY-1, iZ ).template getField<descriptors::VELOCITY>();
778
779 auto u_ppn = tPartner->get(iX+1, iY+1, iZ-1).template getField<descriptors::VELOCITY>();
780 auto u_0pn = tPartner->get(iX, iY+1, iZ-1).template getField<descriptors::VELOCITY>();
781 auto u_npn = tPartner->get(iX-1, iY+1, iZ-1).template getField<descriptors::VELOCITY>();
782 auto u_p0n = tPartner->get(iX+1, iY, iZ-1).template getField<descriptors::VELOCITY>();
783 auto u_00n = tPartner->get(iX, iY, iZ-1).template getField<descriptors::VELOCITY>();
784 auto u_n0n = tPartner->get(iX-1, iY, iZ-1).template getField<descriptors::VELOCITY>();
785 auto u_pnn = tPartner->get(iX+1, iY-1, iZ-1).template getField<descriptors::VELOCITY>();
786 auto u_0nn = tPartner->get(iX, iY-1, iZ-1).template getField<descriptors::VELOCITY>();
787 auto u_nnn = tPartner->get(iX-1, iY-1, iZ-1).template getField<descriptors::VELOCITY>();
788
789 const T *h_ppp = & heatFluxCache[blockLattice.getCellId(iX+1, iY+1, iZ+1)];
790 const T *h_0pp = & heatFluxCache[blockLattice.getCellId(iX, iY+1, iZ+1)];
791 const T *h_npp = & heatFluxCache[blockLattice.getCellId(iX-1, iY+1, iZ+1)];
792 const T *h_p0p = & heatFluxCache[blockLattice.getCellId(iX+1, iY, iZ+1)];
793 const T *h_pnp = & heatFluxCache[blockLattice.getCellId(iX+1, iY-1, iZ+1)];
794 const T *h_00p = & heatFluxCache[blockLattice.getCellId(iX, iY, iZ+1)];
795 const T *h_0np = & heatFluxCache[blockLattice.getCellId(iX, iY-1, iZ+1)];
796 const T *h_n0p = & heatFluxCache[blockLattice.getCellId(iX-1, iY, iZ+1)];
797 const T *h_nnp = & heatFluxCache[blockLattice.getCellId(iX-1, iY-1, iZ+1)];
798
799 const T *h_pp0 = & heatFluxCache[blockLattice.getCellId(iX+1, iY+1, iZ )];
800 const T *h_0p0 = & heatFluxCache[blockLattice.getCellId(iX, iY+1, iZ )];
801 const T *h_np0 = & heatFluxCache[blockLattice.getCellId(iX-1, iY+1, iZ )];
802 const T *h_p00 = & heatFluxCache[blockLattice.getCellId(iX+1, iY, iZ )];
803 const T *h_pn0 = & heatFluxCache[blockLattice.getCellId(iX+1, iY-1, iZ )];
804 const T *h_000 = & heatFluxCache[blockLattice.getCellId(iX, iY, iZ )];
805 const T *h_0n0 = & heatFluxCache[blockLattice.getCellId(iX, iY-1, iZ )];
806 const T *h_n00 = & heatFluxCache[blockLattice.getCellId(iX-1, iY, iZ )];
807 const T *h_nn0 = & heatFluxCache[blockLattice.getCellId(iX-1, iY-1, iZ )];
808
809 const T *h_ppn = & heatFluxCache[blockLattice.getCellId(iX+1, iY+1, iZ-1)];
810 const T *h_0pn = & heatFluxCache[blockLattice.getCellId(iX, iY+1, iZ-1)];
811 const T *h_npn = & heatFluxCache[blockLattice.getCellId(iX-1, iY+1, iZ-1)];
812 const T *h_p0n = & heatFluxCache[blockLattice.getCellId(iX+1, iY, iZ-1)];
813 const T *h_pnn = & heatFluxCache[blockLattice.getCellId(iX+1, iY-1, iZ-1)];
814 const T *h_00n = & heatFluxCache[blockLattice.getCellId(iX, iY, iZ-1)];
815 const T *h_0nn = & heatFluxCache[blockLattice.getCellId(iX, iY-1, iZ-1)];
816 const T *h_n0n = & heatFluxCache[blockLattice.getCellId(iX-1, iY, iZ-1)];
817 const T *h_nnn = & heatFluxCache[blockLattice.getCellId(iX-1, iY-1, iZ-1)];
818
819// cout<<"h_ppn= "<<h_ppp[0]<<endl;
820// cout<<"h_0pn= "<<h_0pp[0]<<endl;
821// cout<<"h_npn= "<<h_npp[0]<<endl;
822// cout<<"h_p0n= "<<h_p0p[0]<<endl;
823// cout<<"h_00n= "<<h_00p[0]<<endl;
824// cout<<"h_n0n= "<<h_n0p[0]<<endl;
825// cout<<"h_pnn= "<<h_pnp[0]<<endl;
826// cout<<"h_0nn= "<<h_0np[0]<<endl;
827// cout<<"h_nnn= "<<h_nnp[0]<<endl;
828
829 //Testfilter h 3D
830 Vector<T,3> filtered_u;
831 T filtered_h = (
832 (
833 ( h_ppp[0] + 2. * h_0pp[0] + h_npp[0]) + 2.*
834 ( h_p0p[0] + 2. * h_00p[0] + h_n0p[0]) +
835 ( h_pnp[0] + 2. * h_0np[0] + h_nnp[0])
836 ) * 0.25 * 0.25 + 2. *
837
838 (
839 ( h_pp0[0] + 2. * h_0p0[0] + h_np0[0]) + 2.*
840 ( h_p00[0] + 2. * h_000[0] + h_n00[0]) +
841 ( h_pn0[0] + 2. * h_0n0[0] + h_nn0[0])
842 ) * 0.25 * 0.25 +
843
844 (
845 ( h_ppn[0] + 2. * h_0pn[0] + h_npn[0]) + 2.*
846 ( h_p0n[0] + 2. * h_00n[0] + h_n0n[0]) +
847 ( h_pnn[0] + 2. * h_0nn[0] + h_nnn[0])
848 ) * 0.25 * 0.25
849 ) * 0.25;
850
851 //Testfilter u 3D
852 filtered_u = (
853 (
854 ( u_ppp + 2. * u_0pp + u_npp) + 2.*
855 ( u_p0p + 2. * u_00p + u_n0p) +
856 ( u_pnp + 2. * u_0np + u_nnp)
857 ) * 0.25 * 0.25 + 2. *
858
859 (
860 ( u_pp0 + 2. * u_0p0 + u_np0) + 2.*
861 ( u_p00 + 2. * u_000 + u_n00) +
862 ( u_pn0 + 2. * u_0n0 + u_nn0)
863 ) * 0.25 * 0.25 +
864
865 (
866 ( u_ppn + 2. * u_0pn + u_npn) + 2.*
867 ( u_p0n + 2. * u_00n + u_n0n) +
868 ( u_pnn + 2. * u_0nn + u_nnn)
869 ) * 0.25 * 0.25
870 ) * 0.25;
871
872 //Ende Testfilter u
873 // filtered_u[i] -= u_00[i];
874
875// cout << "u["<<i<<"]= "<< filtered_u[i] << std::endl;
876// cout << "h[0]= "<< filtered_h[0] << std::endl;
877 const Vector<T,3> filtered_u_reduced = u_000 - filtered_u;
878 const Vector<T,3> filtered_heatFlux_reduced = u_000*h_000[0] - filtered_u*filtered_h;
879 T cutoffKinEnergy = filtered_u_reduced[0] * filtered_u_reduced[0]
880 + filtered_u_reduced[1] * filtered_u_reduced[1]
881 + filtered_u_reduced[2] * filtered_u_reduced[2];
882 T cutoffHeatFlux = filtered_heatFlux_reduced[0] * filtered_heatFlux_reduced[0]
883 + filtered_heatFlux_reduced[1] * filtered_heatFlux_reduced[1]
884 + filtered_heatFlux_reduced[2] * filtered_heatFlux_reduced[2];
885// cout << "filtered_u_reduced= " << filtered_u_reduced << std::endl;
886// cout << "filtered_heatFlux_reduced= " << filtered_heatFlux_reduced << std::endl;
887// cout << "cutoffKinEnergy= " << cutoffKinEnergy << std::endl;
888// cout << "cutoffHeatFlux= " << cutoffHeatFlux << std::endl;
889
890 blockLattice.get(iX,iY,iZ).template setField<descriptors::CUTOFF_KIN_ENERGY>(util::pow(0.5*cutoffKinEnergy, 0.25));
891 tPartner->get(iX,iY,iZ).template setField<descriptors::CUTOFF_HEAT_FLUX>(util::pow(0.5*cutoffHeatFlux, 0.25));
892
893// std::cout<<"cutoffKinEnergy_14= "<<cutoffKinEnergy_14[0]<<std::endl;
894// std::cout<<"cutoffHeatFlux_14= "<<cutoffHeatFlux_14[0]<<std::endl;
895
896
897 // tau coupling
898 auto tauNS = blockLattice.get(iX,iY,iZ).template getField<descriptors::TAU_EFF>();
899 auto tauAD = tPartner->get(iX,iY,iZ).template getField<descriptors::TAU_EFF>();
900
901
902// std::cout<<"tau_NS= "<<tauNS[0]<<std::endl;
903// std::cout<<"tau_AD= "<<tauAD[0]<<std::endl;
904
906 T tau_mol_NS = 1. / blockLattice.get(iX, iY, iZ).getDynamics()->getParameters(blockLattice).template getOrFallback<descriptors::OMEGA>(0);
907 T tau_mol_AD = 1. / tPartner->get(iX, iY, iZ).getDynamics()->getParameters(*tPartner).template getOrFallback<descriptors::OMEGA>(0);
908
909 const T temperature = h_000[0];
910
911 // computation of the bousinessq force
912 auto force = blockLattice.get(iX,iY,iZ).template getField<descriptors::FORCE>();
913 // T temperatureDifference = temperature - T0;
914 // for (unsigned iD = 0; iD < L::d; ++iD) {
915 // force[iD] = forcePrefactor[iD] * temperatureDifference;
916 // }
917
918 auto u = tPartner->get(iX,iY,iZ).template getFieldPointer<descriptors::VELOCITY>();
919 T rho, pi[util::TensorVal<DESCRIPTOR>::n], j[DESCRIPTOR::d];
920 // blockLattice.get(iX,iY).computeAllMomenta(rho, u, pi);
921 rho = blockLattice.get(iX,iY,iZ).computeRho();
922 blockLattice.get(iX,iY,iZ).computeStress(pi);
923
924 int iPi = 0;
925 const int dim = DESCRIPTOR::d;
926 for (int Alpha=0; Alpha<dim; ++Alpha) {
927 for (int Beta=Alpha; Beta<dim; ++Beta) {
928 pi[iPi] += rho/2.*(force[Alpha]*u[Beta] + u[Alpha]*force[Beta]);
929 ++iPi;
930 }
931 }
932 const T piSqr[6] = {pi[0]*pi[0], pi[1]*pi[1], pi[2]*pi[2], pi[3]*pi[3], pi[4]*pi[4], pi[5]*pi[5]};
933 const T PiNeqNormSqr = piSqr[0] + 2.0 * piSqr[1] + 2.0 * piSqr[2] + piSqr[3] + 2.0 * piSqr[4] + piSqr[5];
934 const T PiNeqNorm = util::sqrt(PiNeqNormSqr);
935
936 tPartner->get(iX, iY, iZ).computeJ(j);
937
938 const Vector<T,3> jNeq = {(j[0] - temperature * u[0]), (j[1] - temperature * u[1]), (j[2] - temperature * u[2])};
939
940 const Vector<T,6> t = { 2.0 * jNeq[0] * pi[0], (jNeq[0] + jNeq[1]) * pi[1], (jNeq[0] + jNeq[2]) * pi[2], 2.0 * jNeq[1] * pi[3], (jNeq[1] + jNeq[2]) * pi[4], 2.0 * jNeq[2] * pi[5]};
941 const Vector<T,6> tSqr = {t[0]*t[0], t[1]*t[1], t[2]*t[2], t[3]*t[3], t[4]*t[4], t[5]*t[5]};
942 const T TnormSqr = tSqr[0] + 2.0 * tSqr[1] + 2.0 * tSqr[2] + tSqr[3] + 2.0 * tSqr[4] + tSqr[5];
943 const T Tnorm = util::sqrt(2.0 * 0.25 * TnormSqr);
944
945 auto cutoffKinEnergy_14 = blockLattice.get(iX,iY,iZ).template getField<descriptors::CUTOFF_KIN_ENERGY>();
946 auto cutoffHeatFlux_14 = tPartner->get(iX,iY,iZ).template getField<descriptors::CUTOFF_HEAT_FLUX>();
947
948 const T tau_turb_NS = C_nu * util::sqrt(util::sqrt(2.)/2.) * invCs2 * invCs2 * util::sqrt(PiNeqNorm / rho / tauNS) * cutoffKinEnergy_14;
949 const T tau_turb_AD = C_alpha * invCs2 / rho * util::sqrt(Tnorm * invCs2_g / tauNS / tauAD) * cutoffHeatFlux_14;
950
952 blockLattice.get(iX,iY,iZ).template setField<descriptors::TAU_EFF>(tau_mol_NS + tau_turb_NS);
953 tPartner->get(iX,iY,iZ).template setField<descriptors::TAU_EFF>(tau_mol_AD + tau_turb_AD);
954
955 }
956 }
957 }
958 }
959}
960
961template<typename T, typename DESCRIPTOR>
962void MixedScaleBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>::
963process(BlockLattice<T,DESCRIPTOR>& blockLattice)
964{
965 processSubDomain(blockLattice, x0, x1, y0, y1, z0, z1);
966}
967
969
970template<typename T, typename DESCRIPTOR>
971MixedScaleBoussinesqCouplingGenerator3D<T,DESCRIPTOR>::
972MixedScaleBoussinesqCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
973 T gravity_, T T0_, T deltaTemp_, std::vector<T> dir_, T PrTurb_)
974 : LatticeCouplingGenerator3D<T,DESCRIPTOR>(x0_, x1_, y0_, y1_, z0_, z1_),
975 gravity(gravity_), T0(T0_), deltaTemp(deltaTemp_), dir(dir_), PrTurb(PrTurb_)
976{ }
977
978template<typename T, typename DESCRIPTOR>
979PostProcessor3D<T,DESCRIPTOR>* MixedScaleBoussinesqCouplingGenerator3D<T,DESCRIPTOR>::generate (
980 std::vector<BlockStructureD<3>* > partners) const
981{
982 return new MixedScaleBoussinesqCouplingPostProcessor3D<T,DESCRIPTOR>(
983 this->x0,this->x1,this->y0,this->y1, this-> z0, this->z1, gravity, T0, deltaTemp, dir, PrTurb, partners);
984}
985
986template<typename T, typename DESCRIPTOR>
991
992
993//=====================================================================================
994//============== VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D ===========
995//=====================================================================================
996
997template<typename T, typename DESCRIPTOR, typename POROSITY, typename ADLattice, typename FIELD_A, typename FIELD_B>
998VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D<T,DESCRIPTOR,POROSITY,ADLattice,FIELD_A,FIELD_B>::
999VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, int iC_,
1000 std::vector<BlockStructureD<3>* > partners_,
1001 std::vector<std::reference_wrapper<AdvectionDiffusionForce3D<T,DESCRIPTOR,ADLattice>>> forces_)
1002 : _forces(forces_),
1003 x0(x0_), x1(x1_), y0(y0_), y1(y1_), z0(z0_), z1(z1_), iC(iC_),
1004 _partnerLattice(static_cast<BlockLattice<T,ADLattice> *>(partners_[0])),
1005 _cell(_partnerLattice->get(x0,y0,z0)),
1006 _cellXp(_partnerLattice->get(x0+1,y0,z0)),
1007 _cellXn(_partnerLattice->get(x0-1,y0,z0)),
1008 _cellYp(_partnerLattice->get(x0,y0+1,z0)),
1009 _cellYn(_partnerLattice->get(x0,y0-1,z0)),
1010 _cellZp(_partnerLattice->get(x0,y0,z0+1)),
1011 _cellZn(_partnerLattice->get(x0,y0,z0-1))
1012{
1013 this->getName() = "VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D";
1014}
1015
1016template<typename T, typename DESCRIPTOR, typename POROSITY, typename ADLattice, typename FIELD_A, typename FIELD_B>
1017void VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D<T,DESCRIPTOR,POROSITY,ADLattice,FIELD_A,FIELD_B>::
1018processSubDomain(BlockLattice<T,DESCRIPTOR>& blockLattice,
1019 int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
1020{
1021 auto vel = par ? _cell.template getField<FIELD_B>() : _cell.template getField<FIELD_A>();
1022 auto vel_new = par ? _cell.template getFieldPointer<FIELD_A>() : _cell.template getFieldPointer<FIELD_B>();
1023
1024 auto velXp = par ? _cellXp.template getFieldPointer<FIELD_B>() : _cellXp.template getFieldPointer<FIELD_A>();
1025 auto velXn = par ? _cellXn.template getFieldPointer<FIELD_B>() : _cellXn.template getFieldPointer<FIELD_A>();
1026 auto velYp = par ? _cellYp.template getFieldPointer<FIELD_B>() : _cellYp.template getFieldPointer<FIELD_A>();
1027 auto velYn = par ? _cellYn.template getFieldPointer<FIELD_B>() : _cellYn.template getFieldPointer<FIELD_A>();
1028 auto velZp = par ? _cellZp.template getFieldPointer<FIELD_B>() : _cellZp.template getFieldPointer<FIELD_A>();
1029 auto velZn = par ? _cellZn.template getFieldPointer<FIELD_B>() : _cellZn.template getFieldPointer<FIELD_A>();
1030
1031 auto forceNS = blockLattice.get(x0,y0,z0).template getFieldPointer<descriptors::FORCE>();
1032
1033 auto _cellNSXp = blockLattice.get(x0+1,y0,z0);
1034 auto _cellNSXn = blockLattice.get(x0-1,y0,z0);
1035 auto _cellNSYp = blockLattice.get(x0,y0+1,z0);
1036 auto _cellNSYn = blockLattice.get(x0,y0-1,z0);
1037 auto _cellNSZp = blockLattice.get(x0,y0,z0+1);
1038 auto _cellNSZn = blockLattice.get(x0,y0,z0-1);
1039
1040 int newX0, newX1, newY0, newY1, newZ0, newZ1;
1041 T porosityForce[3];
1042 T pressure;
1043
1044 if ( util::intersect (
1045 x0, x1, y0, y1, z0, z1,
1046 x0_, x1_, y0_, y1_, z0_, z1_,
1047 newX0, newX1, newY0, newY1, newZ0, newZ1 ) ) {
1048
1049 for (int iX=newX0; iX<=newX1; ++iX) {
1050 for (int iY=newY0; iY<=newY1; ++iY) {
1051 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
1052 auto nsCell = blockLattice.get(iX,iY,iZ);
1053 T porosityNS = 1. - _cell.computeRho();
1054 nsCell.template setField<POROSITY>(porosityNS);
1055 }
1056 }
1057 }
1058
1059 for (int iX=newX0; iX<=newX1; ++iX) {
1060 for (int iY=newY0; iY<=newY1; ++iY) {
1061 for (int iZ=newZ0; iZ<=newZ1; ++iZ) {
1062 int latticeR[4] = {iC, iX, iY, iZ};
1063 T velGrad[3] = {0.,0.,0.};
1064 T forceValue[3] = {0.,0.,0.};
1065 T velF[3] = {0.,0.,0.};
1066
1067 auto nsCell = blockLattice.get(iX,iY,iZ);
1068 pressure = nsCell.computeRho() / descriptors::invCs2<T,DESCRIPTOR>();
1069
1070 if (_forces.begin() != _forces.end()) {
1071 // calculating upwind Gradient
1072 // vel contains velocity information on ADlattice
1073 // velGrad contains upwind vel on ADlattice
1074 if (vel[0]<0.) {
1075 velGrad[0] = vel[0]*(velXp[0]-vel[0]);
1076 velGrad[1] = vel[0]*(velXp[1]-vel[1]);
1077 velGrad[2] = vel[0]*(velXp[2]-vel[2]);
1078 }
1079 else {
1080 velGrad[0] = vel[0]*(vel[0]-velXn[0]);
1081 velGrad[1] = vel[0]*(vel[1]-velXn[1]);
1082 velGrad[2] = vel[0]*(vel[2]-velXn[2]);
1083 }
1084 if (vel[1]<0.) {
1085 velGrad[0] += vel[1]*(velYp[0]-vel[0]);
1086 velGrad[1] += vel[1]*(velYp[1]-vel[1]);
1087 velGrad[2] += vel[1]*(velYp[2]-vel[2]);
1088 }
1089 else {
1090 velGrad[0] += vel[1]*(vel[0]-velYn[0]);
1091 velGrad[1] += vel[1]*(vel[1]-velYn[1]);
1092 velGrad[2] += vel[1]*(vel[2]-velYn[2]);
1093 }
1094 if (vel[2]<0.) {
1095 velGrad[0] += vel[2]*(velZp[0]-vel[0]);
1096 velGrad[1] += vel[2]*(velZp[1]-vel[1]);
1097 velGrad[2] += vel[2]*(velZp[2]-vel[2]);
1098 }
1099 else {
1100 velGrad[0] += vel[2]*(vel[0]-velZn[0]);
1101 velGrad[1] += vel[2]*(vel[1]-velZn[1]);
1102 velGrad[2] += vel[2]*(vel[2]-velZn[2]);
1103 }
1104
1105 for (AdvectionDiffusionForce3D<T, DESCRIPTOR, ADLattice>& f : _forces) {
1106 // writes force in forceValues, vel refers to ADlattice
1107 auto adCell = _partnerLattice->get(x0,y0,z0);
1108 f.applyForce(forceValue, &nsCell, &adCell, vel.data(), latticeR);
1109 if (par) {
1110 _cell.template setField<FIELD_B>(vel);
1111 } else {
1112 _cell.template setField<FIELD_A>(vel);
1113 }
1114 }
1115
1116 // compute new particle velocity and opposite fluid force
1117 for (int i=0; i < DESCRIPTOR::d; i++) {
1118 vel_new[i] = vel[i] + forceValue[i] - velGrad[i];
1119 forceNS[i] = -forceValue[i];
1120 }
1121 }
1122 else {
1123
1124 T porXp = _cellNSXp.template getField<POROSITY>();
1125 T porXn = _cellNSXn.template getField<POROSITY>();
1126 T porYp = _cellNSYp.template getField<POROSITY>();
1127 T porYn = _cellNSYn.template getField<POROSITY>();
1128 T porZp = _cellNSZp.template getField<POROSITY>();
1129 T porZn = _cellNSZn.template getField<POROSITY>();
1130
1131 porosityForce[0] = 0.5 * pressure * (porXp - porXn);
1132 porosityForce[1] = 0.5 * pressure * (porYp - porYn);
1133 porosityForce[2] = 0.5 * pressure * (porZp - porZn);
1134
1135 for (int i = 0; i < DESCRIPTOR::d; i++) {
1136 forceNS[i] += porosityForce[i];
1137 }
1138
1139 nsCell.computeU(velF);
1140
1141 for (int i = 0; i < DESCRIPTOR::d; i++) { // set particle velocity to fluid velocity
1142 vel_new[i] = velF[i];
1143 }
1144 }
1145 }
1146 }
1147 }
1148 }
1149 par = !par;
1150}
1151
1152template<typename T, typename DESCRIPTOR, typename POROSITY, typename ADLattice, typename FIELD_A, typename FIELD_B>
1158
1159// LatticeCouplingGenerator for one-way VANS advectionDiffusion coupling with Stokes drag
1160
1161template<typename T, typename DESCRIPTOR, typename POROSITY, typename ADLattice, typename FIELD_A, typename FIELD_B>
1166
1167template<typename T, typename DESCRIPTOR, typename POROSITY, typename ADLattice, typename FIELD_A, typename FIELD_B>
1174
1175template<typename T, typename DESCRIPTOR, typename POROSITY, typename ADLattice, typename FIELD_A, typename FIELD_B>
1180
1181template<typename T, typename DESCRIPTOR, typename POROSITY, typename ADLattice, typename FIELD_A, typename FIELD_B>
1187
1188} // namespace olb
1189
1190#endif
olb::AdvectionDiffusionParticleCouplingGenerator3D::clone
LatticeCouplingGenerator3D< T, DESCRIPTOR > * clone() const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:579
olb::AdvectionDiffusionParticleCouplingGenerator3D::addForce
void addForce(AdvectionDiffusionForce3D< T, DESCRIPTOR, ADLattice > &force)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:585
olb::AdvectionDiffusionParticleCouplingGenerator3D::generate
PostProcessor3D< T, DESCRIPTOR > * generate(std::vector< BlockStructureD< 3 > * > partners) const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:571
olb::AdvectionDiffusionParticleCouplingPostProcessor3D::processSubDomain
void processSubDomain(BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
Execute post-processing step on a sublattice.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:448
olb::AdvectionDiffusionParticleCouplingPostProcessor3D::AdvectionDiffusionParticleCouplingPostProcessor3D
AdvectionDiffusionParticleCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, int iC_, std::vector< BlockStructureD< 3 > * > partners_, std::vector< std::reference_wrapper< AdvectionDiffusionForce3D< T, DESCRIPTOR, ADLattice > > > forces_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:429
olb::AdvectionDiffusionParticleCouplingPostProcessor3D::process
void process(BlockLattice< T, DESCRIPTOR > &blockLattice) override
Execute post-processing step.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:557
olb::BlockLattice
Platform-abstracted block lattice for external access and inter-block interaction.
Definition blockLattice.h:87
olb::BlockLattice::get
Cell< T, DESCRIPTOR > get(CellID iCell)
Get Cell interface for index iCell.
Definition blockLattice.h:171
olb::BlockStructureD
Base of a regular block.
Definition blockStructure.h:54
olb::BlockStructureD::getCellId
CellID getCellId(LatticeR< D > latticeR) const
Get 1D cell ID.
Definition blockStructure.h:125
olb::MixedScaleBoussinesqCouplingGenerator3D::MixedScaleBoussinesqCouplingGenerator3D
MixedScaleBoussinesqCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_, T PrTurb_)
LatticeCouplingGenerator for advectionDiffusion coupling.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:972
olb::MixedScaleBoussinesqCouplingGenerator3D::clone
LatticeCouplingGenerator3D< T, DESCRIPTOR > * clone() const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:987
olb::MixedScaleBoussinesqCouplingGenerator3D::generate
PostProcessor3D< T, DESCRIPTOR > * generate(std::vector< BlockStructureD< 3 > * > partners) const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:979
olb::MixedScaleBoussinesqCouplingPostProcessor3D::MixedScaleBoussinesqCouplingPostProcessor3D
MixedScaleBoussinesqCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_, T PrTurb_, std::vector< BlockStructureD< 3 > * > partners_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:693
olb::MixedScaleBoussinesqCouplingPostProcessor3D::processSubDomain
void processSubDomain(BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
Execute post-processing step on a sublattice.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:720
olb::MixedScaleBoussinesqCouplingPostProcessor3D::process
void process(BlockLattice< T, DESCRIPTOR > &blockLattice) override
Execute post-processing step.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:963
olb::PhaseFieldCouplingGenerator3D::PhaseFieldCouplingGenerator3D
PhaseFieldCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T rho_L, T rho_H, T mu_L, T mu_H, T surface_tension, T interface_thickness)
LatticeCouplingGenerator for advectionDiffusion coupling.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:279
olb::PhaseFieldCouplingGenerator3D::clone
LatticeCouplingGenerator3D< T, DESCRIPTOR > * clone() const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:294
olb::PhaseFieldCouplingGenerator3D::generate
PostProcessor3D< T, DESCRIPTOR > * generate(std::vector< BlockStructureD< 3 > * > partners) const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:286
olb::PhaseFieldCouplingPostProcessor3D::processSubDomain
void processSubDomain(BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
Execute post-processing step on a sublattice.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:164
olb::PhaseFieldCouplingPostProcessor3D::process
void process(BlockLattice< T, DESCRIPTOR > &blockLattice) override
Execute post-processing step.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:270
olb::PhaseFieldCouplingPostProcessor3D::PhaseFieldCouplingPostProcessor3D
PhaseFieldCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T rho_L, T rho_H, T mu_L, T mu_H, T surface_tension, T interface_thickness, std::vector< BlockStructureD< 3 > * > partners_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:151
olb::PorousNavierStokesAdvectionDiffusionCouplingGenerator3D::PorousNavierStokesAdvectionDiffusionCouplingGenerator3D
PorousNavierStokesAdvectionDiffusionCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:666
olb::PorousNavierStokesAdvectionDiffusionCouplingGenerator3D::generate
PostProcessor3D< T, DESCRIPTOR > * generate(std::vector< BlockStructureD< 3 > * > partners) const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:673
olb::PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D::PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D
PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_, std::vector< BlockStructureD< 3 > * > partners_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:598
olb::PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D::processSubDomain
void processSubDomain(BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
Execute post-processing step on a sublattice.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:619
olb::PorousNavierStokesAdvectionDiffusionCouplingPostProcessor3D::process
void process(BlockLattice< T, DESCRIPTOR > &blockLattice) override
Execute post-processing step.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:659
olb::PostProcessor3D::getName
std::string & getName()
read and write access to name
Definition postProcessing.hh:192
olb::SmagorinskyBoussinesqCouplingGenerator3D::clone
LatticeCouplingGenerator3D< T, DESCRIPTOR > * clone() const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:416
olb::SmagorinskyBoussinesqCouplingGenerator3D::SmagorinskyBoussinesqCouplingGenerator3D
SmagorinskyBoussinesqCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_, T PrTurb_, T smagoPrefactor_)
LatticeCouplingGenerator for advectionDiffusion coupling.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:401
olb::SmagorinskyBoussinesqCouplingGenerator3D::generate
PostProcessor3D< T, DESCRIPTOR > * generate(std::vector< BlockStructureD< 3 > * > partners) const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:408
olb::SmagorinskyBoussinesqCouplingPostProcessor3D::process
void process(BlockLattice< T, DESCRIPTOR > &blockLattice) override
Execute post-processing step.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:392
olb::SmagorinskyBoussinesqCouplingPostProcessor3D::SmagorinskyBoussinesqCouplingPostProcessor3D
SmagorinskyBoussinesqCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_, T PrTurb_, T smagoPrefactor_, std::vector< BlockStructureD< 3 > * > partners_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:306
olb::SmagorinskyBoussinesqCouplingPostProcessor3D::processSubDomain
void processSubDomain(BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
Execute post-processing step on a sublattice.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:335
olb::TotalEnthalpyPhaseChangeCouplingGenerator3D::TotalEnthalpyPhaseChangeCouplingGenerator3D
TotalEnthalpyPhaseChangeCouplingGenerator3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_)
LatticeCouplingGenerator for advectionDiffusion coupling.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:124
olb::TotalEnthalpyPhaseChangeCouplingGenerator3D::clone
LatticeCouplingGenerator3D< T, DESCRIPTOR > * clone() const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:139
olb::TotalEnthalpyPhaseChangeCouplingGenerator3D::generate
PostProcessor3D< T, DESCRIPTOR > * generate(std::vector< BlockStructureD< 3 > * > partners) const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:131
olb::TotalEnthalpyPhaseChangeCouplingPostProcessor3D::TotalEnthalpyPhaseChangeCouplingPostProcessor3D
TotalEnthalpyPhaseChangeCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T gravity_, T T0_, T deltaTemp_, std::vector< T > dir_, std::vector< BlockStructureD< 3 > * > partners_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:44
olb::TotalEnthalpyPhaseChangeCouplingPostProcessor3D::process
void process(BlockLattice< T, DESCRIPTOR > &blockLattice) override
Execute post-processing step.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:115
olb::TotalEnthalpyPhaseChangeCouplingPostProcessor3D::processSubDomain
void processSubDomain(BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
Execute post-processing step on a sublattice.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:71
olb::Vector
Plain old scalar vector.
Definition vector.h:47
olb::Vector::data
constexpr const T * data() const any_platform
Definition vector.h:161
olb::VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingGenerator3D::generate
PostProcessor3D< T, DESCRIPTOR > * generate(std::vector< BlockStructureD< 3 > * > partners) const override
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:1168
olb::VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D::processSubDomain
void processSubDomain(BlockLattice< T, DESCRIPTOR > &blockLattice, int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) override
Execute post-processing step on a sublattice.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:1018
olb::VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D::VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D
VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, int iC_, std::vector< BlockStructureD< 3 > * > partners_, std::vector< std::reference_wrapper< AdvectionDiffusionForce3D< T, DESCRIPTOR, ADLattice > > > forces_)
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:999
olb::VolumeAveragedNavierStokesAdvectionDiffusionParticleCouplingPostProcessor3D::process
void process(BlockLattice< T, DESCRIPTOR > &blockLattice) override
Execute post-processing step.
Definition navierStokesAdvectionDiffusionCouplingPostProcessor3D.hh:1154
latticeDescriptors.h
Descriptor for all types of 2D and 3D lattices.
olb::util::sqrt
cpu::simd::Pack< T > sqrt(cpu::simd::Pack< T > value)
Definition pack.h:100
olb::util::min
cpu::simd::Pack< T > min(cpu::simd::Pack< T > rhs, cpu::simd::Pack< T > lhs)
Definition pack.h:124
olb::util::max
cpu::simd::Pack< T > max(cpu::simd::Pack< T > rhs, cpu::simd::Pack< T > lhs)
Definition pack.h:130
olb::util::pow
cpu::simd::Pack< T > pow(cpu::simd::Pack< T > base, cpu::simd::Pack< T > exp)
Definition pack.h:112
olb::util::intersect
bool intersect(int x0, int x1, int y0, int y1, int x0_, int x1_, int y0_, int y1_, int &newX0, int &newX1, int &newY0, int &newY1)
Definition util.h:89
olb
Top level namespace for all of OpenLB.
Definition boundaryPostProcessors2D.h:34
olb::norm
constexpr T norm(const ScalarVector< T, D, IMPL > &a)
Euclidean vector norm.
Definition scalarVector.h:66
olb::ParametersOfDynamicsD
typename ParametersD< typename DYNAMICS::value_t, typename DYNAMICS::descriptor_t >::template include< typename DYNAMICS::parameters > ParametersOfDynamicsD
Deduce ParametersD of DYNAMICS w.r.t. its value type and descriptor.
Definition fieldParametersD.h:190
olb::util::TensorVal
Compute number of elements of a symmetric d-dimensional tensor.
Definition util.h:210
util.h
Set of functions commonly used in LB computations – header file.
Generated on Tue Mar 5 2024 09:57:47 for OpenLB by doxygen 1.10.0