1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
|
#include "lattice.h"
#include "LLBM/collide.h"
#include "LLBM/initialize.h"
#include "LLBM/bounce_back.h"
#include "LLBM/bounce_back_moving_wall.h"
#include "LLBM/equilibrium_density_wall.h"
#include "LLBM/equilibrium_velocity_wall.h"
#include <chrono>
#include <iostream>
#include "pattern/all.h"
using T = SWEEPLB_PRECISION;
using PATTERN = pattern::SWEEPLB_PATTERN<T>;
void simulate(Cuboid cuboid, std::size_t nStep) {
const int nThread = omp_get_max_threads();
Lattice<PATTERN> lattice(cuboid);
LatticeMask<T> bulk_mask(cuboid.volume());
LatticeMask<T> wall_mask(cuboid.volume());
LatticeMask<T> inflow_mask(cuboid.volume());
LatticeMask<T> outflow_mask(cuboid.volume());
cuboid.traverse([&](int iX, int iY, int iZ, std::size_t iCell) {
if ( iY == 0 || iY == cuboid[1]-1
|| iZ == 0 || iZ == cuboid[2]-1) {
wall_mask.set(iCell, true);
} else if (iX == 0) {
inflow_mask.set(iCell, true);
} else if (iX == cuboid[0]-1) {
outflow_mask.set(iCell, true);
} else {
bulk_mask.set(iCell, true);
}
});
bulk_mask.serialize();
wall_mask.serialize();
inflow_mask.serialize();
outflow_mask.serialize();
T tau = 0.56;
T u_inflow = 0.05;
T d_outflow = 1.;
for (std::size_t iStep = 0; iStep < 100; ++iStep) {
lattice.apply(Operator(BgkCollideO(), bulk_mask, tau),
Operator(BounceBackO(), wall_mask),
Operator(EquilibriumVelocityWallO(), inflow_mask, u_inflow, WallNormal<1,0,0>()),
Operator(EquilibriumDensityWallO(), outflow_mask, d_outflow, WallNormal<-1,0,0>()));
lattice.stream();
}
auto start = std::chrono::steady_clock::now();
for (std::size_t iStep = 0; iStep < nStep; ++iStep) {
lattice.apply(Operator(BgkCollideO(), bulk_mask, tau),
Operator(BounceBackO(), wall_mask),
Operator(EquilibriumVelocityWallO(), inflow_mask, u_inflow, WallNormal<1,0,0>()),
Operator(EquilibriumDensityWallO(), outflow_mask, d_outflow, WallNormal<-1,0,0>()));
lattice.stream();
}
auto duration = std::chrono::duration_cast<std::chrono::duration<double>>(
std::chrono::steady_clock::now() - start);
std::cout << cuboid[0] << ", " << cuboid[1] << ", " << cuboid[2]
<< ", " << nStep
<< ", " << nThread
<< ", " << (nStep * lattice.volume()) / (1e6 * duration.count())
<< std::endl;
lattice.write_momenta(bulk_mask, "result.vtk");
}
int main(int argc, char* argv[]) {
const std::size_t nX = atoi(argv[1]);
const std::size_t nY = atoi(argv[2]);
const std::size_t nZ = atoi(argv[3]);
const std::size_t steps = atoi(argv[4]);
simulate({ nX, nY, nZ }, steps);
}
|
