Line data Source code
1 : /********************************************************************/
2 : /* SOFTWARE COPYRIGHT NOTIFICATION */
3 : /* Cardinal */
4 : /* */
5 : /* (c) 2021 UChicago Argonne, LLC */
6 : /* ALL RIGHTS RESERVED */
7 : /* */
8 : /* Prepared by UChicago Argonne, LLC */
9 : /* Under Contract No. DE-AC02-06CH11357 */
10 : /* With the U. S. Department of Energy */
11 : /* */
12 : /* Prepared by Battelle Energy Alliance, LLC */
13 : /* Under Contract No. DE-AC07-05ID14517 */
14 : /* With the U. S. Department of Energy */
15 : /* */
16 : /* See LICENSE for full restrictions */
17 : /********************************************************************/
18 :
19 : #ifdef ENABLE_NEK_COUPLING
20 :
21 : #include "NekMeshDeformation.h"
22 : #include "DisplacedProblem.h"
23 :
24 : registerMooseObject("CardinalApp", NekMeshDeformation);
25 :
26 : InputParameters
27 86 : NekMeshDeformation::validParams()
28 : {
29 86 : auto params = FieldTransferBase::validParams();
30 86 : params.addClassDescription("Reads/writes mesh deformation between NekRS and MOOSE.");
31 86 : return params;
32 0 : }
33 :
34 45 : NekMeshDeformation::NekMeshDeformation(const InputParameters & parameters)
35 45 : : FieldTransferBase(parameters)
36 : {
37 45 : if (_direction == "to_nek" /* && nekrs::hasBlendingSolver() */)
38 : {
39 45 : if (_usrwrk_slot.size() != 3)
40 0 : paramError("usrwrk_slot", "For mesh deformation, 'usrwrk_slot' must be of length 3");
41 :
42 : // all the mesh velocities scale the same way
43 45 : auto d = nekrs::nondimensionalDivisor(field::mesh_velocity_x);
44 45 : auto a = nekrs::nondimensionalAdditive(field::mesh_velocity_x);
45 45 : addExternalVariable(_usrwrk_slot[0], _variable + "_x", a, d);
46 45 : addExternalVariable(_usrwrk_slot[1], _variable + "_y", a, d);
47 90 : addExternalVariable(_usrwrk_slot[2], _variable + "_z", a, d);
48 : }
49 : else
50 : {
51 : // technically, we do not need usrwrk_slot when not using the blending mesh solver,
52 : // since we directly apply the displacements to the entire mesh volume. However,
53 : // the base class is set up to require usrwrk_slot whenever we have displacements
54 : // being applied. Because the mesh deformation is not a widely used feature yet,
55 : // we ignore this annoyance. If someone wants to fix in future with a redesign
56 : // they are welcome to.
57 : // checkUnusedParam(parameters, "usrwrk_slot", "not using the blending mesh solver");
58 : }
59 :
60 45 : if (_direction == "from_nek")
61 0 : paramError("direction",
62 : "The NekMeshDeformation currently only supports transfers 'to_nek'; contact the "
63 : "Cardinal developer team if you require reading of NekRS mesh coordinates.");
64 :
65 45 : if (_app.actionWarehouse().displacedMesh() && !nekrs::hasMovingMesh())
66 1 : mooseWarning(
67 1 : "Your NekRSMesh has 'displacements', but '" + _nek_problem.casename() +
68 : ".par' does not have a\n"
69 : "solver in the [MESH] block! The displacements transferred to NekRS will be unused.");
70 :
71 44 : if (nekrs::hasMovingMesh() && _nek_mesh->exactMirror())
72 0 : mooseError("An exact mesh mirror is not yet implemented for the boundary mesh solver.");
73 :
74 44 : int n_per_surf = _nek_mesh->exactMirror() ? std::pow(_nek_mesh->nekNumQuadraturePoints1D(), 2.0)
75 44 : : _nek_mesh->numVerticesPerSurface();
76 44 : int n_per_vol = _nek_mesh->exactMirror() ? std::pow(_nek_mesh->nekNumQuadraturePoints1D(), 3.0)
77 44 : : _nek_mesh->numVerticesPerVolume();
78 :
79 44 : if (nekrs::hasMovingMesh())
80 : {
81 44 : int n_entries = _nek_mesh->volume() ? n_per_vol : n_per_surf;
82 44 : _displacement_x = (double *)calloc(n_entries, sizeof(double));
83 44 : _displacement_y = (double *)calloc(n_entries, sizeof(double));
84 44 : _displacement_z = (double *)calloc(n_entries, sizeof(double));
85 :
86 44 : if (nekrs::hasBlendingSolver())
87 18 : _mesh_velocity_elem = (double *)calloc(n_entries, sizeof(double));
88 : }
89 :
90 44 : auto boundary = _nek_mesh->boundary();
91 44 : if (!boundary && nekrs::hasBlendingSolver())
92 2 : mooseError("'" + _nek_problem.casename() +
93 : ".par' has a solver in the [MESH] block. This solver uses\n"
94 : "boundary displacement values from MOOSE to move the NekRS mesh. Please indicate\n"
95 : "the 'boundary' for which mesh motion is coupled from MOOSE to NekRS.");
96 :
97 43 : if (nekrs::hasBlendingSolver())
98 : {
99 : bool has_one_mv_bc = false;
100 18 : for (const auto & b : *boundary)
101 : {
102 17 : if (nekrs::isMovingMeshBoundary(b))
103 : {
104 : has_one_mv_bc = true;
105 : break;
106 : }
107 : }
108 :
109 17 : if (!has_one_mv_bc)
110 1 : mooseError("For boundary-coupled moving mesh problems, you need at least one "
111 1 : "boundary in '" +
112 1 : _nek_problem.casename() +
113 : ".par'\nto be of the type 'codedFixedValue'"
114 : " in the [MESH] block.");
115 : }
116 :
117 42 : if (!_nek_mesh->volume() && nekrs::hasUserMeshSolver())
118 1 : mooseError(
119 1 : "'" + _nek_problem.casename() +
120 : ".par' has 'solver = user' in the [MESH] block. With this solver,\n"
121 : "displacement values are sent to every GLL point in NekRS's volume. If you only are "
122 : "building\n"
123 : "a boundary mesh mirror, it's possible that some displacement values could result\n"
124 : "in negative Jacobians if a sideset moves beyond the bounds of an undeformed element.\n"
125 : "To eliminate this possibility, please enable 'volume = true' for NekRSMesh and send a\n"
126 : "whole-domain displacement to NekRS.");
127 :
128 41 : if (nekrs::hasBlendingSolver())
129 16 : _nek_mesh->initializePreviousDisplacements();
130 :
131 41 : if (nekrs::hasUserMeshSolver())
132 25 : _nek_mesh->saveInitialVolMesh();
133 41 : }
134 :
135 40 : NekMeshDeformation::~NekMeshDeformation()
136 : {
137 40 : freePointer(_displacement_x);
138 40 : freePointer(_displacement_y);
139 40 : freePointer(_displacement_z);
140 40 : freePointer(_mesh_velocity_elem);
141 40 : }
142 :
143 : void
144 904 : NekMeshDeformation::sendDataToNek()
145 : {
146 904 : if (nekrs::hasUserMeshSolver())
147 104 : sendVolumeDeformationToNek();
148 800 : else if (nekrs::hasBlendingSolver())
149 800 : sendBoundaryDeformationToNek();
150 : else
151 0 : mooseError("Unhandled mesh solver case in NekMeshDeformation!");
152 :
153 904 : _nek_problem.getDisplacedProblem()->updateMesh();
154 904 : }
155 :
156 : void
157 104 : NekMeshDeformation::sendVolumeDeformationToNek()
158 : {
159 104 : _console << "Sending volume deformation to NekRS" << std::endl;
160 :
161 104 : auto d = nekrs::nondimensionalDivisor(field::x_displacement);
162 104 : auto a = nekrs::nondimensionalAdditive(field::x_displacement);
163 6760 : for (unsigned int e = 0; e < _nek_mesh->numVolumeElems(); e++)
164 : {
165 : // We can only write into the nekRS scratch space if that face is "owned" by the current process
166 6656 : if (nekrs::commRank() != _nek_mesh->volumeCoupling().processor_id(e))
167 5632 : continue;
168 :
169 1024 : _nek_problem.mapVolumeDataToNekVolume(
170 1024 : e, _variable_number[_variable + "_x"], d, a, &_displacement_x);
171 1024 : _nek_problem.writeVolumeDisplacement(
172 1024 : e, _displacement_x, field::x_displacement, &(_nek_mesh->nek_initial_x()));
173 :
174 1024 : _nek_problem.mapVolumeDataToNekVolume(
175 1024 : e, _variable_number[_variable + "_y"], d, a, &_displacement_y);
176 1024 : _nek_problem.writeVolumeDisplacement(
177 1024 : e, _displacement_y, field::y_displacement, &(_nek_mesh->nek_initial_y()));
178 :
179 1024 : _nek_problem.mapVolumeDataToNekVolume(
180 1024 : e, _variable_number[_variable + "_z"], d, a, &_displacement_z);
181 1024 : _nek_problem.writeVolumeDisplacement(
182 1024 : e, _displacement_z, field::z_displacement, &(_nek_mesh->nek_initial_z()));
183 : }
184 104 : }
185 :
186 : void
187 800 : NekMeshDeformation::sendBoundaryDeformationToNek()
188 : {
189 800 : _console << "Sending boundary deformation to NekRS..." << std::endl;
190 :
191 800 : auto d = nekrs::nondimensionalDivisor(field::x_displacement);
192 800 : auto a = nekrs::nondimensionalAdditive(field::x_displacement);
193 800 : if (!_nek_mesh->volume())
194 : {
195 179600 : for (unsigned int e = 0; e < _nek_mesh->numSurfaceElems(); e++)
196 : {
197 : // We can only write into the nekRS scratch space if that face is "owned" by the current
198 : // process
199 179200 : if (nekrs::commRank() != _nek_mesh->boundaryCoupling().processor_id(e))
200 134400 : continue;
201 :
202 44800 : _nek_problem.mapFaceDataToNekFace(
203 44800 : e, _variable_number[_variable + "_x"], d, a, &_displacement_x);
204 44800 : calculateMeshVelocity(e, field::mesh_velocity_x);
205 44800 : _nek_problem.writeBoundarySolution(
206 44800 : _usrwrk_slot[0] * nekrs::fieldOffset(), e, _mesh_velocity_elem);
207 :
208 44800 : _nek_problem.mapFaceDataToNekFace(
209 44800 : e, _variable_number[_variable + "_y"], d, a, &_displacement_y);
210 44800 : calculateMeshVelocity(e, field::mesh_velocity_y);
211 44800 : _nek_problem.writeBoundarySolution(
212 44800 : _usrwrk_slot[1] * nekrs::fieldOffset(), e, _mesh_velocity_elem);
213 :
214 44800 : _nek_problem.mapFaceDataToNekFace(
215 44800 : e, _variable_number[_variable + "_z"], d, a, &_displacement_z);
216 44800 : calculateMeshVelocity(e, field::mesh_velocity_z);
217 44800 : _nek_problem.writeBoundarySolution(
218 44800 : _usrwrk_slot[2] * nekrs::fieldOffset(), e, _mesh_velocity_elem);
219 : }
220 : }
221 : else
222 : {
223 538000 : for (unsigned int e = 0; e < _nek_mesh->numVolumeElems(); ++e)
224 : {
225 : // We can only write into the nekRS scratch space if that face is "owned" by the current
226 : // process
227 537600 : if (nekrs::commRank() != _nek_mesh->volumeCoupling().processor_id(e))
228 403200 : continue;
229 :
230 134400 : _nek_problem.mapFaceDataToNekVolume(
231 134400 : e, _variable_number[_variable + "_x"], d, a, &_displacement_x);
232 134400 : calculateMeshVelocity(e, field::mesh_velocity_x);
233 134400 : _nek_problem.writeVolumeSolution(
234 134400 : _usrwrk_slot[0] * nekrs::fieldOffset(), e, _mesh_velocity_elem);
235 :
236 134400 : _nek_problem.mapFaceDataToNekVolume(
237 134400 : e, _variable_number[_variable + "_y"], d, a, &_displacement_y);
238 134400 : calculateMeshVelocity(e, field::mesh_velocity_y);
239 134400 : _nek_problem.writeVolumeSolution(
240 134400 : _usrwrk_slot[1] * nekrs::fieldOffset(), e, _mesh_velocity_elem);
241 :
242 134400 : _nek_problem.mapFaceDataToNekVolume(
243 134400 : e, _variable_number[_variable + "_z"], d, a, &_displacement_z);
244 134400 : calculateMeshVelocity(e, field::mesh_velocity_z);
245 134400 : _nek_problem.writeVolumeSolution(
246 134400 : _usrwrk_slot[2] * nekrs::fieldOffset(), e, _mesh_velocity_elem);
247 : }
248 : }
249 800 : }
250 :
251 : void
252 537600 : NekMeshDeformation::calculateMeshVelocity(int e, const field::NekWriteEnum & field)
253 : {
254 : int len =
255 537600 : _nek_mesh->volume() ? _nek_mesh->numVerticesPerVolume() : _nek_mesh->numVerticesPerSurface();
256 :
257 537600 : double dt = _nek_problem.transientExecutioner()->getTimeStepper()->getCurrentDT();
258 :
259 : double *displacement = nullptr, *prev_disp = nullptr;
260 : field::NekWriteEnum disp_field;
261 :
262 537600 : switch (field)
263 : {
264 179200 : case field::mesh_velocity_x:
265 179200 : displacement = _displacement_x;
266 : prev_disp = _nek_mesh->prev_disp_x().data();
267 : disp_field = field::x_displacement;
268 179200 : break;
269 179200 : case field::mesh_velocity_y:
270 179200 : displacement = _displacement_y;
271 : prev_disp = _nek_mesh->prev_disp_y().data();
272 : disp_field = field::y_displacement;
273 179200 : break;
274 179200 : case field::mesh_velocity_z:
275 179200 : displacement = _displacement_z;
276 : prev_disp = _nek_mesh->prev_disp_z().data();
277 : disp_field = field::z_displacement;
278 179200 : break;
279 0 : default:
280 0 : mooseError("Unhandled NekWriteEnum in NekRSProblem::calculateMeshVelocity!\n");
281 : }
282 :
283 537600 : auto reference_v = nekrs::nondimensionalDivisor(field::velocity);
284 12633600 : for (int i = 0; i < len; i++)
285 12096000 : _mesh_velocity_elem[i] = (displacement[i] - prev_disp[(e * len) + i]) / dt / reference_v;
286 :
287 537600 : _nek_mesh->updateDisplacement(e, displacement, disp_field);
288 537600 : }
289 :
290 : #endif
|
