Line data Source code
1 : /**********************************************************************/
2 : /* DO NOT MODIFY THIS HEADER */
3 : /* Swift, a Fourier spectral solver for MOOSE */
4 : /* */
5 : /* Copyright 2024 Battelle Energy Alliance, LLC */
6 : /* ALL RIGHTS RESERVED */
7 : /**********************************************************************/
8 :
9 : #include "LBMBounceBack.h"
10 : #include "LatticeBoltzmannProblem.h"
11 : #include "LatticeBoltzmannStencilBase.h"
12 :
13 : using namespace torch::indexing;
14 :
15 : registerMooseObject("SwiftApp", LBMBounceBack);
16 :
17 : InputParameters
18 0 : LBMBounceBack::validParams()
19 : {
20 0 : InputParameters params = LBMBoundaryCondition::validParams();
21 0 : params.addClassDescription("LBMBounceBack object");
22 0 : params.addRequiredParam<TensorInputBufferName>("f_old", "Old state distribution function");
23 0 : params.addParam<bool>(
24 : "exclude_corners_x",
25 0 : false,
26 : "Whether or not apply bounceback in the corners of the domain along x axis");
27 0 : params.addParam<bool>(
28 : "exclude_corners_y",
29 0 : false,
30 : "Whether or not apply bounceback in the corners of the domain along y axis");
31 0 : params.addParam<bool>(
32 : "exclude_corners_z",
33 0 : false,
34 : "Whether or not apply bounceback in the corners of the domain along z axis");
35 0 : return params;
36 0 : }
37 :
38 0 : LBMBounceBack::LBMBounceBack(const InputParameters & parameters)
39 : : LBMBoundaryCondition(parameters),
40 0 : _f_old(_lb_problem.getBufferOld(getParam<TensorInputBufferName>("f_old"), 1)),
41 0 : _exclude_corners_x(getParam<bool>("exclude_corners_x")),
42 0 : _exclude_corners_y(getParam<bool>("exclude_corners_y")),
43 0 : _exclude_corners_z(getParam<bool>("exclude_corners_z"))
44 : {
45 0 : if (_exclude_corners_x)
46 0 : _x_indices = torch::arange(1, _grid_size[0] - 1, MooseTensor::intTensorOptions());
47 : else
48 0 : _x_indices = torch::arange(_grid_size[0], MooseTensor::intTensorOptions());
49 :
50 0 : if (_exclude_corners_y)
51 0 : _y_indices = torch::arange(1, _grid_size[1] - 1, MooseTensor::intTensorOptions());
52 : else
53 0 : _y_indices = torch::arange(_grid_size[1], MooseTensor::intTensorOptions());
54 :
55 0 : if (_exclude_corners_z)
56 0 : _z_indices = torch::arange(1, _grid_size[2] - 1, MooseTensor::intTensorOptions());
57 : else
58 0 : _z_indices = torch::arange(_grid_size[2], MooseTensor::intTensorOptions());
59 :
60 0 : std::vector<torch::Tensor> xyz_mesh = torch::meshgrid({_x_indices, _y_indices, _z_indices});
61 :
62 0 : torch::Tensor flat_x_indices = xyz_mesh[0].reshape(-1);
63 0 : torch::Tensor flat_y_indices = xyz_mesh[1].reshape(-1);
64 0 : torch::Tensor flat_z_indices = xyz_mesh[2].reshape(-1);
65 :
66 0 : _x_indices = flat_x_indices.clone();
67 0 : _y_indices = flat_y_indices.clone();
68 0 : _z_indices = flat_z_indices.clone();
69 :
70 : // for 3D, binary media
71 0 : if (_lb_problem.isBinaryMedia())
72 : {
73 : const torch::Tensor & binary_mesh = _lb_problem.getBinaryMedia();
74 0 : _binary_mesh = binary_mesh.clone();
75 :
76 0 : for (int64_t ic = 1; ic < _stencil._q; ic++)
77 : {
78 0 : int64_t ex = _stencil._ex[ic].item<int64_t>();
79 0 : int64_t ey = _stencil._ey[ic].item<int64_t>();
80 0 : int64_t ez = _stencil._ez[ic].item<int64_t>();
81 0 : torch::Tensor shifted_mesh = torch::roll(binary_mesh, {ex, ey, ez}, {0, 1, 2});
82 0 : torch::Tensor adjacent_to_boundary = (shifted_mesh == 0) & (binary_mesh == 1);
83 0 : _binary_mesh.masked_fill_(adjacent_to_boundary, 2);
84 : }
85 : }
86 0 : }
87 :
88 : void
89 0 : LBMBounceBack::backBoundary()
90 : {
91 0 : for (unsigned int i = 0; i < _stencil._front.size(0); i++)
92 : {
93 0 : const auto & opposite_dir = _stencil._op[_stencil._front[i]];
94 0 : _u.index_put_({_x_indices, _y_indices, _grid_size[2] - 1, opposite_dir},
95 0 : _f_old[0].index({_x_indices, _y_indices, _grid_size[2] - 1, _stencil._front[i]}));
96 : }
97 0 : }
98 :
99 : void
100 0 : LBMBounceBack::frontBoundary()
101 : {
102 0 : for (unsigned int i = 0; i < _stencil._front.size(0); i++)
103 : {
104 0 : const auto & opposite_dir = _stencil._op[_stencil._front[i]];
105 0 : _u.index_put_({_x_indices, _y_indices, 0, _stencil._front[i]},
106 0 : _f_old[0].index({_x_indices, _y_indices, 0, opposite_dir}));
107 : }
108 0 : }
109 :
110 : void
111 0 : LBMBounceBack::leftBoundary()
112 : {
113 0 : for (unsigned int i = 0; i < _stencil._left.size(0); i++)
114 : {
115 0 : const auto & opposite_dir = _stencil._op[_stencil._left[i]];
116 0 : _u.index_put_({0, _y_indices, _z_indices, _stencil._left[i]},
117 0 : _f_old[0].index({0, _y_indices, _z_indices, opposite_dir}));
118 : }
119 0 : }
120 :
121 : void
122 0 : LBMBounceBack::rightBoundary()
123 : {
124 0 : for (unsigned int i = 0; i < _stencil._left.size(0); i++)
125 : {
126 0 : const auto & opposite_dir = _stencil._op[_stencil._left[i]];
127 0 : _u.index_put_({_grid_size[0] - 1, _y_indices, _z_indices, opposite_dir},
128 0 : _f_old[0].index({_grid_size[0] - 1, _y_indices, _z_indices, _stencil._left[i]}));
129 : }
130 0 : }
131 :
132 : void
133 0 : LBMBounceBack::bottomBoundary()
134 : {
135 0 : for (unsigned int i = 0; i < _stencil._bottom.size(0); i++)
136 : {
137 0 : const auto & opposite_dir = _stencil._op[_stencil._bottom[i]];
138 0 : _u.index_put_({_x_indices, 0, _z_indices, _stencil._bottom[i]},
139 0 : _f_old[0].index({_x_indices, 0, _z_indices, opposite_dir}));
140 : }
141 0 : }
142 :
143 : void
144 0 : LBMBounceBack::topBoundary()
145 : {
146 0 : for (unsigned int i = 0; i < _stencil._bottom.size(0); i++)
147 : {
148 0 : const auto & opposite_dir = _stencil._op[_stencil._bottom[i]];
149 0 : _u.index_put_(
150 : {_x_indices, _grid_size[1] - 1, _z_indices, opposite_dir},
151 0 : _f_old[0].index({_x_indices, _grid_size[1] - 1, _z_indices, _stencil._bottom[i]}));
152 : }
153 0 : }
154 :
155 : void
156 0 : LBMBounceBack::wallBoundary()
157 : {
158 0 : if (_domain.getDim() == 3)
159 0 : wallBoundary3D(); // temporary solution to generalization problem
160 : else
161 : // bounce-back
162 0 : _u.index_put_(
163 0 : {_boundary_indices.index({Slice(), 0}),
164 0 : _boundary_indices.index({Slice(), 1}),
165 0 : _boundary_indices.index({Slice(), 2}),
166 0 : _boundary_indices.index({Slice(), 3})},
167 :
168 0 : _f_old[0].index({_boundary_indices.index({Slice(), 0}),
169 0 : _boundary_indices.index({Slice(), 1}),
170 0 : _boundary_indices.index({Slice(), 2}),
171 0 : _stencil._op.index_select(0, _boundary_indices.index({Slice(), 3}))}));
172 0 : }
173 :
174 : void
175 0 : LBMBounceBack::wallBoundary3D()
176 : {
177 0 : if (_lb_problem.getTotalSteps() == 0)
178 : {
179 0 : _boundary_mask = (_binary_mesh.unsqueeze(-1).expand_as(_u) == 2) & (_u == 0);
180 0 : _boundary_mask = _boundary_mask.to(torch::kBool);
181 : }
182 :
183 0 : torch::Tensor f_bounce_back = torch::zeros_like(_u);
184 :
185 0 : for (/* do not use unsigned int */ int ic = 1; ic < _stencil._q; ic++)
186 : {
187 0 : int64_t index = _stencil._op[ic].item<int64_t>();
188 0 : auto lattice_slice = _f_old[0].index({Slice(), Slice(), Slice(), index});
189 0 : auto bounce_back_slice = f_bounce_back.index({Slice(), Slice(), Slice(), ic});
190 : bounce_back_slice.copy_(lattice_slice);
191 : }
192 :
193 0 : _u.index_put_({_boundary_mask}, f_bounce_back.index({_boundary_mask}));
194 0 : }
195 :
196 : void
197 0 : LBMBounceBack::computeBuffer()
198 : {
199 0 : const auto n_old = _f_old.size();
200 0 : if (n_old != 0)
201 : {
202 : // do not overwrite previous
203 0 : _u = _u.clone();
204 :
205 0 : switch (_boundary)
206 : {
207 0 : case Boundary::top:
208 0 : topBoundary();
209 0 : break;
210 0 : case Boundary::bottom:
211 0 : bottomBoundary();
212 0 : break;
213 0 : case Boundary::left:
214 0 : leftBoundary();
215 0 : break;
216 0 : case Boundary::right:
217 0 : rightBoundary();
218 0 : break;
219 0 : case Boundary::front:
220 0 : frontBoundary();
221 0 : break;
222 0 : case Boundary::back:
223 0 : backBoundary();
224 0 : break;
225 0 : case Boundary::wall:
226 0 : wallBoundary3D();
227 0 : break;
228 0 : default:
229 0 : mooseError("Undefined boundary names");
230 : }
231 : }
232 0 : _lb_problem.maskedFillSolids(_u, 0);
233 0 : }
|