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 "LBMFixedZerothOrderBC.h"
10 : #include "LatticeBoltzmannProblem.h"
11 : #include "LatticeBoltzmannStencilBase.h"
12 :
13 : using namespace torch::indexing;
14 :
15 : registerMooseObject("SwiftApp", LBMFixedZerothOrderBC);
16 :
17 : InputParameters
18 0 : LBMFixedZerothOrderBC::validParams()
19 : {
20 0 : InputParameters params = LBMBoundaryCondition::validParams();
21 0 : params.addClassDescription("LBMFixedZerothOrderBC object");
22 0 : params.addRequiredParam<TensorInputBufferName>("f", "Input buffer distribution function");
23 0 : params.addRequiredParam<std::string>("value", "Fixed input value");
24 0 : return params;
25 0 : }
26 :
27 0 : LBMFixedZerothOrderBC::LBMFixedZerothOrderBC(const InputParameters & parameters)
28 : : LBMBoundaryCondition(parameters),
29 0 : _f(getInputBufferByName(getParam<TensorInputBufferName>("f"))),
30 0 : _grid_size(_lb_problem.getGridSize()),
31 0 : _value(_lb_problem.getConstant<Real>(getParam<std::string>("value")))
32 : {
33 0 : }
34 :
35 : void
36 0 : LBMFixedZerothOrderBC::frontBoundary()
37 : {
38 0 : if (_domain.getDim() == 2)
39 0 : mooseError("There is no front boundary in 2 dimensions.");
40 : else
41 0 : mooseError("Front boundary is not implemented, but it can be replaced by any other boundary by "
42 : "rotating the domain.");
43 : }
44 :
45 : void
46 0 : LBMFixedZerothOrderBC::backBoundary()
47 : {
48 0 : if (_domain.getDim() == 2)
49 0 : mooseError("There is no back boundary in 2 dimensions.");
50 : else
51 0 : mooseError("Back boundary is not implemented, but it can be replaced by any other boundary by "
52 : "rotating the domain.");
53 : }
54 :
55 : void
56 0 : LBMFixedZerothOrderBC::leftBoundaryD2Q9()
57 : {
58 : torch::Tensor velocity =
59 0 : 1.0 - (_f.index({0, Slice(), Slice(), 0}) + _f.index({0, Slice(), Slice(), 2}) +
60 0 : _f.index({0, Slice(), Slice(), 4}) +
61 0 : 2 * (_f.index({0, Slice(), Slice(), 3}) + _f.index({0, Slice(), Slice(), 6}) +
62 0 : _f.index({0, Slice(), Slice(), 7}))) /
63 0 : _value;
64 :
65 : // axis aligned direction
66 0 : const auto & opposite_dir = _stencil._op[_stencil._left[0]];
67 0 : _u.index_put_({0, Slice(), Slice(), _stencil._left[0]},
68 0 : _f.index({0, Slice(), Slice(), opposite_dir}) + 2.0 / 3.0 * _value * velocity);
69 :
70 : // other directions
71 0 : for (unsigned int i = 1; i < _stencil._left.size(0); i++)
72 : {
73 0 : const auto & opposite_dir = _stencil._op[_stencil._left[i]];
74 0 : _u.index_put_(
75 0 : {0, Slice(), Slice(), _stencil._left[i]},
76 0 : _f.index({0, Slice(), Slice(), opposite_dir}) -
77 0 : 0.5 * _stencil._ey[_stencil._left[i]] *
78 0 : (_f.index({0, Slice(), Slice(), 2}) - _f.index({0, Slice(), Slice(), 4})) +
79 0 : 1.0 / 6.0 * _value * velocity);
80 : }
81 0 : }
82 :
83 : void
84 0 : LBMFixedZerothOrderBC::leftBoundary()
85 : {
86 0 : if (_stencil._q == 9)
87 0 : leftBoundaryD2Q9();
88 : else
89 : {
90 : torch::Tensor velocity =
91 0 : 1.0 - (torch::sum(_f.index({0, Slice(), Slice(), -_stencil._neutral_x}), -1) +
92 0 : 2 * torch::sum(_f.index({0, Slice(), Slice(), _stencil._right}), -1)) /
93 0 : _value;
94 :
95 0 : _u.index_put_({0, Slice(), Slice(), _stencil._left[0]},
96 0 : _f.index({0, Slice(), Slice(), _stencil._right[0]}) +
97 0 : 2.0 * _stencil._weights[_stencil._left[0]] / _lb_problem._cs2 * _value *
98 : velocity);
99 :
100 0 : for (unsigned int i = 1; i < _stencil._left.size(0); i++)
101 : {
102 0 : _u.index_put_({0, Slice(), Slice(), _stencil._left[i]},
103 0 : _f.index({0, Slice(), Slice(), _stencil._right[i]}) +
104 0 : 2.0 * _stencil._weights[_stencil._left[i]] / _lb_problem._cs2 * _value *
105 : velocity);
106 : }
107 : }
108 0 : }
109 :
110 : void
111 0 : LBMFixedZerothOrderBC::rightBoundaryD2Q9()
112 : {
113 0 : torch::Tensor velocity = (_f.index({_grid_size[0] - 1, Slice(), Slice(), 0}) +
114 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), 2}) +
115 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), 4}) +
116 0 : 2 * (_f.index({_grid_size[0] - 1, Slice(), Slice(), 1}) +
117 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), 5}) +
118 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), 8}))) /
119 0 : _value -
120 0 : 1.0;
121 :
122 : // axis aligned direction
123 0 : const auto & opposite_dir = _stencil._op[_stencil._left[0]];
124 0 : _u.index_put_({_grid_size[0] - 1, Slice(), Slice(), opposite_dir},
125 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), _stencil._left[0]}) -
126 0 : 2.0 / 3.0 * _value * velocity);
127 :
128 : // other directions
129 0 : for (unsigned int i = 1; i < _stencil._left.size(0); i++)
130 : {
131 0 : const auto & opposite_dir = _stencil._op[_stencil._left[i]];
132 0 : _u.index_put_({_grid_size[0] - 1, Slice(), Slice(), opposite_dir},
133 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), _stencil._left[i]}) +
134 0 : 0.5 * _stencil._ey[opposite_dir] *
135 0 : (_f.index({_grid_size[0] - 1, Slice(), Slice(), 4}) -
136 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), 2})) -
137 0 : 1.0 / 6.0 * _value * velocity);
138 : }
139 0 : }
140 :
141 : void
142 0 : LBMFixedZerothOrderBC::rightBoundary()
143 : {
144 0 : if (_stencil._q == 9)
145 0 : rightBoundaryD2Q9();
146 : else
147 : {
148 : torch::Tensor velocity =
149 0 : (torch::sum(_f.index({_grid_size[0] - 1, Slice(), Slice(), -_stencil._neutral_x}), -1) +
150 0 : 2 * torch::sum(_f.index({_grid_size[0] - 1, Slice(), Slice(), _stencil._left}), -1)) /
151 0 : _value -
152 0 : 1.0;
153 :
154 0 : _u.index_put_({_grid_size[0] - 1, Slice(), Slice(), _stencil._right[0]},
155 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), _stencil._left[0]}) -
156 0 : 2.0 * _stencil._weights[_stencil._right[0]] / _lb_problem._cs2 * _value *
157 : velocity);
158 :
159 0 : for (unsigned int i = 1; i < _stencil._right.size(0); i++)
160 : {
161 0 : _u.index_put_({_grid_size[0] - 1, Slice(), Slice(), _stencil._right[i]},
162 0 : _f.index({_grid_size[0] - 1, Slice(), Slice(), _stencil._left[i]}) -
163 0 : 2.0 * _stencil._weights[_stencil._right[i]] / _lb_problem._cs2 * _value *
164 : velocity);
165 : }
166 : }
167 0 : }
168 :
169 : void
170 0 : LBMFixedZerothOrderBC::bottomBoundaryD2Q9()
171 : {
172 : torch::Tensor velocity =
173 0 : 1.0 - (_f.index({Slice(), 0, Slice(), 0}) + _f.index({Slice(), 0, Slice(), 1}) +
174 0 : _f.index({Slice(), 0, Slice(), 3}) +
175 0 : 2 * (_f.index({Slice(), 0, Slice(), 4}) + _f.index({Slice(), 0, Slice(), 7}) +
176 0 : _f.index({Slice(), 0, Slice(), 8}))) /
177 0 : _value;
178 :
179 : // axis aligned direction
180 0 : const auto & opposite_dir = _stencil._op[_stencil._bottom[0]];
181 0 : _u.index_put_({Slice(), 0, Slice(), _stencil._bottom[0]},
182 0 : _f.index({Slice(), 0, Slice(), opposite_dir}) + 2.0 / 3.0 * _value * velocity);
183 :
184 : // other directions
185 0 : for (unsigned int i = 1; i < _stencil._bottom.size(0); i++)
186 : {
187 0 : const auto & opposite_dir = _stencil._op[_stencil._bottom[i]];
188 0 : _u.index_put_(
189 0 : {Slice(), 0, Slice(), _stencil._bottom[i]},
190 0 : _f.index({Slice(), 0, Slice(), opposite_dir}) -
191 0 : 0.5 * _stencil._ex[_stencil._bottom[i]] *
192 0 : (_f.index({Slice(), 0, Slice(), 1}) - _f.index({Slice(), 0, Slice(), 3})) +
193 0 : 1.0 / 6.0 * _value * velocity);
194 : }
195 0 : }
196 :
197 : void
198 0 : LBMFixedZerothOrderBC::bottomBoundary()
199 : {
200 : // TBD
201 0 : if (_stencil._q == 9)
202 0 : bottomBoundaryD2Q9();
203 : else
204 0 : mooseError("Bottom boundary is not implemented, but it can be replaced by any other boundary "
205 : "by rotating the domain");
206 0 : }
207 :
208 : void
209 0 : LBMFixedZerothOrderBC::topBoundaryD2Q9()
210 : {
211 0 : torch::Tensor velocity = (_f.index({Slice(), _grid_size[1] - 1, Slice(), 0}) +
212 0 : _f.index({Slice(), _grid_size[1] - 1, Slice(), 1}) +
213 0 : _f.index({Slice(), _grid_size[1] - 1, Slice(), 3}) +
214 0 : 2 * (_f.index({Slice(), _grid_size[1] - 1, Slice(), 2}) +
215 0 : _f.index({Slice(), _grid_size[1] - 1, Slice(), 5}) +
216 0 : _f.index({Slice(), _grid_size[1] - 1, Slice(), 6}))) /
217 0 : _value -
218 0 : 1.0;
219 :
220 : // axis aligned direction
221 0 : const auto & opposite_dir = _stencil._op[_stencil._bottom[0]];
222 0 : _u.index_put_({Slice(), _grid_size[1] - 1, Slice(), opposite_dir},
223 0 : _f.index({Slice(), _grid_size[1] - 1, Slice(), _stencil._bottom[0]}) -
224 0 : 2.0 / 3.0 * _value * velocity);
225 :
226 : // other directions
227 0 : for (unsigned int i = 1; i < _stencil._bottom.size(0); i++)
228 : {
229 0 : const auto & opposite_dir = _stencil._op[_stencil._bottom[i]];
230 0 : _u.index_put_({Slice(), _grid_size[1] - 1, Slice(), opposite_dir},
231 0 : _f.index({Slice(), _grid_size[1] - 1, Slice(), _stencil._bottom[i]}) +
232 0 : 0.5 * _stencil._ex[opposite_dir] *
233 0 : (_f.index({Slice(), _grid_size[1] - 1, Slice(), 3}) -
234 0 : _f.index({Slice(), _grid_size[1] - 1, Slice(), 1})) -
235 0 : 1.0 / 6.0 * _value * velocity);
236 : }
237 0 : }
238 :
239 : void
240 0 : LBMFixedZerothOrderBC::topBoundary()
241 : {
242 : // TBD
243 0 : if (_stencil._q == 9)
244 0 : topBoundaryD2Q9();
245 : else
246 0 : mooseError("Top boundary is not implemented, but it can be replaced by any other boundary by "
247 : "rotating the domain");
248 0 : }
249 :
250 : void
251 0 : LBMFixedZerothOrderBC::computeBuffer()
252 : {
253 0 : switch (_boundary)
254 : {
255 0 : case Boundary::top:
256 0 : topBoundary();
257 0 : break;
258 0 : case Boundary::bottom:
259 0 : bottomBoundary();
260 0 : break;
261 0 : case Boundary::left:
262 0 : leftBoundary();
263 0 : break;
264 0 : case Boundary::right:
265 0 : rightBoundary();
266 0 : break;
267 0 : case Boundary::front:
268 0 : frontBoundary();
269 0 : break;
270 0 : case Boundary::back:
271 0 : backBoundary();
272 0 : break;
273 0 : case Boundary::wall:
274 0 : wallBoundary();
275 0 : break;
276 0 : default:
277 0 : mooseError("Undefined boundary names");
278 : }
279 0 : _lb_problem.maskedFillSolids(_u, 0);
280 0 : }
|
