Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "MultiAppProjectionTransfer.h"
11 :
12 : // MOOSE includes
13 : #include "AddVariableAction.h"
14 : #include "FEProblem.h"
15 : #include "MooseMesh.h"
16 : #include "MooseVariableFE.h"
17 : #include "SystemBase.h"
18 : #include "MooseAppCoordTransform.h"
19 :
20 : #include "libmesh/dof_map.h"
21 : #include "libmesh/linear_implicit_system.h"
22 : #include "libmesh/mesh_function.h"
23 : #include "libmesh/mesh_tools.h"
24 : #include "libmesh/numeric_vector.h"
25 : #include "libmesh/parallel_algebra.h"
26 : #include "libmesh/quadrature_gauss.h"
27 : #include "libmesh/sparse_matrix.h"
28 : #include "libmesh/string_to_enum.h"
29 : #include "libmesh/default_coupling.h"
30 :
31 : // TIMPI includes
32 : #include "timpi/parallel_sync.h"
33 :
34 : using namespace libMesh;
35 :
36 : void
37 544 : assemble_l2(EquationSystems & es, const std::string & system_name)
38 : {
39 : MultiAppProjectionTransfer * transfer =
40 544 : es.parameters.get<MultiAppProjectionTransfer *>("transfer");
41 544 : transfer->assembleL2(es, system_name);
42 544 : }
43 :
44 : registerMooseObject("MooseApp", MultiAppProjectionTransfer);
45 :
46 : InputParameters
47 14817 : MultiAppProjectionTransfer::validParams()
48 : {
49 14817 : InputParameters params = MultiAppConservativeTransfer::validParams();
50 14817 : params.addClassDescription(
51 : "Perform a projection between a master and sub-application mesh of a field variable.");
52 :
53 14817 : MooseEnum proj_type("l2", "l2");
54 14817 : params.addParam<MooseEnum>("proj_type", proj_type, "The type of the projection.");
55 :
56 44451 : params.addParam<bool>("fixed_meshes",
57 29634 : false,
58 : "Set to true when the meshes are not changing (ie, "
59 : "no movement or adaptivity). This will cache some "
60 : "information to speed up the transfer.");
61 :
62 : // Need one layer of ghosting
63 14817 : params.addRelationshipManager("ElementSideNeighborLayers",
64 : Moose::RelationshipManagerType::GEOMETRIC |
65 : Moose::RelationshipManagerType::ALGEBRAIC);
66 14817 : MultiAppTransfer::addBBoxFactorParam(params);
67 29634 : return params;
68 14817 : }
69 :
70 276 : MultiAppProjectionTransfer::MultiAppProjectionTransfer(const InputParameters & parameters)
71 : : MultiAppConservativeTransfer(parameters),
72 276 : _proj_type(getParam<MooseEnum>("proj_type")),
73 276 : _compute_matrix(true),
74 276 : _fixed_meshes(getParam<bool>("fixed_meshes")),
75 828 : _qps_cached(false)
76 : {
77 276 : if (_to_var_names.size() != 1)
78 0 : paramError("variable", " Support single to-variable only ");
79 :
80 276 : if (_from_var_names.size() != 1)
81 0 : paramError("source_variable", " Support single from-variable only ");
82 276 : }
83 :
84 : void
85 276 : MultiAppProjectionTransfer::initialSetup()
86 : {
87 276 : MultiAppConservativeTransfer::initialSetup();
88 :
89 276 : _proj_sys.resize(_to_problems.size(), NULL);
90 :
91 596 : for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
92 : {
93 320 : FEProblemBase & to_problem = *_to_problems[i_to];
94 320 : EquationSystems & to_es = to_problem.es();
95 :
96 : // Add the projection system.
97 : FEType fe_type = to_problem
98 320 : .getVariable(0,
99 : _to_var_name,
100 : Moose::VarKindType::VAR_ANY,
101 : Moose::VarFieldType::VAR_FIELD_STANDARD)
102 320 : .feType();
103 :
104 320 : LinearImplicitSystem & proj_sys = to_es.add_system<LinearImplicitSystem>("proj-sys-" + name());
105 :
106 640 : proj_sys.get_dof_map().add_coupling_functor(
107 320 : proj_sys.get_dof_map().default_coupling(),
108 : false); // The false keeps it from getting added to the mesh
109 :
110 320 : _proj_var_num = proj_sys.add_variable("var", fe_type);
111 320 : proj_sys.attach_assemble_function(assemble_l2);
112 320 : _proj_sys[i_to] = &proj_sys;
113 :
114 : // Prevent the projection system from being written to checkpoint
115 : // files. In the event of a recover or restart, we'll read the checkpoint
116 : // before this initialSetup method is called. As a result, we'll find
117 : // systems in the checkpoint (the projection systems) that we don't know
118 : // what to do with, and there will be a crash. We could fix this by making
119 : // the systems in the constructor, except we don't know how many sub apps
120 : // there are at the time of construction. So instead, we'll just nuke the
121 : // projection system and rebuild it from scratch every recover/restart.
122 320 : proj_sys.hide_output() = true;
123 :
124 : // Reinitialize EquationSystems since we added a system.
125 320 : to_es.reinit();
126 : }
127 276 : }
128 :
129 : void
130 544 : MultiAppProjectionTransfer::assembleL2(EquationSystems & es, const std::string & system_name)
131 : {
132 : // Get the system and mesh from the input arguments.
133 544 : LinearImplicitSystem & system = es.get_system<LinearImplicitSystem>(system_name);
134 544 : MeshBase & to_mesh = es.get_mesh();
135 :
136 : // Get the meshfunction evaluations and the map that was stashed in the es.
137 544 : std::vector<Real> & final_evals = *es.parameters.get<std::vector<Real> *>("final_evals");
138 : std::map<dof_id_type, unsigned int> & element_map =
139 544 : *es.parameters.get<std::map<dof_id_type, unsigned int> *>("element_map");
140 :
141 : // Setup system vectors and matrices.
142 544 : FEType fe_type = system.variable_type(0);
143 544 : std::unique_ptr<FEBase> fe(FEBase::build(to_mesh.mesh_dimension(), fe_type));
144 544 : QGauss qrule(to_mesh.mesh_dimension(), fe_type.default_quadrature_order());
145 544 : fe->attach_quadrature_rule(&qrule);
146 544 : const DofMap & dof_map = system.get_dof_map();
147 544 : DenseMatrix<Number> Ke;
148 544 : DenseVector<Number> Fe;
149 544 : std::vector<dof_id_type> dof_indices;
150 544 : const std::vector<Real> & JxW = fe->get_JxW();
151 544 : const std::vector<std::vector<Real>> & phi = fe->get_phi();
152 544 : auto & system_matrix = system.get_system_matrix();
153 :
154 35184 : for (const auto & elem : to_mesh.active_local_element_ptr_range())
155 : {
156 17320 : fe->reinit(elem);
157 :
158 17320 : dof_map.dof_indices(elem, dof_indices);
159 17320 : Ke.resize(dof_indices.size(), dof_indices.size());
160 17320 : Fe.resize(dof_indices.size());
161 :
162 72528 : for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
163 : {
164 55208 : Real meshfun_eval = 0.;
165 55208 : if (element_map.find(elem->id()) != element_map.end())
166 : {
167 : // We have evaluations for this element.
168 44568 : meshfun_eval = final_evals[element_map[elem->id()] + qp];
169 : }
170 :
171 : // Now compute the element matrix and RHS contributions.
172 263312 : for (unsigned int i = 0; i < phi.size(); i++)
173 : {
174 : // RHS
175 208104 : Fe(i) += JxW[qp] * (meshfun_eval * phi[i][qp]);
176 :
177 208104 : if (_compute_matrix)
178 1046160 : for (unsigned int j = 0; j < phi.size(); j++)
179 : {
180 : // The matrix contribution
181 838056 : Ke(i, j) += JxW[qp] * (phi[i][qp] * phi[j][qp]);
182 : }
183 : }
184 55208 : dof_map.constrain_element_matrix_and_vector(Ke, Fe, dof_indices);
185 :
186 55208 : if (_compute_matrix)
187 55208 : system_matrix.add_matrix(Ke, dof_indices);
188 55208 : system.rhs->add_vector(Fe, dof_indices);
189 : }
190 544 : }
191 544 : }
192 :
193 : void
194 443 : MultiAppProjectionTransfer::execute()
195 : {
196 443 : TIME_SECTION(
197 : "MultiAppProjectionTransfer::execute()", 5, "Transferring variables through projection");
198 :
199 : ////////////////////
200 : // We are going to project the solutions by solving some linear systems. In
201 : // order to assemble the systems, we need to evaluate the "from" domain
202 : // solutions at quadrature points in the "to" domain. Some parallel
203 : // communication is necessary because each processor doesn't necessarily have
204 : // all the "from" information it needs to set its "to" values. We don't want
205 : // to use a bunch of big all-to-all broadcasts, so we'll use bounding boxes to
206 : // figure out which processors have the information we need and only
207 : // communicate with those processors.
208 : //
209 : // Each processor will
210 : // 1. Check its local quadrature points in the "to" domains to see which
211 : // "from" domains they might be in.
212 : // 2. Send quadrature points to the processors with "from" domains that might
213 : // contain those points.
214 : // 3. Recieve quadrature points from other processors, evaluate its mesh
215 : // functions at those points, and send the values back to the proper
216 : // processor
217 : // 4. Recieve mesh function evaluations from all relevant processors and
218 : // decide which one to use at every quadrature point (the lowest global app
219 : // index always wins)
220 : // 5. And use the mesh function evaluations to assemble and solve an L2
221 : // projection system on its local elements.
222 : ////////////////////
223 :
224 : ////////////////////
225 : // For every combination of global "from" problem and local "to" problem, find
226 : // which "from" bounding boxes overlap with which "to" elements. Keep track
227 : // of which processors own bounding boxes that overlap with which elements.
228 : // Build vectors of quadrature points to send to other processors for mesh
229 : // function evaluations.
230 : ////////////////////
231 :
232 : // Get the bounding boxes for the "from" domains.
233 443 : std::vector<BoundingBox> bboxes = getFromBoundingBoxes();
234 :
235 : // Figure out how many "from" domains each processor owns.
236 443 : std::vector<unsigned int> froms_per_proc = getFromsPerProc();
237 :
238 443 : std::map<processor_id_type, std::vector<Point>> outgoing_qps;
239 : std::map<processor_id_type, std::map<std::pair<unsigned int, unsigned int>, unsigned int>>
240 443 : element_index_map;
241 : // element_index_map[i_to, element_id] = index
242 : // outgoing_qps[index] is the first quadrature point in element
243 :
244 443 : if (!_qps_cached)
245 : {
246 907 : for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
247 : {
248 : // Indexing into the coordinate transforms vector
249 : const auto to_global_num =
250 504 : _current_direction == FROM_MULTIAPP ? 0 : _to_local2global_map[i_to];
251 504 : MeshBase & to_mesh = _to_meshes[i_to]->getMesh();
252 :
253 504 : LinearImplicitSystem & system = *_proj_sys[i_to];
254 :
255 504 : FEType fe_type = system.variable_type(0);
256 504 : std::unique_ptr<FEBase> fe(FEBase::build(to_mesh.mesh_dimension(), fe_type));
257 504 : QGauss qrule(to_mesh.mesh_dimension(), fe_type.default_quadrature_order());
258 504 : fe->attach_quadrature_rule(&qrule);
259 504 : const std::vector<Point> & xyz = fe->get_xyz();
260 :
261 504 : unsigned int from0 = 0;
262 1230 : for (processor_id_type i_proc = 0; i_proc < n_processors();
263 726 : from0 += froms_per_proc[i_proc], i_proc++)
264 : {
265 726 : for (const auto & elem :
266 45712 : as_range(to_mesh.local_elements_begin(), to_mesh.local_elements_end()))
267 : {
268 22130 : fe->reinit(elem);
269 :
270 22130 : bool qp_hit = false;
271 45700 : for (unsigned int i_from = 0; i_from < froms_per_proc[i_proc] && !qp_hit; i_from++)
272 : {
273 69067 : for (unsigned int qp = 0; qp < qrule.n_points() && !qp_hit; qp++)
274 : {
275 45497 : Point qpt = xyz[qp];
276 45497 : if (bboxes[from0 + i_from].contains_point((*_to_transforms[to_global_num])(qpt)))
277 13177 : qp_hit = true;
278 : }
279 : }
280 :
281 22130 : if (qp_hit)
282 : {
283 : // The selected processor's bounding box contains at least one
284 : // quadrature point from this element. Send all qps from this element
285 : // and remember where they are in the array using the map.
286 13177 : std::pair<unsigned int, unsigned int> key(i_to, elem->id());
287 13177 : element_index_map[i_proc][key] = outgoing_qps[i_proc].size();
288 56643 : for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
289 : {
290 43466 : Point qpt = xyz[qp];
291 43466 : outgoing_qps[i_proc].push_back((*_to_transforms[to_global_num])(qpt));
292 : }
293 : }
294 726 : }
295 : }
296 504 : }
297 :
298 403 : if (_fixed_meshes)
299 48 : _cached_index_map = element_index_map;
300 : }
301 : else
302 : {
303 40 : element_index_map = _cached_index_map;
304 : }
305 :
306 : ////////////////////
307 : // Request quadrature point evaluations from other processors and handle
308 : // requests sent to this processor.
309 : ////////////////////
310 :
311 : // Get the local bounding boxes.
312 443 : std::vector<BoundingBox> local_bboxes(froms_per_proc[processor_id()]);
313 : {
314 : // Find the index to the first of this processor's local bounding boxes.
315 443 : unsigned int local_start = 0;
316 566 : for (processor_id_type i_proc = 0; i_proc < n_processors() && i_proc != processor_id();
317 : i_proc++)
318 123 : local_start += froms_per_proc[i_proc];
319 :
320 : // Extract the local bounding boxes.
321 906 : for (unsigned int i_from = 0; i_from < froms_per_proc[processor_id()]; i_from++)
322 463 : local_bboxes[i_from] = bboxes[local_start + i_from];
323 : }
324 :
325 : // Setup the local mesh functions.
326 443 : std::vector<MeshFunction> local_meshfuns;
327 906 : for (unsigned int i_from = 0; i_from < _from_problems.size(); i_from++)
328 : {
329 463 : FEProblemBase & from_problem = *_from_problems[i_from];
330 463 : MooseVariableFEBase & from_var = from_problem.getVariable(
331 : 0, _from_var_name, Moose::VarKindType::VAR_ANY, Moose::VarFieldType::VAR_FIELD_STANDARD);
332 463 : System & from_sys = from_var.sys().system();
333 463 : unsigned int from_var_num = from_sys.variable_number(from_var.name());
334 :
335 926 : local_meshfuns.emplace_back(
336 463 : from_problem.es(), *from_sys.current_local_solution, from_sys.get_dof_map(), from_var_num);
337 463 : local_meshfuns.back().init();
338 463 : local_meshfuns.back().enable_out_of_mesh_mode(OutOfMeshValue);
339 : }
340 :
341 : // Recieve quadrature points from other processors, evaluate mesh frunctions
342 : // at those points, and send the values back.
343 443 : std::map<processor_id_type, std::vector<std::pair<Real, unsigned int>>> outgoing_evals_ids;
344 :
345 : // If there is no cached data, we need to do communication
346 : // Quadrature points I will receive from remote processors
347 443 : std::map<processor_id_type, std::vector<Point>> incoming_qps;
348 443 : if (!_qps_cached)
349 : {
350 486 : auto qps_action_functor = [&incoming_qps](processor_id_type pid, const std::vector<Point> & qps)
351 : {
352 : // Quadrature points from processor 'pid'
353 486 : auto & incoming_qps_from_pid = incoming_qps[pid];
354 : // Store data for late use
355 486 : incoming_qps_from_pid.reserve(incoming_qps_from_pid.size() + qps.size());
356 486 : std::copy(qps.begin(), qps.end(), std::back_inserter(incoming_qps_from_pid));
357 889 : };
358 :
359 403 : Parallel::push_parallel_vector_data(comm(), outgoing_qps, qps_action_functor);
360 : }
361 :
362 : // Cache data
363 443 : if (!_qps_cached)
364 403 : _cached_qps = incoming_qps;
365 :
366 969 : for (auto & qps : _cached_qps)
367 : {
368 526 : const processor_id_type pid = qps.first;
369 :
370 1052 : outgoing_evals_ids[pid].resize(qps.second.size(),
371 526 : std::make_pair(OutOfMeshValue, libMesh::invalid_uint));
372 :
373 47492 : for (unsigned int qp = 0; qp < qps.second.size(); qp++)
374 : {
375 46966 : Point qpt = qps.second[qp];
376 :
377 : // Loop until we've found the lowest-ranked app that actually contains
378 : // the quadrature point.
379 94832 : for (unsigned int i_from = 0; i_from < _from_problems.size(); i_from++)
380 : {
381 47866 : if (local_bboxes[i_from].contains_point(qpt))
382 : {
383 : const auto from_global_num =
384 46966 : _current_direction == TO_MULTIAPP ? 0 : _from_local2global_map[i_from];
385 46966 : outgoing_evals_ids[pid][qp].first =
386 46966 : (local_meshfuns[i_from])(_from_transforms[from_global_num]->mapBack(qpt));
387 46966 : if (_current_direction == FROM_MULTIAPP)
388 10592 : outgoing_evals_ids[pid][qp].second = from_global_num;
389 : }
390 : }
391 : }
392 : }
393 :
394 : ////////////////////
395 : // Gather all of the qp evaluations and pick out the best ones for each qp.
396 : ////////////////////
397 :
398 : // Values back from remote processors for my local quadrature points
399 443 : std::map<processor_id_type, std::vector<std::pair<Real, unsigned int>>> incoming_evals_ids;
400 :
401 : auto evals_action_functor =
402 526 : [&incoming_evals_ids](processor_id_type pid,
403 526 : const std::vector<std::pair<Real, unsigned int>> & evals)
404 : {
405 : // evals for processor 'pid'
406 526 : auto & incoming_evals_ids_for_pid = incoming_evals_ids[pid];
407 : // Copy evals for late use
408 526 : incoming_evals_ids_for_pid.reserve(incoming_evals_ids_for_pid.size() + evals.size());
409 526 : std::copy(evals.begin(), evals.end(), std::back_inserter(incoming_evals_ids_for_pid));
410 969 : };
411 :
412 443 : Parallel::push_parallel_vector_data(comm(), outgoing_evals_ids, evals_action_functor);
413 :
414 443 : std::vector<std::vector<Real>> final_evals(_to_problems.size());
415 443 : std::vector<std::map<dof_id_type, unsigned int>> trimmed_element_maps(_to_problems.size());
416 :
417 987 : for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
418 : {
419 544 : MeshBase & to_mesh = _to_meshes[i_to]->getMesh();
420 544 : LinearImplicitSystem & system = *_proj_sys[i_to];
421 :
422 544 : FEType fe_type = system.variable_type(0);
423 544 : std::unique_ptr<FEBase> fe(FEBase::build(to_mesh.mesh_dimension(), fe_type));
424 544 : QGauss qrule(to_mesh.mesh_dimension(), fe_type.default_quadrature_order());
425 :
426 35184 : for (const auto & elem : to_mesh.active_local_element_ptr_range())
427 : {
428 17320 : qrule.init(*elem);
429 :
430 17320 : bool element_is_evaled = false;
431 17320 : std::vector<Real> evals(qrule.n_points(), 0.);
432 :
433 72528 : for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
434 : {
435 55208 : unsigned int lowest_app_rank = libMesh::invalid_uint;
436 123510 : for (auto & values_ids : incoming_evals_ids)
437 : {
438 : // Current processor id
439 68302 : const processor_id_type pid = values_ids.first;
440 :
441 : // Ignore the selected processor if the element wasn't found in it's
442 : // bounding box.
443 : std::map<std::pair<unsigned int, unsigned int>, unsigned int> & map =
444 68302 : element_index_map[pid];
445 68302 : std::pair<unsigned int, unsigned int> key(i_to, elem->id());
446 68302 : if (map.find(key) == map.end())
447 22142 : continue;
448 46966 : unsigned int qp0 = map[key];
449 :
450 : // Ignore the selected processor if it's app has a higher rank than the
451 : // previously found lowest app rank.
452 46966 : if (_current_direction == FROM_MULTIAPP)
453 10592 : if (values_ids.second[qp0 + qp].second >= lowest_app_rank)
454 0 : continue;
455 :
456 : // Ignore the selected processor if the qp was actually outside the
457 : // processor's subapp's mesh.
458 46966 : if (values_ids.second[qp0 + qp].first == OutOfMeshValue)
459 806 : continue;
460 :
461 : // This is the best meshfunction evaluation so far, save it.
462 46160 : element_is_evaled = true;
463 46160 : evals[qp] = values_ids.second[qp0 + qp].first;
464 : }
465 : }
466 :
467 : // If we found good evaluations for any of the qps in this element, save
468 : // those evaluations for later.
469 17320 : if (element_is_evaled)
470 : {
471 13904 : trimmed_element_maps[i_to][elem->id()] = final_evals[i_to].size();
472 58472 : for (unsigned int qp = 0; qp < qrule.n_points(); qp++)
473 44568 : final_evals[i_to].push_back(evals[qp]);
474 : }
475 17864 : }
476 544 : }
477 :
478 : ////////////////////
479 : // We now have just one or zero mesh function values at all of our local
480 : // quadrature points. Stash those values (and a map linking them to element
481 : // ids) in the equation systems parameters and project the solution.
482 : ////////////////////
483 :
484 987 : for (unsigned int i_to = 0; i_to < _to_problems.size(); i_to++)
485 : {
486 544 : _to_es[i_to]->parameters.set<std::vector<Real> *>("final_evals") = &final_evals[i_to];
487 544 : _to_es[i_to]->parameters.set<std::map<dof_id_type, unsigned int> *>("element_map") =
488 544 : &trimmed_element_maps[i_to];
489 544 : projectSolution(i_to);
490 544 : _to_es[i_to]->parameters.set<std::vector<Real> *>("final_evals") = NULL;
491 544 : _to_es[i_to]->parameters.set<std::map<dof_id_type, unsigned int> *>("element_map") = NULL;
492 : }
493 :
494 443 : if (_fixed_meshes)
495 88 : _qps_cached = true;
496 :
497 443 : postExecute();
498 443 : }
499 :
500 : void
501 544 : MultiAppProjectionTransfer::projectSolution(unsigned int i_to)
502 : {
503 544 : FEProblemBase & to_problem = *_to_problems[i_to];
504 544 : EquationSystems & proj_es = to_problem.es();
505 544 : LinearImplicitSystem & ls = *_proj_sys[i_to];
506 : // activate the current transfer
507 544 : proj_es.parameters.set<MultiAppProjectionTransfer *>("transfer") = this;
508 :
509 : // TODO: specify solver params in an input file
510 : // solver tolerance
511 544 : Real tol = proj_es.parameters.get<Real>("linear solver tolerance");
512 544 : proj_es.parameters.set<Real>("linear solver tolerance") = 1e-10; // set our tolerance
513 : // solve it
514 544 : ls.solve();
515 544 : proj_es.parameters.set<Real>("linear solver tolerance") = tol; // restore the original tolerance
516 :
517 : // copy projected solution into target es
518 544 : MeshBase & to_mesh = proj_es.get_mesh();
519 :
520 544 : MooseVariableFEBase & to_var = to_problem.getVariable(
521 : 0, _to_var_name, Moose::VarKindType::VAR_ANY, Moose::VarFieldType::VAR_FIELD_STANDARD);
522 544 : System & to_sys = to_var.sys().system();
523 544 : NumericVector<Number> * to_solution = to_sys.solution.get();
524 :
525 46792 : for (const auto & node : to_mesh.local_node_ptr_range())
526 : {
527 39738 : for (unsigned int comp = 0; comp < node->n_comp(to_sys.number(), to_var.number()); comp++)
528 : {
529 16614 : const dof_id_type proj_index = node->dof_number(ls.number(), _proj_var_num, comp);
530 16614 : const dof_id_type to_index = node->dof_number(to_sys.number(), to_var.number(), comp);
531 16614 : to_solution->set(to_index, (*ls.solution)(proj_index));
532 : }
533 544 : }
534 35184 : for (const auto & elem : to_mesh.active_local_element_ptr_range())
535 22776 : for (unsigned int comp = 0; comp < elem->n_comp(to_sys.number(), to_var.number()); comp++)
536 : {
537 5456 : const dof_id_type proj_index = elem->dof_number(ls.number(), _proj_var_num, comp);
538 5456 : const dof_id_type to_index = elem->dof_number(to_sys.number(), to_var.number(), comp);
539 5456 : to_solution->set(to_index, (*ls.solution)(proj_index));
540 544 : }
541 :
542 544 : to_solution->close();
543 544 : to_sys.update();
544 544 : }
|
