VirtualFluids 0.2.0
Parallel CFD LBM Solver
Loading...
Searching...
No Matches
SlipTurbulentViscosityCompressible.cu
Go to the documentation of this file.
1//=======================================================================================
2// ____ ____ __ ______ __________ __ __ __ __
3// \ \ | | | | | _ \ |___ ___| | | | | / \ | |
4// \ \ | | | | | |_) | | | | | | | / \ | |
5// \ \ | | | | | _ / | | | | | | / /\ \ | |
6// \ \ | | | | | | \ \ | | | \__/ | / ____ \ | |____
7// \ \ | | |__| |__| \__\ |__| \________/ /__/ \__\ |_______|
8// \ \ | | ________________________________________________________________
9// \ \ | | | ______________________________________________________________|
10// \ \| | | | __ __ __ __ ______ _______
11// \ | | |_____ | | | | | | | | | _ \ / _____)
12// \ | | _____| | | | | | | | | | | \ \ \_______
13// \ | | | | |_____ | \_/ | | | | |_/ / _____ |
14// \ _____| |__| |________| \_______/ |__| |______/ (_______/
15//
16// This file is part of VirtualFluids. VirtualFluids is free software: you can
17// redistribute it and/or modify it under the terms of the GNU General Public
18// License as published by the Free Software Foundation, either version 3 of
19// the License, or (at your option) any later version.
20//
21// VirtualFluids is distributed in the hope that it will be useful, but WITHOUT
22// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
23// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
24// for more details.
25//
26// SPDX-License-Identifier: GPL-3.0-or-later
27// SPDX-FileCopyrightText: Copyright © VirtualFluids Project contributors, see AUTHORS.md in root folder
28//
29//! \addtogroup gpu_BoundaryConditions BoundaryConditions
30//! \ingroup gpu_core core
31//! \{
32//! \author Martin Schoenherr, Anna Wellmann
33//======================================================================================
34#include "Calculation/Calculation.h"
35#include "lbm/constants/D3Q27.h"
36#include "basics/constants/NumericConstants.h"
37#include "Utilities/KernelUtilities.h"
38#include "cuda_helper/CudaIndexCalculation.h"
39
40using namespace vf::basics::constant;
41using namespace vf::lbm::dir;
42
43namespace vf::gpu {
44
45__global__ void SlipTurbulentViscosityCompressible_Device(
46 real* distributions,
47 int* subgridDistanceIndices,
48 real* subgridDistances,
49 unsigned int numberOfBCnodes,
50 real omega,
51 unsigned int* neighborX,
52 unsigned int* neighborY,
53 unsigned int* neighborZ,
54 real* turbViscosity,
55 unsigned long long numberOfLBnodes,
56 bool isEvenTimestep)
57{
58 //! The slip boundary condition is executed in the following steps
59 //!
60
61 ////////////////////////////////////////////////////////////////////////////////
62 //! - Get node index coordinates from threadIdx, blockIdx, blockDim and gridDim.
63 //!
64 const unsigned nodeIndex = vf::cuda::get1DIndexFrom2DBlock();
65
66 if(nodeIndex < numberOfBCnodes)
67 {
68 //////////////////////////////////////////////////////////////////////////
69 //! - Read distributions: style of reading and writing the distributions from/to stored arrays dependent on timestep is based on the esoteric twist algorithm \ref
70 //! <a href="https://doi.org/10.3390/computation5020019"><b>[ M. Geier et al. (2017), DOI:10.3390/computation5020019 ]</b></a>
71 //!
72 Distributions27 dist;
73 getPointersToDistributions(dist, distributions, numberOfLBnodes, isEvenTimestep);
74
75 ////////////////////////////////////////////////////////////////////////////////
76 //! - Set local subgrid distances (q's)
77 //!
78 SubgridDistances27 subgridD;
79 getPointersToSubgridDistances(subgridD, subgridDistances, numberOfBCnodes);
80
81 ////////////////////////////////////////////////////////////////////////////////
82 //! - Set neighbor indices (necessary for indirect addressing)
83 //!
84 unsigned int indexOfBCnode = subgridDistanceIndices[nodeIndex];
85 unsigned int kzero= indexOfBCnode;
86 unsigned int ke = indexOfBCnode;
87 unsigned int kw = neighborX[indexOfBCnode];
88 unsigned int kn = indexOfBCnode;
89 unsigned int ks = neighborY[indexOfBCnode];
90 unsigned int kt = indexOfBCnode;
91 unsigned int kb = neighborZ[indexOfBCnode];
92 unsigned int ksw = neighborY[kw];
93 unsigned int kne = indexOfBCnode;
94 unsigned int kse = ks;
95 unsigned int knw = kw;
96 unsigned int kbw = neighborZ[kw];
97 unsigned int kte = indexOfBCnode;
98 unsigned int kbe = kb;
99 unsigned int ktw = kw;
100 unsigned int kbs = neighborZ[ks];
101 unsigned int ktn = indexOfBCnode;
102 unsigned int kbn = kb;
103 unsigned int kts = ks;
104 unsigned int ktse = ks;
105 unsigned int kbnw = kbw;
106 unsigned int ktnw = kw;
107 unsigned int kbse = kbs;
108 unsigned int ktsw = ksw;
109 unsigned int kbne = kb;
110 unsigned int ktne = indexOfBCnode;
111 unsigned int kbsw = neighborZ[ksw];
112
113 ////////////////////////////////////////////////////////////////////////////////
114 //! - Set local distributions
115 //!
116 real f_W = (dist.f[dP00])[ke ];
117 real f_E = (dist.f[dM00])[kw ];
118 real f_S = (dist.f[d0P0])[kn ];
119 real f_N = (dist.f[d0M0])[ks ];
120 real f_B = (dist.f[d00P])[kt ];
121 real f_T = (dist.f[d00M])[kb ];
122 real f_SW = (dist.f[dPP0])[kne ];
123 real f_NE = (dist.f[dMM0])[ksw ];
124 real f_NW = (dist.f[dPM0])[kse ];
125 real f_SE = (dist.f[dMP0])[knw ];
126 real f_BW = (dist.f[dP0P])[kte ];
127 real f_TE = (dist.f[dM0M])[kbw ];
128 real f_TW = (dist.f[dP0M])[kbe ];
129 real f_BE = (dist.f[dM0P])[ktw ];
130 real f_BS = (dist.f[d0PP])[ktn ];
131 real f_TN = (dist.f[d0MM])[kbs ];
132 real f_TS = (dist.f[d0PM])[kbn ];
133 real f_BN = (dist.f[d0MP])[kts ];
134 real f_BSW = (dist.f[dPPP])[ktne ];
135 real f_BNE = (dist.f[dMMP])[ktsw ];
136 real f_BNW = (dist.f[dPMP])[ktse ];
137 real f_BSE = (dist.f[dMPP])[ktnw ];
138 real f_TSW = (dist.f[dPPM])[kbne ];
139 real f_TNE = (dist.f[dMMM])[kbsw ];
140 real f_TNW = (dist.f[dPMM])[kbse ];
141 real f_TSE = (dist.f[dMPM])[kbnw ];
142
143 ////////////////////////////////////////////////////////////////////////////////
144 //! - Calculate macroscopic quantities
145 //!
146 real drho = f_TSE + f_TNW + f_TNE + f_TSW + f_BSE + f_BNW + f_BNE + f_BSW +
147 f_BN + f_TS + f_TN + f_BS + f_BE + f_TW + f_TE + f_BW + f_SE + f_NW + f_NE + f_SW +
148 f_T + f_B + f_N + f_S + f_E + f_W + ((dist.f[d000])[kzero]);
149
150 real vx1 = (((f_TSE - f_BNW) - (f_TNW - f_BSE)) + ((f_TNE - f_BSW) - (f_TSW - f_BNE)) +
151 ((f_BE - f_TW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_NE - f_SW)) +
152 (f_E - f_W)) / (c1o1 + drho);
153
154 real vx2 = ((-(f_TSE - f_BNW) + (f_TNW - f_BSE)) + ((f_TNE - f_BSW) - (f_TSW - f_BNE)) +
155 ((f_BN - f_TS) + (f_TN - f_BS)) + (-(f_SE - f_NW) + (f_NE - f_SW)) +
156 (f_N - f_S)) / (c1o1 + drho);
157
158 real vx3 = (((f_TSE - f_BNW) + (f_TNW - f_BSE)) + ((f_TNE - f_BSW) + (f_TSW - f_BNE)) +
159 (-(f_BN - f_TS) + (f_TN - f_BS)) + ((f_TE - f_BW) - (f_BE - f_TW)) +
160 (f_T - f_B)) / (c1o1 + drho);
161
162 real cu_sq = c3o2 * (vx1 * vx1 + vx2 * vx2 + vx3 * vx3) * (c1o1 + drho);
163
164 ////////////////////////////////////////////////////////////////////////////////
165 //! - change the pointer to write the results in the correct array
166 //!
167 getPointersToDistributions(dist, distributions, numberOfLBnodes, !isEvenTimestep);
168
169 ////////////////////////////////////////////////////////////////////////////////
170 //! - compute local relaxation rate
171 //!
172 real om_turb = omega / (c1o1 + c3o1* omega* max(c0o1, turbViscosity[indexOfBCnode]) );
173
174 ////////////////////////////////////////////////////////////////////////////////
175 //! - Multiply the local velocities by the slipLength
176 //!
177 real slipLength = c1o1;
178 real VeloX = slipLength*vx1;
179 real VeloY = slipLength*vx2;
180 real VeloZ = slipLength*vx3;
181
182 ////////////////////////////////////////////////////////////////////////////////
183 //! - Update distributions with subgrid distance (q) between zero and one
184 //!
185 real feq, q, velocityLB, velocityBC;
186
187 bool x = false;
188 bool y = false;
189 bool z = false;
190
191 q = (subgridD.q[dP00])[nodeIndex];
192 if (q>=c0o1 && q<=c1o1) // only update distribution for q between zero and one
193 {
194 VeloX = c0o1;
195 x = true;
196
197 velocityLB = vx1;
198 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c2o27);
199 velocityBC = VeloX;
200 (dist.f[dM00])[kw] = getInterpolatedDistributionForVeloBC(q, f_E, f_W, feq, om_turb, velocityBC, c2o27);
201 }
202
203 q = (subgridD.q[dM00])[nodeIndex];
204 if (q>=c0o1 && q<=c1o1)
205 {
206 VeloX = c0o1;
207 x = true;
208
209 velocityLB = -vx1;
210 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c2o27);
211 velocityBC = -VeloX;
212 (dist.f[dP00])[ke] = getInterpolatedDistributionForVeloBC(q, f_W, f_E, feq, om_turb, velocityBC, c2o27);
213 }
214
215 q = (subgridD.q[d0P0])[nodeIndex];
216 if (q>=c0o1 && q<=c1o1)
217 {
218 VeloY = c0o1;
219 y = true;
220
221 velocityLB = vx2;
222 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c2o27);
223 velocityBC = VeloY;
224 (dist.f[d0M0])[ks] = getInterpolatedDistributionForVeloBC(q, f_N, f_S, feq, om_turb, velocityBC, c2o27);
225 }
226
227 q = (subgridD.q[d0M0])[nodeIndex];
228 if (q>=c0o1 && q<=c1o1)
229 {
230 VeloY = c0o1;
231 y = true;
232
233 velocityLB = -vx2;
234 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c2o27);
235 velocityBC = -VeloY;
236 (dist.f[d0P0])[kn] = getInterpolatedDistributionForVeloBC(q, f_S, f_N, feq, om_turb, velocityBC, c2o27);
237 }
238
239 q = (subgridD.q[d00P])[nodeIndex];
240 if (q>=c0o1 && q<=c1o1)
241 {
242 VeloZ = c0o1;
243 z = true;
244
245 velocityLB = vx3;
246 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c2o27);
247 velocityBC = VeloZ;
248 (dist.f[d00M])[kb] = getInterpolatedDistributionForVeloBC(q, f_T, f_B, feq, om_turb, velocityBC, c2o27);
249 }
250
251 q = (subgridD.q[d00M])[nodeIndex];
252 if (q>=c0o1 && q<=c1o1)
253 {
254 VeloZ = c0o1;
255 z = true;
256
257 velocityLB = -vx3;
258 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c2o27);
259 velocityBC = -VeloZ;
260 (dist.f[d00P])[kt] = getInterpolatedDistributionForVeloBC(q, f_B, f_T, feq, om_turb, velocityBC, c2o27);
261 }
262
263 q = (subgridD.q[dPP0])[nodeIndex];
264 if (q>=c0o1 && q<=c1o1)
265 {
266 VeloX = slipLength*vx1;
267 VeloY = slipLength*vx2;
268 if (x == true) VeloX = c0o1;
269 if (y == true) VeloY = c0o1;
270
271 velocityLB = vx1 + vx2;
272 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
273 velocityBC = VeloX + VeloY;
274 (dist.f[dMM0])[ksw] = getInterpolatedDistributionForVeloBC(q, f_NE, f_SW, feq, om_turb, velocityBC, c1o54);
275 }
276
277 q = (subgridD.q[dMM0])[nodeIndex];
278 if (q>=c0o1 && q<=c1o1)
279 {
280 VeloX = slipLength*vx1;
281 VeloY = slipLength*vx2;
282 if (x == true) VeloX = c0o1;
283 if (y == true) VeloY = c0o1;
284
285 velocityLB = -vx1 - vx2;
286 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
287 velocityBC = -VeloX - VeloY;
288 (dist.f[dPP0])[kne] = getInterpolatedDistributionForVeloBC(q, f_SW, f_NE, feq, om_turb, velocityBC, c1o54);
289 }
290
291 q = (subgridD.q[dPM0])[nodeIndex];
292 if (q>=c0o1 && q<=c1o1)
293 {
294 VeloX = slipLength*vx1;
295 VeloY = slipLength*vx2;
296 if (x == true) VeloX = c0o1;
297 if (y == true) VeloY = c0o1;
298
299 velocityLB = vx1 - vx2;
300 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
301 velocityBC = VeloX - VeloY;
302 (dist.f[dMP0])[knw] = getInterpolatedDistributionForVeloBC(q, f_SE, f_NW, feq, om_turb, velocityBC, c1o54);
303 }
304
305 q = (subgridD.q[dMP0])[nodeIndex];
306 if (q>=c0o1 && q<=c1o1)
307 {
308 VeloX = slipLength*vx1;
309 VeloY = slipLength*vx2;
310 if (x == true) VeloX = c0o1;
311 if (y == true) VeloY = c0o1;
312
313 velocityLB = -vx1 + vx2;
314 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
315 velocityBC = -VeloX + VeloY;
316 (dist.f[dPM0])[kse] = getInterpolatedDistributionForVeloBC(q, f_NW, f_SE, feq, om_turb, velocityBC, c1o54);
317 }
318
319 q = (subgridD.q[dP0P])[nodeIndex];
320 if (q>=c0o1 && q<=c1o1)
321 {
322 VeloX = slipLength*vx1;
323 VeloZ = slipLength*vx3;
324 if (x == true) VeloX = c0o1;
325 if (z == true) VeloZ = c0o1;
326
327 velocityLB = vx1 + vx3;
328 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
329 velocityBC = VeloX + VeloZ;
330 (dist.f[dM0M])[kbw] = getInterpolatedDistributionForVeloBC(q, f_TE, f_BW, feq, om_turb, velocityBC, c1o54);
331 }
332
333 q = (subgridD.q[dM0M])[nodeIndex];
334 if (q>=c0o1 && q<=c1o1)
335 {
336 VeloX = slipLength*vx1;
337 VeloZ = slipLength*vx3;
338 if (x == true) VeloX = c0o1;
339 if (z == true) VeloZ = c0o1;
340
341 velocityLB = -vx1 - vx3;
342 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
343 velocityBC = -VeloX - VeloZ;
344 (dist.f[dP0P])[kte] = getInterpolatedDistributionForVeloBC(q, f_BW, f_TE, feq, om_turb, velocityBC, c1o54);
345 }
346
347 q = (subgridD.q[dP0M])[nodeIndex];
348 if (q>=c0o1 && q<=c1o1)
349 {
350 VeloX = slipLength*vx1;
351 VeloZ = slipLength*vx3;
352 if (x == true) VeloX = c0o1;
353 if (z == true) VeloZ = c0o1;
354
355 velocityLB = vx1 - vx3;
356 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
357 velocityBC = VeloX - VeloZ;
358 (dist.f[dM0P])[ktw] = getInterpolatedDistributionForVeloBC(q, f_BE, f_TW, feq, om_turb, velocityBC, c1o54);
359 }
360
361 q = (subgridD.q[dM0P])[nodeIndex];
362 if (q>=c0o1 && q<=c1o1)
363 {
364 VeloX = slipLength*vx1;
365 VeloZ = slipLength*vx3;
366 if (x == true) VeloX = c0o1;
367 if (z == true) VeloZ = c0o1;
368
369 velocityLB = -vx1 + vx3;
370 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
371 velocityBC = -VeloX + VeloZ;
372 (dist.f[dP0M])[kbe] = getInterpolatedDistributionForVeloBC(q, f_TW, f_BE, feq, om_turb, velocityBC, c1o54);
373 }
374
375 q = (subgridD.q[d0PP])[nodeIndex];
376 if (q>=c0o1 && q<=c1o1)
377 {
378 VeloY = slipLength*vx2;
379 VeloZ = slipLength*vx3;
380 if (y == true) VeloY = c0o1;
381 if (z == true) VeloZ = c0o1;
382
383 velocityLB = vx2 + vx3;
384 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
385 velocityBC = VeloY + VeloZ;
386 (dist.f[d0MM])[kbs] = getInterpolatedDistributionForVeloBC(q, f_TN, f_BS, feq, om_turb, velocityBC, c1o54);
387 }
388
389 q = (subgridD.q[d0MM])[nodeIndex];
390 if (q>=c0o1 && q<=c1o1)
391 {
392 VeloY = slipLength*vx2;
393 VeloZ = slipLength*vx3;
394 if (y == true) VeloY = c0o1;
395 if (z == true) VeloZ = c0o1;
396
397 velocityLB = -vx2 - vx3;
398 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
399 velocityBC = -VeloY - VeloZ;
400 (dist.f[d0PP])[ktn] = getInterpolatedDistributionForVeloBC(q, f_BS, f_TN, feq, om_turb, velocityBC, c1o54);
401 }
402
403
404 q = (subgridD.q[d0PM])[nodeIndex];
405 if (q>=c0o1 && q<=c1o1)
406 {
407 VeloY = slipLength*vx2;
408 VeloZ = slipLength*vx3;
409 if (y == true) VeloY = c0o1;
410 if (z == true) VeloZ = c0o1;
411
412 velocityLB = vx2 - vx3;
413 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
414 velocityBC = VeloY - VeloZ;
415 (dist.f[d0MP])[kts] = getInterpolatedDistributionForVeloBC(q, f_BN, f_TS, feq, om_turb, velocityBC, c1o54);
416 }
417
418 q = (subgridD.q[d0MP])[nodeIndex];
419 if (q>=c0o1 && q<=c1o1)
420 {
421 VeloY = slipLength*vx2;
422 VeloZ = slipLength*vx3;
423 if (y == true) VeloY = c0o1;
424 if (z == true) VeloZ = c0o1;
425
426 velocityLB = -vx2 + vx3;
427 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o54);
428 velocityBC = -VeloY + VeloZ;
429 (dist.f[d0PM])[kbn] = getInterpolatedDistributionForVeloBC(q, f_TS, f_BN, feq, om_turb, velocityBC, c1o54);
430 }
431
432 q = (subgridD.q[dPPP])[nodeIndex];
433 if (q>=c0o1 && q<=c1o1)
434 {
435 VeloX = slipLength*vx1;
436 VeloY = slipLength*vx2;
437 VeloZ = slipLength*vx3;
438 if (x == true) VeloX = c0o1;
439 if (y == true) VeloY = c0o1;
440 if (z == true) VeloZ = c0o1;
441 velocityLB = vx1 + vx2 + vx3;
442 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
443 velocityBC = VeloX + VeloY + VeloZ;
444 (dist.f[dMMM])[kbsw] = getInterpolatedDistributionForVeloBC(q, f_TNE, f_BSW, feq, om_turb, velocityBC, c1o216);
445 }
446
447 q = (subgridD.q[dMMM])[nodeIndex];
448 if (q>=c0o1 && q<=c1o1)
449 {
450 VeloX = slipLength*vx1;
451 VeloY = slipLength*vx2;
452 VeloZ = slipLength*vx3;
453 if (x == true) VeloX = c0o1;
454 if (y == true) VeloY = c0o1;
455 if (z == true) VeloZ = c0o1;
456 velocityLB = -vx1 - vx2 - vx3;
457 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
458 velocityBC = -VeloX - VeloY - VeloZ;
459 (dist.f[dPPP])[ktne] = getInterpolatedDistributionForVeloBC(q, f_BSW, f_TNE, feq, om_turb, velocityBC, c1o216);
460 }
461
462
463 q = (subgridD.q[dPPM])[nodeIndex];
464 if (q>=c0o1 && q<=c1o1)
465 {
466 VeloX = slipLength*vx1;
467 VeloY = slipLength*vx2;
468 VeloZ = slipLength*vx3;
469 if (x == true) VeloX = c0o1;
470 if (y == true) VeloY = c0o1;
471 if (z == true) VeloZ = c0o1;
472 velocityLB = vx1 + vx2 - vx3;
473 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
474 velocityBC = VeloX + VeloY - VeloZ;
475 (dist.f[dMMP])[ktsw] = getInterpolatedDistributionForVeloBC(q, f_BNE, f_TSW, feq, om_turb, velocityBC, c1o216);
476 }
477
478 q = (subgridD.q[dMMP])[nodeIndex];
479 if (q>=c0o1 && q<=c1o1)
480 {
481 VeloX = slipLength*vx1;
482 VeloY = slipLength*vx2;
483 VeloZ = slipLength*vx3;
484 if (x == true) VeloX = c0o1;
485 if (y == true) VeloY = c0o1;
486 if (z == true) VeloZ = c0o1;
487 velocityLB = -vx1 - vx2 + vx3;
488 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
489 velocityBC = -VeloX - VeloY + VeloZ;
490 (dist.f[dPPM])[kbne] = getInterpolatedDistributionForVeloBC(q, f_TSW, f_BNE, feq, om_turb, velocityBC, c1o216);
491 }
492
493 q = (subgridD.q[dPMP])[nodeIndex];
494 if (q>=c0o1 && q<=c1o1)
495 {
496 VeloX = slipLength*vx1;
497 VeloY = slipLength*vx2;
498 VeloZ = slipLength*vx3;
499 if (x == true) VeloX = c0o1;
500 if (y == true) VeloY = c0o1;
501 if (z == true) VeloZ = c0o1;
502 velocityLB = vx1 - vx2 + vx3;
503 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
504 velocityBC = VeloX - VeloY + VeloZ;
505 (dist.f[dMPM])[kbnw] = getInterpolatedDistributionForVeloBC(q, f_TSE, f_BNW, feq, om_turb, velocityBC, c1o216);
506 }
507
508 q = (subgridD.q[dMPM])[nodeIndex];
509 if (q>=c0o1 && q<=c1o1)
510 {
511 VeloX = slipLength*vx1;
512 VeloY = slipLength*vx2;
513 VeloZ = slipLength*vx3;
514 if (x == true) VeloX = c0o1;
515 if (y == true) VeloY = c0o1;
516 if (z == true) VeloZ = c0o1;
517 velocityLB = -vx1 + vx2 - vx3;
518 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
519 velocityBC = -VeloX + VeloY - VeloZ;
520 (dist.f[dPMP])[ktse] = getInterpolatedDistributionForVeloBC(q, f_BNW, f_TSE, feq, om_turb, velocityBC, c1o216);
521 }
522
523 q = (subgridD.q[dPMM])[nodeIndex];
524 if (q>=c0o1 && q<=c1o1)
525 {
526 VeloX = slipLength*vx1;
527 VeloY = slipLength*vx2;
528 VeloZ = slipLength*vx3;
529 if (x == true) VeloX = c0o1;
530 if (y == true) VeloY = c0o1;
531 if (z == true) VeloZ = c0o1;
532 velocityLB = vx1 - vx2 - vx3;
533 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
534 velocityBC = VeloX - VeloY - VeloZ;
535 (dist.f[dMPP])[ktnw] = getInterpolatedDistributionForVeloBC(q, f_BSE, f_TNW, feq, om_turb, velocityBC, c1o216);
536 }
537
538 q = (subgridD.q[dMPP])[nodeIndex];
539 if (q>=c0o1 && q<=c1o1)
540 {
541 VeloX = slipLength*vx1;
542 VeloY = slipLength*vx2;
543 VeloZ = slipLength*vx3;
544 if (x == true) VeloX = c0o1;
545 if (y == true) VeloY = c0o1;
546 if (z == true) VeloZ = c0o1;
547 velocityLB = -vx1 + vx2 + vx3;
548 feq = getEquilibriumForBC(drho, velocityLB, cu_sq, c1o216);
549 velocityBC = -VeloX + VeloY + VeloZ;
550 (dist.f[dPMM])[kbse] = getInterpolatedDistributionForVeloBC(q, f_TNW, f_BSE, feq, om_turb, velocityBC, c1o216);
551 }
552 }
553}
554
555}
556
557//! \}