v0.13.1
VecManager.cpp
Go to the documentation of this file.
1/** \file VecManager.cpp
2 * \brief Implementation of VecManager
3 *
4 * MoFEM is free software: you can redistribute it and/or modify it under
5 * the terms of the GNU Lesser General Public License as published by the
6 * Free Software Foundation, either version 3 of the License, or (at your
7 * option) any later version.
8 *
9 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
10 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
12 * License for more details.
13 *
14 * You should have received a copy of the GNU Lesser General Public
15 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
16 */
17
18namespace MoFEM {
19
21VecManager::query_interface(boost::typeindex::type_index type_index,
22 UnknownInterface **iface) const {
23 *iface = const_cast<VecManager *>(this);
24 return 0;
25}
26
28 : cOre(const_cast<MoFEM::Core &>(core)), dEbug(false) {
29 if (!LogManager::checkIfChannelExist("VECWORLD")) {
30 auto core_log = logging::core::get();
31
32 core_log->add_sink(
34 core_log->add_sink(
36 core_log->add_sink(
38
39 LogManager::setLog("VECWORLD");
40 LogManager::setLog("VECSYNC");
41 LogManager::setLog("VECSELF");
42
43 MOFEM_LOG_TAG("VECWORLD", "VecMng");
44 MOFEM_LOG_TAG("VECSYNC", "VecMng");
45 MOFEM_LOG_TAG("VECSELF", "VecMng");
46 }
47
48 MOFEM_LOG("VECWORLD", Sev::noisy) << "Vector manager created";
49}
50
52 Vec *V) const {
53 const MoFEM::Interface &m_field = cOre;
54 const Problem *problem_ptr;
56 CHKERR m_field.get_problem(name, &problem_ptr);
57 DofIdx nb_local_dofs, nb_ghost_dofs;
58 switch (rc) {
59 case ROW:
60 nb_local_dofs = problem_ptr->getNbLocalDofsRow();
61 nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
62 break;
63 case COL:
64 nb_local_dofs = problem_ptr->getNbLocalDofsCol();
65 nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
66 break;
67 default:
68 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
69 }
70 CHKERR ::VecCreateSeq(PETSC_COMM_SELF, nb_local_dofs + nb_ghost_dofs, V);
72}
73
75 Vec *V) const {
76 const MoFEM::Interface &m_field = cOre;
77 const Problem *problem_ptr;
79 CHKERR m_field.get_problem(name, &problem_ptr);
80 DofIdx nb_dofs, nb_local_dofs, nb_ghost_dofs;
82 switch (rc) {
83 case ROW:
84 nb_dofs = problem_ptr->getNbDofsRow();
85 nb_local_dofs = problem_ptr->getNbLocalDofsRow();
86 nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
87 dofs = const_cast<NumeredDofEntityByLocalIdx *>(
88 &problem_ptr->numeredRowDofsPtr->get<PetscLocalIdx_mi_tag>());
89 break;
90 case COL:
91 nb_dofs = problem_ptr->getNbDofsCol();
92 nb_local_dofs = problem_ptr->getNbLocalDofsCol();
93 nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
94 dofs = const_cast<NumeredDofEntityByLocalIdx *>(
95 &problem_ptr->numeredColDofsPtr->get<PetscLocalIdx_mi_tag>());
96 break;
97 default:
98 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
99 }
100 // get ghost dofs
101 auto miit = dofs->lower_bound(nb_local_dofs);
102 auto hi_miit = dofs->upper_bound(nb_local_dofs + nb_ghost_dofs);
103 int count = std::distance(miit, hi_miit);
104 if (count != nb_ghost_dofs) {
105 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
106 }
107 std::vector<DofIdx> ghost_idx(count);
108 for (auto vit = ghost_idx.begin(); miit != hi_miit; ++miit, ++vit) {
109 *vit = (*miit)->petscGloablDofIdx;
110
111#ifndef NDEBUG
112 if (PetscUnlikely(*vit > nb_dofs || *vit < 0)) {
114 MOFEM_LOG_ATTRIBUTES("VECSELF",
116 MOFEM_LOG("VECSELF", Sev::error) << "Problem with DOF " << **miit;
117 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong ghost index");
118 }
119#endif
120
121 }
122 CHKERR ::VecCreateGhost(m_field.get_comm(), nb_local_dofs, nb_dofs,
123 nb_ghost_dofs, &ghost_idx[0], V);
125}
126
128 SmartPetscObj<Vec> &v_ptr) const {
130 Vec vec;
131 CHKERR vecCreateGhost(name, rc, &vec);
132 v_ptr.reset(vec, false);
134}
135
137 Vec xin, const std::string x_problem, const std::string x_field_name,
138 RowColData x_rc, Vec yin, const std::string y_problem,
139 const std::string y_field_name, RowColData y_rc, VecScatter *newctx) const {
140 const MoFEM::Interface &m_field = cOre;
142 std::vector<int> idx(0), idy(0);
143 CHKERR m_field.getInterface<ISManager>()
144 ->isCreateFromProblemFieldToOtherProblemField(
145 x_problem, x_field_name, x_rc, y_problem, y_field_name, y_rc, idx,
146 idy);
147 IS ix, iy;
148 CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
149 PETSC_USE_POINTER, &ix);
150 CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
151 PETSC_USE_POINTER, &iy);
152 CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
153 CHKERR ISDestroy(&ix);
154 CHKERR ISDestroy(&iy);
156}
157
159 Vec xin, const std::string x_problem, RowColData x_rc, Vec yin,
160 const std::string y_problem, RowColData y_rc, VecScatter *newctx) const {
161 const MoFEM::Interface &m_field = cOre;
163 std::vector<int> idx(0), idy(0);
164 CHKERR m_field.getInterface<ISManager>()->isCreateFromProblemToOtherProblem(
165 x_problem, x_rc, y_problem, y_rc, idx, idy);
166 IS ix, iy;
167 CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
168 PETSC_USE_POINTER, &ix);
169 CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
170 PETSC_USE_POINTER, &iy);
171 if (dEbug) {
172 ISView(ix, PETSC_VIEWER_STDOUT_WORLD);
173 ISView(iy, PETSC_VIEWER_STDOUT_WORLD);
174 }
175 CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
176 CHKERR ISDestroy(&ix);
177 CHKERR ISDestroy(&iy);
179}
180
182VecManager::vecScatterCreate(Vec xin, const std::string x_problem,
183 const std::string x_field_name, RowColData x_rc,
184 Vec yin, const std::string y_problem,
185 const std::string y_field_name, RowColData y_rc,
186 SmartPetscObj<VecScatter> &smart_newctx) const {
188 VecScatter newctx;
189 CHKERR vecScatterCreate(xin, x_problem, x_field_name, x_rc, yin, y_problem,
190 y_field_name, y_rc, &newctx);
191 smart_newctx = SmartPetscObj<VecScatter>(newctx);
193}
194
196VecManager::vecScatterCreate(Vec xin, const std::string x_problem,
197 RowColData x_rc, Vec yin,
198 const std::string y_problem, RowColData y_rc,
199 SmartPetscObj<VecScatter> &smart_newctx) const {
201 VecScatter newctx;
202 CHKERR vecScatterCreate(xin, x_problem, x_rc, yin, y_problem, y_rc, &newctx);
203 smart_newctx = SmartPetscObj<VecScatter>(newctx);
205}
206
208 RowColData rc, Vec V,
209 InsertMode mode,
210 ScatterMode scatter_mode) const {
213 DofIdx nb_local_dofs, nb_ghost_dofs;
214 switch (rc) {
215 case ROW:
216 nb_local_dofs = problem_ptr->getNbLocalDofsRow();
217 nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
218 dofs = const_cast<NumeredDofEntityByLocalIdx *>(
219 &problem_ptr->numeredRowDofsPtr->get<PetscLocalIdx_mi_tag>());
220 break;
221 case COL:
222 nb_local_dofs = problem_ptr->getNbLocalDofsCol();
223 nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
224 dofs = const_cast<NumeredDofEntityByLocalIdx *>(
225 &problem_ptr->numeredColDofsPtr->get<PetscLocalIdx_mi_tag>());
226 break;
227 default:
228 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
229 }
230 Vec Vlocal;
231 CHKERR VecGhostGetLocalForm(V, &Vlocal);
232 int size;
233 CHKERR VecGetLocalSize(V, &size);
234 if (size != nb_local_dofs) {
235 SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
236 "data inconsistency: check ghost vector, problem <%s> with nb. of "
237 "local nodes %d != %d",
238 problem_ptr->getName().c_str(), size, nb_local_dofs);
239 }
240 CHKERR VecGetLocalSize(Vlocal, &size);
241 if (size != nb_local_dofs + nb_ghost_dofs) {
242 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
243 "data inconsistency: check ghost vector, problem with nb. of "
244 "ghost nodes %d != ",
245 size, nb_local_dofs + nb_ghost_dofs);
246 }
247 NumeredDofEntityByLocalIdx::iterator miit = dofs->lower_bound(0);
248 NumeredDofEntityByLocalIdx::iterator hi_miit =
249 dofs->upper_bound(nb_local_dofs + nb_ghost_dofs);
250 DofIdx ii = 0;
251 switch (scatter_mode) {
252 case SCATTER_FORWARD: {
253 PetscScalar *array;
254 CHKERR VecGetArray(Vlocal, &array);
255 switch (mode) {
256 case INSERT_VALUES:
257 for (; miit != hi_miit; ++miit, ++ii)
258 array[ii] = (*miit)->getFieldData();
259 break;
260 case ADD_VALUES:
261 for (; miit != hi_miit; ++miit, ++ii)
262 array[ii] += (*miit)->getFieldData();
263 break;
264 default:
265 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
266 }
267 CHKERR VecRestoreArray(Vlocal, &array);
268 }; break;
269 case SCATTER_REVERSE: {
270 const PetscScalar *array;
271 CHKERR VecGetArrayRead(Vlocal, &array);
272 switch (mode) {
273 case INSERT_VALUES:
274 for (; miit != hi_miit; ++miit, ++ii) {
275 // std::cerr << *miit << std::endl;
276 // std::cerr << array[ii] << std::endl;
277 (*miit)->getFieldData() = array[ii];
278 }
279 break;
280 case ADD_VALUES:
281 for (; miit != hi_miit; ++miit, ++ii)
282 (*miit)->getFieldData() += array[ii];
283 break;
284 default:
285 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
286 }
287 CHKERR VecRestoreArrayRead(Vlocal, &array);
288 }; break;
289 default:
290 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
291 }
292 CHKERR VecGhostRestoreLocalForm(V, &Vlocal);
294}
295
297 RowColData rc, Vec V,
298 InsertMode mode,
299 ScatterMode scatter_mode) const {
300 const MoFEM::Interface &m_field = cOre;
301 const Problem *problem_ptr;
303 CHKERR m_field.get_problem(name, &problem_ptr);
304 CHKERR setLocalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
306}
307
310 Vec V, InsertMode mode,
311 ScatterMode scatter_mode) const {
314 DofsByGlobalIdx;
315 DofsByGlobalIdx *dofs;
316 DofIdx nb_dofs;
317 switch (rc) {
318 case ROW:
319 nb_dofs = problem_ptr->getNbDofsRow();
320 dofs = const_cast<DofsByGlobalIdx *>(
321 &problem_ptr->numeredRowDofsPtr->get<PetscGlobalIdx_mi_tag>());
322 break;
323 case COL:
324 nb_dofs = problem_ptr->getNbDofsCol();
325 dofs = const_cast<DofsByGlobalIdx *>(
326 &problem_ptr->numeredColDofsPtr->get<PetscGlobalIdx_mi_tag>());
327 break;
328 default:
329 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
330 "Function works only for ROWs and COLs");
331 }
332 DofsByGlobalIdx::iterator miit = dofs->lower_bound(0);
333 DofsByGlobalIdx::iterator hi_miit = dofs->upper_bound(nb_dofs);
334 auto comm = PetscObjectComm((PetscObject)V);
335 switch (scatter_mode) {
336 case SCATTER_REVERSE: {
337 VecScatter ctx;
338 Vec V_glob;
339 CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
340 CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
341 CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
342 int size;
343 CHKERR VecGetSize(V_glob, &size);
344 if (size != nb_dofs) {
345 SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
346 "Size of vector is inconsistent with size of problem. You could "
347 "use wrong vector with wrong problem, or you created vector "
348 "before you remove DOFs from problem.");
349 }
350 PetscScalar *array;
351 CHKERR VecGetArray(V_glob, &array);
352 switch (mode) {
353 case INSERT_VALUES:
354 for (; miit != hi_miit; miit++)
355 (*miit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
356 break;
357 case ADD_VALUES:
358 for (; miit != hi_miit; miit++)
359 (*miit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
360 break;
361 default:
362 SETERRQ(comm, MOFEM_NOT_IMPLEMENTED, "not implemented");
363 }
364 CHKERR VecRestoreArray(V_glob, &array);
365 CHKERR VecScatterDestroy(&ctx);
366 CHKERR VecDestroy(&V_glob);
367 break;
368 }
369 default:
370 SETERRQ(comm, MOFEM_NOT_IMPLEMENTED, "not implemented");
371 }
373}
374
376VecManager::setGlobalGhostVector(const std::string name, RowColData rc, Vec V,
377 InsertMode mode,
378 ScatterMode scatter_mode) const {
379 const MoFEM::Interface &m_field = cOre;
380 const Problem *problem_ptr;
382 CHKERR m_field.get_problem(name, &problem_ptr);
383 CHKERR setGlobalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
385}
386
387template <int MODE> struct SetOtherLocalGhostVector {
388 template <typename A0, typename A1, typename A2, typename A3, typename A4>
389 inline MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit,
390 A4 &cpy_bit_number,
391 const std::string cpy_field_name) {
393 for (; miit != hi_miit; miit++) {
394 if (miit->get()->getHasLocalIndex()) {
395 const auto uid = DofEntity::getUniqueIdCalculate(
396 (*miit)->getEntDofIdx(), FieldEntity::getLocalUniqueIdCalculate(
397 cpy_bit_number, (*miit)->getEnt()));
398 auto diiiit = dofs_ptr->template get<Unique_mi_tag>().find(uid);
399 if (diiiit == dofs_ptr->template get<Unique_mi_tag>().end()) {
401 MOFEM_LOG_ATTRIBUTES("VECSELF",
403 MOFEM_LOG("VECSELF", Sev::error)
404 << "Problem finding DOFs in the copy field";
405 MOFEM_LOG("VECSELF", Sev::error) << "Field DOF: "<< (**miit);
406 MOFEM_LOG("VECSELF", Sev::error)
407 << "Copy field name: " << cpy_field_name;
408 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
409 "Automatic creation of entity and dof not implemented");
410 }
411 if constexpr (MODE == INSERT_VALUES)
412 (*diiiit)->getFieldData() = array[(*miit)->getPetscLocalDofIdx()];
413 else if constexpr (MODE == ADD_VALUES)
414 (*diiiit)->getFieldData() += array[(*miit)->getPetscLocalDofIdx()];
415 }
416 }
418 }
419};
420
422 const Problem *problem_ptr, const std::string field_name,
423 const std::string cpy_field_name, RowColData rc, Vec V, InsertMode mode,
424 ScatterMode scatter_mode) const {
425 const MoFEM::Interface &m_field = cOre;
426 auto fields_ptr = m_field.get_fields();
427 auto dofs_ptr = m_field.get_dofs();
430 DofsByUId *dofs;
431
432 switch (rc) {
433 case ROW:
434 dofs = const_cast<DofsByUId *>(
435 &problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>());
436 break;
437 case COL:
438 dofs = const_cast<DofsByUId *>(
439 &problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>());
440 break;
441 default:
442 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
443 }
444
445 auto cpy_fit = fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
446 if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end())
447 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
448 "cpy field < %s > not found, (top tip: check spelling)",
449 cpy_field_name.c_str());
450 const auto cpy_bit_number = (*cpy_fit)->getBitNumber();
451
452 auto miit = dofs->lower_bound(
454 if (miit == dofs->end()) {
455 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
456 "cpy field < %s > not found, (top tip: check spelling)",
457 field_name.c_str());
458 }
459 auto hi_miit = dofs->upper_bound(
461
462 if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
463 SETERRQ4(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
464 "fields have to have same space (%s) %s != (%s) %s",
465 (*miit)->getName().c_str(), FieldSpaceNames[(*miit)->getSpace()],
466 cpy_field_name.c_str(), FieldSpaceNames[(*cpy_fit)->getSpace()]);
467 }
468 if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
469 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
470 "fields have to have same rank");
471 }
472
473 switch (scatter_mode) {
474 case SCATTER_REVERSE: {
475
476 PetscScalar *array;
477 CHKERR VecGetArray(V, &array);
478 if (mode == INSERT_VALUES)
480 dofs_ptr, array, miit, hi_miit, cpy_bit_number, cpy_field_name);
481 else if (mode == ADD_VALUES)
483 dofs_ptr, array, miit, hi_miit, cpy_bit_number, cpy_field_name);
484 else
485 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
486 "Operation mode not implemented");
487 CHKERR VecRestoreArray(V, &array);
488 } break;
489 case SCATTER_FORWARD: {
490 for (; miit != hi_miit; miit++) {
491 if (!miit->get()->getHasLocalIndex())
492 continue;
493 const auto uid = DofEntity::getUniqueIdCalculate(
494 (*miit)->getEntDofIdx(), FieldEntity::getLocalUniqueIdCalculate(
495 cpy_bit_number, (*miit)->getEnt()));
496 auto diiiit = dofs_ptr->get<Unique_mi_tag>().find(uid);
497 if (diiiit == dofs_ptr->get<Unique_mi_tag>().end()) {
498 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
499 "no data to fill the vector (top tip: you want scatter forward "
500 "or scatter reverse?)");
501 }
502 CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
503 (*diiiit)->getFieldData(), mode);
504 }
505 CHKERR VecAssemblyBegin(V);
506 CHKERR VecAssemblyEnd(V);
507 } break;
508 default:
509 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
510 }
512}
513
515 const std::string name, const std::string field_name,
516 const std::string cpy_field_name, RowColData rc, Vec V, InsertMode mode,
517 ScatterMode scatter_mode) const {
518 const MoFEM::Interface &m_field = cOre;
519 const Problem *problem_ptr;
521 CHKERR m_field.get_problem(name, &problem_ptr);
522 CHKERR setOtherLocalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
523 V, mode, scatter_mode);
525}
526
527template <int MODE> struct SetOtherGlobalGhostVector {
528 template <typename A0, typename A1, typename A2, typename A3, typename A4>
529 inline MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit,
530 A4 &cpy_bit_number) {
532 for (; miit != hi_miit; miit++) {
533 const auto uid = DofEntity::getUniqueIdCalculate(
534 (*miit)->getEntDofIdx(), FieldEntity::getLocalUniqueIdCalculate(
535 cpy_bit_number, (*miit)->getEnt()));
536 auto diiiit = dofs_ptr->template get<Unique_mi_tag>().find(uid);
537#ifndef NDEBUG
538 if (diiiit == dofs_ptr->template get<Unique_mi_tag>().end()) {
539 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
540 "Automatic creation of entity and dof not implemented");
541 }
542#endif
543 if constexpr (MODE == INSERT_VALUES)
544 (*diiiit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
545 else if constexpr (MODE == ADD_VALUES)
546 (*diiiit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
547 }
549 }
550};
551
553 const Problem *problem_ptr, const std::string field_name,
554 const std::string cpy_field_name, RowColData rc, Vec V, InsertMode mode,
555 ScatterMode scatter_mode) const {
556 const MoFEM::Interface &m_field = cOre;
557 auto fields_ptr = m_field.get_fields();
558 auto dofs_ptr = m_field.get_dofs();
559 auto *field_ents = m_field.get_field_ents();
562 DofIdx nb_dofs;
563 switch (rc) {
564 case ROW:
565 nb_dofs = problem_ptr->getNbDofsRow();
566 dofs = const_cast<NumeredDofEntityByUId *>(
567 &problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>());
568 break;
569 case COL:
570 nb_dofs = problem_ptr->getNbDofsCol();
571 dofs = const_cast<NumeredDofEntityByUId *>(
572 &problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>());
573 break;
574 default:
575 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
576 }
577 auto cpy_fit = fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
578 if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
579 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
580 "cpy field < %s > not found, (top tip: check spelling)",
581 cpy_field_name.c_str());
582 }
583 const auto cpy_bit_number = (*cpy_fit)->getBitNumber();
584
585 auto miit = dofs->lower_bound(
587 if (miit == dofs->end()) {
588 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
589 "problem field < %s > not found, (top tip: check spelling)",
590 field_name.c_str());
591 }
592 auto hi_miit = dofs->upper_bound(
594 if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
595 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
596 "fields have to have same space");
597 }
598 if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
599 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
600 "fields have to have same rank");
601 }
602 switch (scatter_mode) {
603 case SCATTER_REVERSE: {
604 Vec V_glob;
605 VecScatter ctx;
606 CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
607 CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
608 CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
609 int size;
610 CHKERR VecGetSize(V_glob, &size);
611 if (size != nb_dofs)
612 SETERRQ(
613 PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
614 "data inconsistency: nb. of dofs and declared nb. dofs in database");
615 PetscScalar *array;
616 CHKERR VecGetArray(V_glob, &array);
617
618 if (mode == INSERT_VALUES)
620 dofs_ptr, array, miit, hi_miit, cpy_bit_number);
621 else if (mode == ADD_VALUES)
622 CHKERR SetOtherGlobalGhostVector<ADD_VALUES>()(dofs_ptr, array, miit,
623 hi_miit, cpy_bit_number);
624 else
625 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
626 "Operation mode not implemented");
627
628 CHKERR VecRestoreArray(V_glob, &array);
629 CHKERR VecDestroy(&V_glob);
630 CHKERR VecScatterDestroy(&ctx);
631 } break;
632 case SCATTER_FORWARD: {
633 for (; miit != hi_miit; miit++) {
634 const auto uid = DofEntity::getUniqueIdCalculate(
635 (*miit)->getEntDofIdx(), FieldEntity::getLocalUniqueIdCalculate(
636 cpy_bit_number, (*miit)->getEnt()));
637 auto diiiit = dofs_ptr->get<Unique_mi_tag>().find(uid);
638 if (diiiit == dofs_ptr->get<Unique_mi_tag>().end()) {
639 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
640 "No data to fill the vector (top tip: you want scatter forward "
641 "or scatter reverse?)");
642 }
643 CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
644 (*diiiit)->getFieldData(), mode);
645 }
646 CHKERR VecAssemblyBegin(V);
647 CHKERR VecAssemblyEnd(V);
648 } break;
649 default:
650 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
651 }
653}
654
656 const std::string name, const std::string field_name,
657 const std::string cpy_field_name, RowColData rc, Vec V, InsertMode mode,
658 ScatterMode scatter_mode) const {
659 const MoFEM::Interface &m_field = cOre;
660 const Problem *problem_ptr;
662 CHKERR m_field.get_problem(name, &problem_ptr);
663 CHKERR setOtherGlobalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
664 V, mode, scatter_mode);
666}
667
668} // namespace MoFEM
RowColData
RowColData.
Definition: definitions.h:136
@ COL
Definition: definitions.h:136
@ ROW
Definition: definitions.h:136
static const char *const FieldSpaceNames[]
Definition: definitions.h:105
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_NOT_FOUND
Definition: definitions.h:46
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
NumeredDofEntity_multiIndex::index< Unique_mi_tag >::type NumeredDofEntityByUId
Numbered DoF multi-index by UId.
NumeredDofEntity_multiIndex::index< PetscLocalIdx_mi_tag >::type NumeredDofEntityByLocalIdx
Numbered DoF multi-index by local index.
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
virtual const FieldEntity_multiIndex * get_field_ents() const =0
Get the field ents object.
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:342
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
Definition: LogManager.hpp:299
#define MOFEM_LOG_FUNCTION()
Set scope.
Definition: LogManager.hpp:328
MoFEMErrorCode vecCreateGhost(const std::string name, RowColData rc, Vec *V) const
create ghost vector for problem (collective)
Definition: VecManager.cpp:74
MoFEMErrorCode setOtherLocalGhostVector(const Problem *problem_ptr, const std::string field_name, const std::string cpy_field_name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
Copy vector to field which is not part of the problem.
Definition: VecManager.cpp:421
MoFEMErrorCode setGlobalGhostVector(const Problem *problem_ptr, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
set values of vector from/to mesh database (collective)
Definition: VecManager.cpp:309
MoFEMErrorCode vecScatterCreate(Vec xin, const std::string x_problem, const std::string x_field_name, RowColData x_rc, Vec yin, const std::string y_problem, const std::string y_field_name, RowColData y_rc, VecScatter *newctx) const
create scatter for vectors form one to another problem (collective)
Definition: VecManager.cpp:136
MoFEMErrorCode vecCreateSeq(const std::string name, RowColData rc, Vec *V) const
create local vector for problem
Definition: VecManager.cpp:51
MoFEMErrorCode setLocalGhostVector(const Problem *problem_ptr, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
set values of vector from/to meshdatabase
Definition: VecManager.cpp:207
MoFEMErrorCode setOtherGlobalGhostVector(const Problem *problem_ptr, const std::string field_name, const std::string cpy_field_name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
Copy vector to field which is not part of the problem (collective)
Definition: VecManager.cpp:552
const FTensor::Tensor2< T, Dim, Dim > Vec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
int DofIdx
Index of DOF.
Definition: Types.hpp:29
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
constexpr auto field_name
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
virtual MPI_Comm & get_comm() const =0
Core (interface) class.
Definition: Core.hpp:92
Deprecated interface functions.
static UId getUniqueIdCalculate(const DofIdx dof, UId ent_uid)
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
MultiIndex Tag for field name.
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:33
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:327
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:374
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:319
keeps basic data about problem
DofIdx getNbLocalDofsRow() const
DofIdx getNbDofsRow() const
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem
DofIdx getNbGhostDofsCol() const
DofIdx getNbDofsCol() const
DofIdx getNbLocalDofsCol() const
DofIdx getNbGhostDofsRow() const
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredColDofsPtr
store DOFs on columns for this problem
MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit, A4 &cpy_bit_number)
Definition: VecManager.cpp:529
MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit, A4 &cpy_bit_number, const std::string cpy_field_name)
Definition: VecManager.cpp:389
base class for all interface classes
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:33
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: VecManager.cpp:21
const MoFEM::Interface & cOre
Definition: VecManager.hpp:38
VecManager(const MoFEM::Core &core)
Definition: VecManager.cpp:27