v0.9.2
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 VecManager::vecScatterCreate(Vec xin, const std::string &x_problem,
160  const std::string &x_field_name, RowColData x_rc,
161  Vec yin, const std::string &y_problem,
162  const std::string &y_field_name, RowColData y_rc,
163  SmartPetscObj<VecScatter> &smart_newctx) const {
165  VecScatter newctx;
166  CHKERR vecScatterCreate(xin, x_problem, x_field_name, x_rc, yin, y_problem,
167  y_field_name, y_rc, &newctx);
168  smart_newctx = newctx;
170 }
171 
173 VecManager::vecScatterCreate(Vec xin, const std::string &x_problem,
174  RowColData x_rc, Vec yin,
175  const std::string &y_problem, RowColData y_rc,
176  SmartPetscObj<VecScatter> &smart_newctx) const {
178  VecScatter newctx;
179  CHKERR vecScatterCreate(xin, x_problem, x_rc, yin, y_problem, y_rc, &newctx);
180  smart_newctx = newctx;
182 }
183 
185  RowColData rc, Vec V,
186  InsertMode mode,
187  ScatterMode scatter_mode) const {
190  DofIdx nb_local_dofs, nb_ghost_dofs;
191  switch (rc) {
192  case ROW:
193  nb_local_dofs = problem_ptr->getNbLocalDofsRow();
194  nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
195  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
196  &problem_ptr->numeredDofsRows->get<PetscLocalIdx_mi_tag>());
197  break;
198  case COL:
199  nb_local_dofs = problem_ptr->getNbLocalDofsCol();
200  nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
201  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
202  &problem_ptr->numeredDofsCols->get<PetscLocalIdx_mi_tag>());
203  break;
204  default:
205  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
206  }
207  Vec Vlocal;
208  CHKERR VecGhostGetLocalForm(V, &Vlocal);
209  int size;
210  CHKERR VecGetLocalSize(V, &size);
211  if (size != nb_local_dofs) {
212  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
213  "data inconsistency: check ghost vector, problem with nb. of "
214  "local nodes %d != %d",
215  size, nb_local_dofs);
216  }
217  CHKERR VecGetLocalSize(Vlocal, &size);
218  if (size != nb_local_dofs + nb_ghost_dofs) {
219  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
220  "data inconsistency: check ghost vector, problem with nb. of "
221  "ghost nodes %d != ",
222  size, nb_local_dofs + nb_ghost_dofs);
223  }
224  NumeredDofEntityByLocalIdx::iterator miit = dofs->lower_bound(0);
225  NumeredDofEntityByLocalIdx::iterator hi_miit =
226  dofs->upper_bound(nb_local_dofs + nb_ghost_dofs);
227  DofIdx ii = 0;
228  switch (scatter_mode) {
229  case SCATTER_FORWARD: {
230  PetscScalar *array;
231  CHKERR VecGetArray(Vlocal, &array);
232  switch (mode) {
233  case INSERT_VALUES:
234  for (; miit != hi_miit; ++miit, ++ii)
235  array[ii] = (*miit)->getFieldData();
236  break;
237  case ADD_VALUES:
238  for (; miit != hi_miit; ++miit, ++ii)
239  array[ii] += (*miit)->getFieldData();
240  break;
241  default:
242  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
243  }
244  CHKERR VecRestoreArray(Vlocal, &array);
245  }; break;
246  case SCATTER_REVERSE: {
247  const PetscScalar *array;
248  CHKERR VecGetArrayRead(Vlocal, &array);
249  switch (mode) {
250  case INSERT_VALUES:
251  for (; miit != hi_miit; ++miit, ++ii) {
252  // std::cerr << *miit << std::endl;
253  // std::cerr << array[ii] << std::endl;
254  (*miit)->getFieldData() = array[ii];
255  }
256  break;
257  case ADD_VALUES:
258  for (; miit != hi_miit; ++miit, ++ii)
259  (*miit)->getFieldData() += array[ii];
260  break;
261  default:
262  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
263  }
264  CHKERR VecRestoreArrayRead(Vlocal, &array);
265  }; break;
266  default:
267  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
268  }
269  CHKERR VecGhostRestoreLocalForm(V, &Vlocal);
271 }
272 
274  RowColData rc, Vec V,
275  InsertMode mode,
276  ScatterMode scatter_mode) const {
277  const MoFEM::Interface &m_field = cOre;
278  const Problem *problem_ptr;
280  CHKERR m_field.get_problem(name, &problem_ptr);
281  CHKERR setLocalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
283 }
284 
287  Vec V, InsertMode mode,
288  ScatterMode scatter_mode) const {
291  DofsByGlobalIdx;
292  DofsByGlobalIdx *dofs;
293  DofIdx nb_dofs;
294  switch (rc) {
295  case ROW:
296  nb_dofs = problem_ptr->getNbDofsRow();
297  dofs = const_cast<DofsByGlobalIdx *>(
298  &problem_ptr->numeredDofsRows->get<PetscGlobalIdx_mi_tag>());
299  break;
300  case COL:
301  nb_dofs = problem_ptr->getNbDofsCol();
302  dofs = const_cast<DofsByGlobalIdx *>(
303  &problem_ptr->numeredDofsCols->get<PetscGlobalIdx_mi_tag>());
304  break;
305  default:
306  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
307  }
308  DofsByGlobalIdx::iterator miit = dofs->lower_bound(0);
309  DofsByGlobalIdx::iterator hi_miit = dofs->upper_bound(nb_dofs);
310  switch (scatter_mode) {
311  case SCATTER_REVERSE: {
312  VecScatter ctx;
313  Vec V_glob;
314  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
315  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
316  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
317  int size;
318  CHKERR VecGetSize(V_glob, &size);
319  if (size != nb_dofs) {
320  SETERRQ(
321  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
322  "data inconsistency: nb. of dofs and declared nb. dofs in database");
323  }
324  PetscScalar *array;
325  CHKERR VecGetArray(V_glob, &array);
326  switch (mode) {
327  case INSERT_VALUES:
328  for (; miit != hi_miit; miit++)
329  (*miit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
330  break;
331  case ADD_VALUES:
332  for (; miit != hi_miit; miit++)
333  (*miit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
334  break;
335  default:
336  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
337  }
338  CHKERR VecRestoreArray(V_glob, &array);
339  CHKERR VecScatterDestroy(&ctx);
340  CHKERR VecDestroy(&V_glob);
341  break;
342  }
343  default:
344  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
345  }
347 }
348 
350 VecManager::setGlobalGhostVector(const std::string &name, RowColData rc, Vec V,
351  InsertMode mode,
352  ScatterMode scatter_mode) const {
353  const MoFEM::Interface &m_field = cOre;
354  const Problem *problem_ptr;
356  CHKERR m_field.get_problem(name, &problem_ptr);
357  CHKERR setGlobalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
359 }
360 
361 template <int MODE> struct SetOtherLocalGhostVector {
362  template <typename A0, typename A1, typename A2, typename A3, typename A4>
363  inline MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit,
364  A4 &cpy_field_name) {
366  for (; miit != hi_miit; miit++) {
367  if (miit->get()->getHasLocalIndex()) {
368  auto diiiit =
369  dofs_ptr
370  ->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
371  .find(boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
372  (*miit)->getEntDofIdx()));
373  if (diiiit ==
374  dofs_ptr
375  ->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
376  .end()) {
377  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
378  "Automatic creation of entity and dof not implemented");
379  }
380  if (MODE == INSERT_VALUES)
381  (*diiiit)->getFieldData() = array[(*miit)->getPetscLocalDofIdx()];
382  else if (MODE == ADD_VALUES)
383  (*diiiit)->getFieldData() += array[(*miit)->getPetscLocalDofIdx()];
384  }
385  }
387  }
388 };
389 
391  const Problem *problem_ptr, const std::string &field_name,
392  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
393  ScatterMode scatter_mode) const {
394  const MoFEM::Interface &m_field = cOre;
395  auto fields_ptr = m_field.get_fields();
396  auto dofs_ptr = m_field.get_dofs();
399  DofsByName *dofs;
400  switch (rc) {
401  case ROW:
402  dofs = const_cast<DofsByName *>(
403  &problem_ptr->numeredDofsRows->get<FieldName_mi_tag>());
404  break;
405  case COL:
406  dofs = const_cast<DofsByName *>(
407  &problem_ptr->numeredDofsCols->get<FieldName_mi_tag>());
408  break;
409  default:
410  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
411  }
412  Field_multiIndex::index<FieldName_mi_tag>::type::iterator cpy_fit =
413  fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
414  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
415  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
416  "cpy field < %s > not found, (top tip: check spelling)",
417  cpy_field_name.c_str());
418  }
419  DofsByName::iterator miit = dofs->lower_bound(field_name);
420  if (miit == dofs->end()) {
421  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
422  "cpy field < %s > not found, (top tip: check spelling)",
423  field_name.c_str());
424  }
425  DofsByName::iterator hi_miit = dofs->upper_bound(field_name);
426  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
427  SETERRQ4(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
428  "fields have to have same space (%s) %s != (%s) %s",
429  (*miit)->getName().c_str(), FieldSpaceNames[(*miit)->getSpace()],
430  cpy_field_name.c_str(), FieldSpaceNames[(*cpy_fit)->getSpace()]);
431  }
432  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
433  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
434  "fields have to have same rank");
435  }
436 
437  switch (scatter_mode) {
438  case SCATTER_REVERSE: {
439 
440  PetscScalar *array;
441  CHKERR VecGetArray(V, &array);
442  if (mode == INSERT_VALUES)
443  CHKERR SetOtherLocalGhostVector<INSERT_VALUES>()(dofs_ptr, array, miit,
444  hi_miit, cpy_field_name);
445  else if (mode == ADD_VALUES)
446  CHKERR SetOtherLocalGhostVector<ADD_VALUES>()(dofs_ptr, array, miit,
447  hi_miit, cpy_field_name);
448  else
449  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
450  "Operation mode not implemented");
451  CHKERR VecRestoreArray(V, &array);
452  } break;
453  case SCATTER_FORWARD: {
454  for (; miit != hi_miit; miit++) {
455  if (!miit->get()->getHasLocalIndex())
456  continue;
457  DofEntity_multiIndex::index<
458  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator diiiit;
459  diiiit =
460  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
461  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
462  (*miit)->getEntDofIdx()));
463  if (diiiit ==
464  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
465  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
466  "no data to fill the vector (top tip: you want scatter forward "
467  "or scatter reverse?)");
468  }
469  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
470  (*diiiit)->getFieldData(), mode);
471  }
472  CHKERR VecAssemblyBegin(V);
473  CHKERR VecAssemblyEnd(V);
474  } break;
475  default:
476  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
477  }
479 }
480 
482  const std::string &name, const std::string &field_name,
483  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
484  ScatterMode scatter_mode) const {
485  const MoFEM::Interface &m_field = cOre;
486  const Problem *problem_ptr;
488  CHKERR m_field.get_problem(name, &problem_ptr);
489  CHKERR setOtherLocalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
490  V, mode, scatter_mode);
492 }
493 
494 template <int MODE> struct SetOtherGlobalGhostVector {
495  template <typename A0, typename A1, typename A2, typename A3, typename A4>
496  inline MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit,
497  A4 &cpy_field_name) {
499  for (; miit != hi_miit; miit++) {
500  auto diiiit =
501  dofs_ptr->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
502  .find(boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
503  (*miit)->getEntDofIdx()));
504  if (diiiit ==
505  dofs_ptr->template get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>()
506  .end()) {
507  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
508  "Automatic creation of entity and dof not implemented");
509  }
510  if (MODE == INSERT_VALUES)
511  (*diiiit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
512  else if (MODE == ADD_VALUES)
513  (*diiiit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
514  }
516  }
517 };
518 
520  const Problem *problem_ptr, const std::string &field_name,
521  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
522  ScatterMode scatter_mode) const {
523  const MoFEM::Interface &m_field = cOre;
524  auto fields_ptr = m_field.get_fields();
525  auto dofs_ptr = m_field.get_dofs();
526  auto *field_ents = m_field.get_field_ents();
528  typedef NumeredDofEntityByFieldName DofsByName;
529  DofsByName *dofs;
530  DofIdx nb_dofs;
531  switch (rc) {
532  case ROW:
533  nb_dofs = problem_ptr->getNbDofsRow();
534  dofs = const_cast<DofsByName *>(
535  &problem_ptr->numeredDofsRows->get<FieldName_mi_tag>());
536  break;
537  case COL:
538  nb_dofs = problem_ptr->getNbDofsCol();
539  dofs = const_cast<DofsByName *>(
540  &problem_ptr->numeredDofsCols->get<FieldName_mi_tag>());
541  break;
542  default:
543  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
544  }
545  auto cpy_fit = fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
546  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
547  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
548  "cpy field < %s > not found, (top tip: check spelling)",
549  cpy_field_name.c_str());
550  }
551  auto miit = dofs->lower_bound(field_name);
552  if (miit == dofs->end()) {
553  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
554  "problem field < %s > not found, (top tip: check spelling)",
555  field_name.c_str());
556  }
557  auto hi_miit = dofs->upper_bound(field_name);
558  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
559  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
560  "fields have to have same space");
561  }
562  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
563  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
564  "fields have to have same rank");
565  }
566  switch (scatter_mode) {
567  case SCATTER_REVERSE: {
568  Vec V_glob;
569  VecScatter ctx;
570  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
571  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
572  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
573  int size;
574  CHKERR VecGetSize(V_glob, &size);
575  if (size != nb_dofs)
576  SETERRQ(
577  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
578  "data inconsistency: nb. of dofs and declared nb. dofs in database");
579  PetscScalar *array;
580  CHKERR VecGetArray(V_glob, &array);
581 
582  if (mode == INSERT_VALUES)
584  dofs_ptr, array, miit, hi_miit, cpy_field_name);
585  else if (mode == ADD_VALUES)
586  CHKERR SetOtherGlobalGhostVector<ADD_VALUES>()(dofs_ptr, array, miit,
587  hi_miit, cpy_field_name);
588  else
589  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
590  "Operation mode not implemented");
591 
592  CHKERR VecRestoreArray(V_glob, &array);
593  CHKERR VecDestroy(&V_glob);
594  CHKERR VecScatterDestroy(&ctx);
595  } break;
596  case SCATTER_FORWARD: {
597  for (; miit != hi_miit; miit++) {
598  auto diiiit =
599  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
600  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
601  (*miit)->getEntDofIdx()));
602  if (diiiit ==
603  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
604  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
605  "No data to fill the vector (top tip: you want scatter forward "
606  "or scatter reverse?)");
607  }
608  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
609  (*diiiit)->getFieldData(), mode);
610  }
611  CHKERR VecAssemblyBegin(V);
612  CHKERR VecAssemblyEnd(V);
613  } break;
614  default:
615  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
616  }
618 }
619 
621  const std::string &name, const std::string &field_name,
622  const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode,
623  ScatterMode scatter_mode) const {
624  const MoFEM::Interface &m_field = cOre;
625  const Problem *problem_ptr;
627  CHKERR m_field.get_problem(name, &problem_ptr);
628  CHKERR setOtherGlobalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
629  V, mode, scatter_mode);
631 }
632 
633 } // 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:390
NumeredDofEntity_multiIndex::index< FieldName_mi_tag >::type NumeredDofEntityByFieldName
Numbered DoF multi-index by field name.
DofIdx getNbDofsRow() const
Deprecated interface functions.
MoFEM interface unique ID.
DofIdx getNbGhostDofsCol() const
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
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:507
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
virtual const FieldEntity_multiIndex * get_field_ents() const =0
Get the field ents object.
base class for all interface classes
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
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:514
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:192
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
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:184
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:286
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.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
~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:496
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsRows
store DOFs on rows for this problem
static const char *const FieldSpaceNames[]
Definition: definitions.h:184
#define CHKERR
Inline error check.
Definition: definitions.h:602
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
virtual const Problem * get_problem(const std::string &problem_name) const =0
Get the problem object.
int DofIdx
Index of DOF.
Definition: Types.hpp:29
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
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:519
MoFEMErrorCode operator()(A0 dofs_ptr, A1 array, A2 miit, A3 hi_miit, A4 &cpy_field_name)
Definition: VecManager.cpp:363
DofIdx getNbLocalDofsCol() const
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsCols
store DOFs on columns for this problem