v0.8.23
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 
18 namespace MoFEM {
19 
21  UnknownInterface **iface) const {
23  *iface = NULL;
24  if (uuid == IDD_MOFEMVEC) {
25  *iface = const_cast<VecManager *>(this);
27  }
28  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
30 }
31 
33  : cOre(const_cast<MoFEM::Core &>(core)), dEbug(false) {}
35 
37  Vec *V) const {
38  const MoFEM::Interface &m_field = cOre;
39  const Problem *problem_ptr;
41  CHKERR m_field.get_problem(name, &problem_ptr);
42  DofIdx nb_local_dofs, nb_ghost_dofs;
43  switch (rc) {
44  case ROW:
45  nb_local_dofs = problem_ptr->getNbLocalDofsRow();
46  nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
47  break;
48  case COL:
49  nb_local_dofs = problem_ptr->getNbLocalDofsCol();
50  nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
51  break;
52  default:
53  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
54  }
55  CHKERR ::VecCreateSeq(PETSC_COMM_SELF, nb_local_dofs + nb_ghost_dofs, V);
57 }
58 
60  RowColData rc, Vec *V) const {
61  const MoFEM::Interface &m_field = cOre;
62  const Problem *problem_ptr;
64  CHKERR m_field.get_problem(name, &problem_ptr);
65  DofIdx nb_dofs, nb_local_dofs, nb_ghost_dofs;
67  switch (rc) {
68  case ROW:
69  nb_dofs = problem_ptr->getNbDofsRow();
70  nb_local_dofs = problem_ptr->getNbLocalDofsRow();
71  nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
72  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
73  &problem_ptr->numeredDofsRows->get<PetscLocalIdx_mi_tag>());
74  break;
75  case COL:
76  nb_dofs = problem_ptr->getNbDofsCol();
77  nb_local_dofs = problem_ptr->getNbLocalDofsCol();
78  nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
79  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
80  &problem_ptr->numeredDofsCols->get<PetscLocalIdx_mi_tag>());
81  break;
82  default:
83  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
84  }
85  // get ghost dofs
86  auto miit = dofs->lower_bound(nb_local_dofs);
87  auto hi_miit = dofs->upper_bound(nb_local_dofs + nb_ghost_dofs);
88  int count = std::distance(miit, hi_miit);
89  if (count != nb_ghost_dofs) {
90  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
91  }
92  std::vector<DofIdx> ghost_idx(count);
93  for (auto vit = ghost_idx.begin(); miit != hi_miit; ++miit, ++vit) {
94  *vit = (*miit)->petscGloablDofIdx;
95  }
96  CHKERR ::VecCreateGhost(PETSC_COMM_WORLD, nb_local_dofs, nb_dofs,
97  nb_ghost_dofs, &ghost_idx[0], V);
99 }
100 
102  RowColData rc,
103  SmartPetscObj<Vec> &v_ptr) const {
105  Vec vec;
106  CHKERR vecCreateGhost(name, rc, &vec);
107  v_ptr.reset(vec, false);
109 }
110 
112 VecManager::vecScatterCreate(Vec xin, const std::string &x_problem,
113  const std::string &x_field_name, RowColData x_rc,
114  Vec yin, const std::string &y_problem,
115  const std::string &y_field_name, RowColData y_rc,
116  VecScatter *newctx) const {
117  const MoFEM::Interface &m_field = cOre;
119  std::vector<int> idx(0), idy(0);
120  CHKERR m_field.getInterface<ISManager>()
121  ->isCreateFromProblemFieldToOtherProblemField(
122  x_problem, x_field_name, x_rc, y_problem, y_field_name, y_rc, idx,
123  idy);
124  IS ix, iy;
125  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
126  PETSC_USE_POINTER, &ix);
127  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
128  PETSC_USE_POINTER, &iy);
129  CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
130  CHKERR ISDestroy(&ix);
131  CHKERR ISDestroy(&iy);
133 }
134 
136  Vec xin, const std::string &x_problem, RowColData x_rc, Vec yin,
137  const std::string &y_problem, RowColData y_rc, VecScatter *newctx) const {
138  const MoFEM::Interface &m_field = cOre;
140  std::vector<int> idx(0), idy(0);
141  CHKERR m_field.getInterface<ISManager>()->isCreateFromProblemToOtherProblem(
142  x_problem, x_rc, y_problem, y_rc, idx, idy);
143  IS ix, iy;
144  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
145  PETSC_USE_POINTER, &ix);
146  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
147  PETSC_USE_POINTER, &iy);
148  if (dEbug) {
149  ISView(ix, PETSC_VIEWER_STDOUT_WORLD);
150  ISView(iy, PETSC_VIEWER_STDOUT_WORLD);
151  }
152  CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
153  CHKERR ISDestroy(&ix);
154  CHKERR ISDestroy(&iy);
156 }
157 
159  RowColData rc, Vec V,
160  InsertMode mode,
161  ScatterMode scatter_mode) const {
164  DofIdx nb_local_dofs, nb_ghost_dofs;
165  switch (rc) {
166  case ROW:
167  nb_local_dofs = problem_ptr->getNbLocalDofsRow();
168  nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
169  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
170  &problem_ptr->numeredDofsRows->get<PetscLocalIdx_mi_tag>());
171  break;
172  case COL:
173  nb_local_dofs = problem_ptr->getNbLocalDofsCol();
174  nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
175  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
176  &problem_ptr->numeredDofsCols->get<PetscLocalIdx_mi_tag>());
177  break;
178  default:
179  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
180  }
181  Vec Vlocal;
182  CHKERR VecGhostGetLocalForm(V, &Vlocal);
183  int size;
184  CHKERR VecGetLocalSize(V, &size);
185  if (size != nb_local_dofs) {
186  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
187  "data inconsistency: check ghost vector, problem with nb. of "
188  "local nodes %d != %d",
189  size, nb_local_dofs);
190  }
191  CHKERR VecGetLocalSize(Vlocal, &size);
192  if (size != nb_local_dofs + nb_ghost_dofs) {
193  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
194  "data inconsistency: check ghost vector, problem with nb. of "
195  "ghost nodes %d != ",
196  size, nb_local_dofs + nb_ghost_dofs);
197  }
198  NumeredDofEntityByLocalIdx::iterator miit = dofs->lower_bound(0);
199  NumeredDofEntityByLocalIdx::iterator hi_miit =
200  dofs->upper_bound(nb_local_dofs + nb_ghost_dofs);
201  DofIdx ii = 0;
202  switch (scatter_mode) {
203  case SCATTER_FORWARD: {
204  PetscScalar *array;
205  CHKERR VecGetArray(Vlocal, &array);
206  switch (mode) {
207  case INSERT_VALUES:
208  for (; miit != hi_miit; ++miit, ++ii)
209  array[ii] = (*miit)->getFieldData();
210  break;
211  case ADD_VALUES:
212  for (; miit != hi_miit; ++miit, ++ii)
213  array[ii] += (*miit)->getFieldData();
214  break;
215  default:
216  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
217  }
218  CHKERR VecRestoreArray(Vlocal, &array);
219  }; break;
220  case SCATTER_REVERSE: {
221  const PetscScalar *array;
222  CHKERR VecGetArrayRead(Vlocal, &array);
223  switch (mode) {
224  case INSERT_VALUES:
225  for (; miit != hi_miit; ++miit, ++ii) {
226  // std::cerr << *miit << std::endl;
227  // std::cerr << array[ii] << std::endl;
228  (*miit)->getFieldData() = array[ii];
229  }
230  break;
231  case ADD_VALUES:
232  for (; miit != hi_miit; ++miit, ++ii)
233  (*miit)->getFieldData() += array[ii];
234  break;
235  default:
236  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
237  }
238  CHKERR VecRestoreArrayRead(Vlocal, &array);
239  }; break;
240  default:
241  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
242  }
243  CHKERR VecGhostRestoreLocalForm(V, &Vlocal);
245 }
246 
248  RowColData rc, Vec V,
249  InsertMode mode,
250  ScatterMode scatter_mode) const {
251  const MoFEM::Interface &m_field = cOre;
252  const Problem *problem_ptr;
254  CHKERR m_field.get_problem(name, &problem_ptr);
255  CHKERR setLocalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
257 }
258 
261  Vec V, InsertMode mode,
262  ScatterMode scatter_mode) const {
264  typedef NumeredDofEntity_multiIndex::index<PetscGlobalIdx_mi_tag>::type
265  DofsByGlobalIdx;
266  DofsByGlobalIdx *dofs;
267  DofIdx nb_dofs;
268  switch (rc) {
269  case ROW:
270  nb_dofs = problem_ptr->getNbDofsRow();
271  dofs = const_cast<DofsByGlobalIdx *>(
272  &problem_ptr->numeredDofsRows->get<PetscGlobalIdx_mi_tag>());
273  break;
274  case COL:
275  nb_dofs = problem_ptr->getNbDofsCol();
276  dofs = const_cast<DofsByGlobalIdx *>(
277  &problem_ptr->numeredDofsCols->get<PetscGlobalIdx_mi_tag>());
278  break;
279  default:
280  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
281  }
282  DofsByGlobalIdx::iterator miit = dofs->lower_bound(0);
283  DofsByGlobalIdx::iterator hi_miit = dofs->upper_bound(nb_dofs);
284  switch (scatter_mode) {
285  case SCATTER_REVERSE: {
286  VecScatter ctx;
287  Vec V_glob;
288  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
289  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
290  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
291  int size;
292  CHKERR VecGetSize(V_glob, &size);
293  if (size != nb_dofs) {
294  SETERRQ(
295  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
296  "data inconsistency: nb. of dofs and declared nb. dofs in database");
297  }
298  PetscScalar *array;
299  CHKERR VecGetArray(V_glob, &array);
300  switch (mode) {
301  case INSERT_VALUES:
302  for (; miit != hi_miit; miit++)
303  (*miit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
304  break;
305  case ADD_VALUES:
306  for (; miit != hi_miit; miit++)
307  (*miit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
308  break;
309  default:
310  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
311  }
312  CHKERR VecRestoreArray(V_glob, &array);
313  CHKERR VecScatterDestroy(&ctx);
314  CHKERR VecDestroy(&V_glob);
315  break;
316  }
317  default:
318  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
319  }
321 }
322 
324 VecManager::setGlobalGhostVector(const std::string &name, RowColData rc, Vec V,
325  InsertMode mode,
326  ScatterMode scatter_mode) const {
327  const MoFEM::Interface &m_field = cOre;
328  const Problem *problem_ptr;
330  CHKERR m_field.get_problem(name, &problem_ptr);
331  CHKERR setGlobalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
333 }
334 
335 template <int MODE> struct SetOtherLocalGhostVector {
336  template <typename A0, typename A1, typename A2, typename A3, typename A4>
337  inline MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit,
338  A4 &cpy_field_name) {
340  for (; miit != hi_miit; miit++) {
341  if (miit->get()->getHasLocalIndex()) {
342  auto diiiit =
343  dofs_ptr
344  ->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
345  .find(boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
346  (*miit)->getEntDofIdx()));
347  if (diiiit ==
348  dofs_ptr
349  ->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
350  .end()) {
351  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
352  "Automatic creation of entity and dof not implemented");
353  }
354  if (MODE == INSERT_VALUES)
355  (*diiiit)->getFieldData() = array[(*miit)->getPetscLocalDofIdx()];
356  else if (MODE == ADD_VALUES)
357  (*diiiit)->getFieldData() += array[(*miit)->getPetscLocalDofIdx()];
358  }
359  }
361  }
362 };
363 
365  const Problem *problem_ptr, const std::string &field_name,
366  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
367  ScatterMode scatter_mode) const {
368  const MoFEM::Interface &m_field = cOre;
369  const Field_multiIndex *fields_ptr;
370  const DofEntity_multiIndex *dofs_ptr;
372  CHKERR m_field.get_fields(&fields_ptr);
373  CHKERR m_field.get_dofs(&dofs_ptr);
374  typedef NumeredDofEntity_multiIndex::index<FieldName_mi_tag>::type DofsByName;
375  DofsByName *dofs;
376  switch (rc) {
377  case ROW:
378  dofs = const_cast<DofsByName *>(
379  &problem_ptr->numeredDofsRows->get<FieldName_mi_tag>());
380  break;
381  case COL:
382  dofs = const_cast<DofsByName *>(
383  &problem_ptr->numeredDofsCols->get<FieldName_mi_tag>());
384  break;
385  default:
386  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
387  }
388  Field_multiIndex::index<FieldName_mi_tag>::type::iterator cpy_fit =
389  fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
390  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
391  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
392  "cpy field < %s > not found, (top tip: check spelling)",
393  cpy_field_name.c_str());
394  }
395  DofsByName::iterator miit = dofs->lower_bound(field_name);
396  if (miit == dofs->end()) {
397  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
398  "cpy field < %s > not found, (top tip: check spelling)",
399  field_name.c_str());
400  }
401  DofsByName::iterator hi_miit = dofs->upper_bound(field_name);
402  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
403  SETERRQ4(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
404  "fields have to have same space (%s) %s != (%s) %s",
405  (*miit)->getName().c_str(), FieldSpaceNames[(*miit)->getSpace()],
406  cpy_field_name.c_str(), FieldSpaceNames[(*cpy_fit)->getSpace()]);
407  }
408  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
409  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
410  "fields have to have same rank");
411  }
412 
413  switch (scatter_mode) {
414  case SCATTER_REVERSE: {
415 
416  PetscScalar *array;
417  CHKERR VecGetArray(V, &array);
418  if (mode == INSERT_VALUES)
419  CHKERR SetOtherLocalGhostVector<INSERT_VALUES>()(dofs_ptr, array, miit,
420  hi_miit, cpy_field_name);
421  else if (mode == ADD_VALUES)
422  CHKERR SetOtherLocalGhostVector<ADD_VALUES>()(dofs_ptr, array, miit,
423  hi_miit, cpy_field_name);
424  else
425  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
426  "Operation mode not implemented");
427  CHKERR VecRestoreArray(V, &array);
428  } break;
429  case SCATTER_FORWARD: {
430  for (; miit != hi_miit; miit++) {
431  if (!miit->get()->getHasLocalIndex())
432  continue;
433  DofEntity_multiIndex::index<
434  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator diiiit;
435  diiiit =
436  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
437  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
438  (*miit)->getEntDofIdx()));
439  if (diiiit ==
440  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
441  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
442  "no data to fill the vector (top tip: you want scatter forward "
443  "or scatter reverse?)");
444  }
445  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
446  (*diiiit)->getFieldData(), mode);
447  }
448  CHKERR VecAssemblyBegin(V);
449  CHKERR VecAssemblyEnd(V);
450  } break;
451  default:
452  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
453  }
455 }
456 
458  const std::string &name, const std::string &field_name,
459  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
460  ScatterMode scatter_mode) const {
461  const MoFEM::Interface &m_field = cOre;
462  const Problem *problem_ptr;
464  CHKERR m_field.get_problem(name, &problem_ptr);
465  CHKERR setOtherLocalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
466  V, mode, scatter_mode);
468 }
469 
470 template <int MODE> struct SetOtherGlobalGhostVector {
471  template <typename A0, typename A1, typename A2, typename A3, typename A4>
472  inline MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit,
473  A4 &cpy_field_name) {
475  for (; miit != hi_miit; miit++) {
476  auto diiiit =
477  dofs_ptr->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
478  .find(boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
479  (*miit)->getEntDofIdx()));
480  if (diiiit ==
481  dofs_ptr->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
482  .end()) {
483  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
484  "Automatic creation of entity and dof not implemented");
485  }
486  if (MODE == INSERT_VALUES)
487  (*diiiit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
488  else if (MODE == ADD_VALUES)
489  (*diiiit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
490  }
492  }
493 };
494 
496  const Problem *problem_ptr, const std::string &field_name,
497  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
498  ScatterMode scatter_mode) const {
499  const MoFEM::Interface &m_field = cOre;
500  const Field_multiIndex *fields_ptr;
501  const DofEntity_multiIndex *dofs_ptr;
502  const FieldEntity_multiIndex *field_ents;
504  CHKERR m_field.get_fields(&fields_ptr);
505  CHKERR m_field.get_dofs(&dofs_ptr);
506  CHKERR m_field.get_field_ents(&field_ents);
507  typedef NumeredDofEntityByFieldName DofsByName;
508  DofsByName *dofs;
509  DofIdx nb_dofs;
510  switch (rc) {
511  case ROW:
512  nb_dofs = problem_ptr->getNbDofsRow();
513  dofs = const_cast<DofsByName *>(
514  &problem_ptr->numeredDofsRows->get<FieldName_mi_tag>());
515  break;
516  case COL:
517  nb_dofs = problem_ptr->getNbDofsCol();
518  dofs = const_cast<DofsByName *>(
519  &problem_ptr->numeredDofsCols->get<FieldName_mi_tag>());
520  break;
521  default:
522  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
523  }
524  auto cpy_fit = fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
525  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
526  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
527  "cpy field < %s > not found, (top tip: check spelling)",
528  cpy_field_name.c_str());
529  }
530  auto miit = dofs->lower_bound(field_name);
531  if (miit == dofs->end()) {
532  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
533  "problem field < %s > not found, (top tip: check spelling)",
534  field_name.c_str());
535  }
536  auto hi_miit = dofs->upper_bound(field_name);
537  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
538  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
539  "fields have to have same space");
540  }
541  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
542  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
543  "fields have to have same rank");
544  }
545  switch (scatter_mode) {
546  case SCATTER_REVERSE: {
547  Vec V_glob;
548  VecScatter ctx;
549  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
550  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
551  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
552  int size;
553  CHKERR VecGetSize(V_glob, &size);
554  if (size != nb_dofs)
555  SETERRQ(
556  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
557  "data inconsistency: nb. of dofs and declared nb. dofs in database");
558  PetscScalar *array;
559  CHKERR VecGetArray(V_glob, &array);
560 
561  if (mode == INSERT_VALUES)
563  dofs_ptr, array, miit, hi_miit, cpy_field_name);
564  else if (mode == ADD_VALUES)
565  CHKERR SetOtherGlobalGhostVector<ADD_VALUES>()(dofs_ptr, array, miit,
566  hi_miit, cpy_field_name);
567  else
568  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
569  "Operation mode not implemented");
570 
571  CHKERR VecRestoreArray(V_glob, &array);
572  CHKERR VecDestroy(&V_glob);
573  CHKERR VecScatterDestroy(&ctx);
574  } break;
575  case SCATTER_FORWARD: {
576  for (; miit != hi_miit; miit++) {
577  auto diiiit =
578  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
579  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
580  (*miit)->getEntDofIdx()));
581  if (diiiit ==
582  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
583  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
584  "No data to fill the vector (top tip: you want scatter forward "
585  "or scatter reverse?)");
586  }
587  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
588  (*diiiit)->getFieldData(), mode);
589  }
590  CHKERR VecAssemblyBegin(V);
591  CHKERR VecAssemblyEnd(V);
592  } break;
593  default:
594  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
595  }
597 }
598 
600  const std::string &name, const std::string &field_name,
601  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
602  ScatterMode scatter_mode) const {
603  const MoFEM::Interface &m_field = cOre;
604  const Problem *problem_ptr;
606  CHKERR m_field.get_problem(name, &problem_ptr);
607  CHKERR setOtherGlobalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
608  V, mode, scatter_mode);
610 }
611 
612 } // namespace MoFEM
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:364
NumeredDofEntity_multiIndex::index< FieldName_mi_tag >::type NumeredDofEntityByFieldName
Numbered DoF multi-index by field name.
DofIdx getNbDofsRow() const
MoFEM interface unique ID.
DofIdx getNbGhostDofsCol() const
Section manager is used to create indexes and sectionsFIXME: ISManager is not properly testsed by ato...
Definition: ISManager.hpp:36
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
base class for all interface classes
multi_index_container< boost::shared_ptr< FieldEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< FieldEntity, UId, &FieldEntity::globalUId > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< FieldEntity::interface_type_Field, boost::string_ref, &FieldEntity::getNameRef > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntity, EntityHandle, &FieldEntity::getEnt > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< FieldEntity, const_mem_fun< FieldEntity::interface_type_Field, boost::string_ref, &FieldEntity::getNameRef >, const_mem_fun< FieldEntity, EntityHandle, &FieldEntity::getEnt > > > > > FieldEntity_multiIndex
MultiIndex container keeps FieldEntity.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
multi_index_container< boost::shared_ptr< Field >, indexed_by< hashed_unique< tag< BitFieldId_mi_tag >, const_mem_fun< Field, const BitFieldId &, &Field::getId >, HashBit< BitFieldId >, EqBit< BitFieldId > >, ordered_unique< tag< Meshset_mi_tag >, member< Field, EntityHandle, &Field::meshSet > >, ordered_unique< tag< FieldName_mi_tag >, const_mem_fun< Field, boost::string_ref, &Field::getNameRef > >, ordered_non_unique< tag< BitFieldId_space_mi_tag >, const_mem_fun< Field, FieldSpace, &Field::getSpace > > > > Field_multiIndex
Field_multiIndex for Field.
DofIdx getNbLocalDofsRow() const
static const MOFEMuuid IDD_MOFEMVEC
Definition: VecManager.hpp:26
Core (interface) class.
Definition: Core.hpp:50
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
NumeredDofEntity_multiIndex::index< PetscLocalIdx_mi_tag >::type NumeredDofEntityByLocalIdx
Numbered DoF multi-index by local index.
keeps basic data about problemThis is low level structure with information about problem,...
DofIdx getNbGhostDofsRow() const
RowColData
RowColData.
Definition: definitions.h:186
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: VecManager.cpp:20
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
const MoFEM::Interface & cOre
Definition: VecManager.hpp:41
virtual MoFEMErrorCode get_dofs(const DofEntity_multiIndex **dofs_ptr) const =0
Get dofs multi index.
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:158
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)collective - need tu be run on all processors ...
Definition: VecManager.cpp:260
MoFEMErrorCode vecCreateGhost(const std::string &name, RowColData rc, Vec *V) const
create ghost vector for problem (collective)collective - need to be run on all processors in communic...
Definition: VecManager.cpp:59
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)User specify what name of field on...
Definition: VecManager.cpp:112
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< DofEntity, UId, &DofEntity::globalUId > >, ordered_non_unique< tag< Composite_Name_And_Ent_And_EntDofIdx_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt >, const_mem_fun< DofEntity, DofIdx, &DofEntity::getEntDofIdx > > >, ordered_non_unique< tag< Unique_Ent_mi_tag >, const_mem_fun< DofEntity, const UId &, &DofEntity::getEntGlobalUniqueId > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity::interface_type_Field, boost::string_ref, &DofEntity::getNameRef >, const_mem_fun< DofEntity::interface_type_RefEntity, EntityType, &DofEntity::getEntType > > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
~VecManager()
Destructor.
Definition: VecManager.cpp:34
MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit, A4 &cpy_field_name)
Definition: VecManager.cpp:472
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsRows
store DOFs on rows for this problem
static const char *const FieldSpaceNames[]
Definition: definitions.h:178
#define CHKERR
Inline error check.
Definition: definitions.h:596
virtual MoFEMErrorCode get_field_ents(const FieldEntity_multiIndex **field_ents) const =0
Get field multi index.
DofIdx getNbDofsCol() const
MoFEMErrorCode vecCreateSeq(const std::string &name, RowColData rc, Vec *V) const
create local vector for problem
Definition: VecManager.cpp:36
MultiIndex Tag for field name.
VecManager(const MoFEM::Core &core)
Definition: VecManager.cpp:32
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
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)collective - need tu be run on all ...
Definition: VecManager.cpp:495
virtual MoFEMErrorCode get_fields(const Field_multiIndex **fields_ptr) const =0
Get fields multi-index from database.
MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit, A4 &cpy_field_name)
Definition: VecManager.cpp:337
DofIdx getNbLocalDofsCol() const
intrusive_ptr for managing petsc objects
Definition: AuxPETSc.hpp:124
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsCols
store DOFs on columns for this problem