v0.13.0
Files | Functions
Vectors (Vec)

Creating and scattering vectors on the mesh for given problem. More...

Collaboration diagram for Vectors (Vec):

Files

file  VecManager.hpp
 Interface managing vectors.
 

Functions

MoFEMErrorCode MoFEM::VecManager::vecCreateSeq (const std::string name, RowColData rc, Vec *V) const
 create local vector for problem More...
 
MoFEMErrorCode MoFEM::VecManager::vecCreateGhost (const std::string name, RowColData rc, Vec *V) const
 create ghost vector for problem (collective) More...
 
MoFEMErrorCode MoFEM::VecManager::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) More...
 
MoFEMErrorCode MoFEM::VecManager::vecScatterCreate (Vec xin, const std::string x_problem, RowColData x_rc, Vec yin, const std::string y_problem, RowColData y_rc, VecScatter *newctx) const
 create scatter for vectors from one to another problem (collective) More...
 
MoFEMErrorCode MoFEM::VecManager::setLocalGhostVector (const Problem *problem_ptr, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
 set values of vector from/to meshdatabase More...
 
MoFEMErrorCode MoFEM::VecManager::setLocalGhostVector (const std::string name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
 set values of vector from/to meshdatabase More...
 
MoFEMErrorCode MoFEM::VecManager::setGlobalGhostVector (const Problem *problem_ptr, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
 set values of vector from/to mesh database (collective) More...
 
MoFEMErrorCode MoFEM::VecManager::setGlobalGhostVector (const std::string name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
 set values of vector from/to mesh database (collective) More...
 
MoFEMErrorCode MoFEM::VecManager::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. More...
 
MoFEMErrorCode MoFEM::VecManager::setOtherLocalGhostVector (const std::string name, 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. More...
 
MoFEMErrorCode MoFEM::VecManager::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) More...
 
MoFEMErrorCode MoFEM::VecManager::setOtherGlobalGhostVector (const std::string name, 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) More...
 

Detailed Description

Creating and scattering vectors on the mesh for given problem.

Function Documentation

◆ setGlobalGhostVector() [1/2]

MoFEMErrorCode MoFEM::VecManager::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 in communicator

Parameters
pointerto porblem struture
RowColDatafor row or column (i.e. Row,Col)
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities form V vector.

Definition at line 298 of file VecManager.cpp.

300  {
303  DofsByGlobalIdx;
304  DofsByGlobalIdx *dofs;
305  DofIdx nb_dofs;
306  switch (rc) {
307  case ROW:
308  nb_dofs = problem_ptr->getNbDofsRow();
309  dofs = const_cast<DofsByGlobalIdx *>(
310  &problem_ptr->numeredRowDofsPtr->get<PetscGlobalIdx_mi_tag>());
311  break;
312  case COL:
313  nb_dofs = problem_ptr->getNbDofsCol();
314  dofs = const_cast<DofsByGlobalIdx *>(
315  &problem_ptr->numeredColDofsPtr->get<PetscGlobalIdx_mi_tag>());
316  break;
317  default:
318  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
319  "Function works only for ROWs and COLs");
320  }
321  DofsByGlobalIdx::iterator miit = dofs->lower_bound(0);
322  DofsByGlobalIdx::iterator hi_miit = dofs->upper_bound(nb_dofs);
323  auto comm = PetscObjectComm((PetscObject)V);
324  switch (scatter_mode) {
325  case SCATTER_REVERSE: {
326  VecScatter ctx;
327  Vec V_glob;
328  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
329  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
330  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
331  int size;
332  CHKERR VecGetSize(V_glob, &size);
333  if (size != nb_dofs) {
334  SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
335  "Size of vector is inconsistent with size of problem. You could "
336  "use wrong vector with wrong problem, or you created vector "
337  "before you remove DOFs from problem.");
338  }
339  PetscScalar *array;
340  CHKERR VecGetArray(V_glob, &array);
341  switch (mode) {
342  case INSERT_VALUES:
343  for (; miit != hi_miit; miit++)
344  (*miit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
345  break;
346  case ADD_VALUES:
347  for (; miit != hi_miit; miit++)
348  (*miit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
349  break;
350  default:
351  SETERRQ(comm, MOFEM_NOT_IMPLEMENTED, "not implemented");
352  }
353  CHKERR VecRestoreArray(V_glob, &array);
354  CHKERR VecScatterDestroy(&ctx);
355  CHKERR VecDestroy(&V_glob);
356  break;
357  }
358  default:
359  SETERRQ(comm, MOFEM_NOT_IMPLEMENTED, "not implemented");
360  }
362 }
@ COL
Definition: definitions.h:136
@ ROW
Definition: definitions.h:136
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ 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
const FTensor::Tensor2< T, Dim, Dim > Vec
int DofIdx
Index of DOF.
Definition: Types.hpp:29

◆ setGlobalGhostVector() [2/2]

MoFEMErrorCode MoFEM::VecManager::setGlobalGhostVector ( const std::string  name,
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 in communicator

Parameters
nameof the problem
RowColDatafor row or column (i.e. Row,Col)
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities form V vector.

Definition at line 365 of file VecManager.cpp.

367  {
368  const MoFEM::Interface &m_field = cOre;
369  const Problem *problem_ptr;
371  CHKERR m_field.get_problem(name, &problem_ptr);
372  CHKERR setGlobalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
374 }
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
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:298
Deprecated interface functions.
const MoFEM::Interface & cOre
Definition: VecManager.hpp:38

◆ setLocalGhostVector() [1/2]

MoFEMErrorCode MoFEM::VecManager::setLocalGhostVector ( const Problem problem_ptr,
RowColData  rc,
Vec  V,
InsertMode  mode,
ScatterMode  scatter_mode 
) const

set values of vector from/to meshdatabase

Parameters
pointerto problem struture
RowColDatafor row or column:e (i.e. Row,Col)
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities from V vector.

SCATTER_FORWARD set vector V from data field entities

Definition at line 196 of file VecManager.cpp.

199  {
202  DofIdx nb_local_dofs, nb_ghost_dofs;
203  switch (rc) {
204  case ROW:
205  nb_local_dofs = problem_ptr->getNbLocalDofsRow();
206  nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
207  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
208  &problem_ptr->numeredRowDofsPtr->get<PetscLocalIdx_mi_tag>());
209  break;
210  case COL:
211  nb_local_dofs = problem_ptr->getNbLocalDofsCol();
212  nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
213  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
214  &problem_ptr->numeredColDofsPtr->get<PetscLocalIdx_mi_tag>());
215  break;
216  default:
217  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
218  }
219  Vec Vlocal;
220  CHKERR VecGhostGetLocalForm(V, &Vlocal);
221  int size;
222  CHKERR VecGetLocalSize(V, &size);
223  if (size != nb_local_dofs) {
224  SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
225  "data inconsistency: check ghost vector, problem <%s> with nb. of "
226  "local nodes %d != %d",
227  problem_ptr->getName().c_str(), size, nb_local_dofs);
228  }
229  CHKERR VecGetLocalSize(Vlocal, &size);
230  if (size != nb_local_dofs + nb_ghost_dofs) {
231  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
232  "data inconsistency: check ghost vector, problem with nb. of "
233  "ghost nodes %d != ",
234  size, nb_local_dofs + nb_ghost_dofs);
235  }
236  NumeredDofEntityByLocalIdx::iterator miit = dofs->lower_bound(0);
237  NumeredDofEntityByLocalIdx::iterator hi_miit =
238  dofs->upper_bound(nb_local_dofs + nb_ghost_dofs);
239  DofIdx ii = 0;
240  switch (scatter_mode) {
241  case SCATTER_FORWARD: {
242  PetscScalar *array;
243  CHKERR VecGetArray(Vlocal, &array);
244  switch (mode) {
245  case INSERT_VALUES:
246  for (; miit != hi_miit; ++miit, ++ii)
247  array[ii] = (*miit)->getFieldData();
248  break;
249  case ADD_VALUES:
250  for (; miit != hi_miit; ++miit, ++ii)
251  array[ii] += (*miit)->getFieldData();
252  break;
253  default:
254  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
255  }
256  CHKERR VecRestoreArray(Vlocal, &array);
257  }; break;
258  case SCATTER_REVERSE: {
259  const PetscScalar *array;
260  CHKERR VecGetArrayRead(Vlocal, &array);
261  switch (mode) {
262  case INSERT_VALUES:
263  for (; miit != hi_miit; ++miit, ++ii) {
264  // std::cerr << *miit << std::endl;
265  // std::cerr << array[ii] << std::endl;
266  (*miit)->getFieldData() = array[ii];
267  }
268  break;
269  case ADD_VALUES:
270  for (; miit != hi_miit; ++miit, ++ii)
271  (*miit)->getFieldData() += array[ii];
272  break;
273  default:
274  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
275  }
276  CHKERR VecRestoreArrayRead(Vlocal, &array);
277  }; break;
278  default:
279  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
280  }
281  CHKERR VecGhostRestoreLocalForm(V, &Vlocal);
283 }
NumeredDofEntity_multiIndex::index< PetscLocalIdx_mi_tag >::type NumeredDofEntityByLocalIdx
Numbered DoF multi-index by local index.

◆ setLocalGhostVector() [2/2]

MoFEMErrorCode MoFEM::VecManager::setLocalGhostVector ( const std::string  name,
RowColData  rc,
Vec  V,
InsertMode  mode,
ScatterMode  scatter_mode 
) const

set values of vector from/to meshdatabase

Parameters
nameof the problem
RowColDatafor row or column:e (i.e. Row,Col)
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities from V vector.

SCATTER_FORWARD set vector V from data field entities

Definition at line 285 of file VecManager.cpp.

288  {
289  const MoFEM::Interface &m_field = cOre;
290  const Problem *problem_ptr;
292  CHKERR m_field.get_problem(name, &problem_ptr);
293  CHKERR setLocalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
295 }
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:196

◆ setOtherGlobalGhostVector() [1/2]

MoFEMErrorCode MoFEM::VecManager::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 processors in communicator

Parameters
problem_ptrpointer to problem
field_namefield name used for indexing petsc vectors used in the problem
cpy_fieldfield name where data from vector are stored
RowColDatafor row or column
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities form V vector.

Definition at line 541 of file VecManager.cpp.

544  {
545  const MoFEM::Interface &m_field = cOre;
546  auto fields_ptr = m_field.get_fields();
547  auto dofs_ptr = m_field.get_dofs();
548  auto *field_ents = m_field.get_field_ents();
550  NumeredDofEntityByUId *dofs;
551  DofIdx nb_dofs;
552  switch (rc) {
553  case ROW:
554  nb_dofs = problem_ptr->getNbDofsRow();
555  dofs = const_cast<NumeredDofEntityByUId *>(
556  &problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>());
557  break;
558  case COL:
559  nb_dofs = problem_ptr->getNbDofsCol();
560  dofs = const_cast<NumeredDofEntityByUId *>(
561  &problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>());
562  break;
563  default:
564  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
565  }
566  auto cpy_fit = fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
567  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
568  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
569  "cpy field < %s > not found, (top tip: check spelling)",
570  cpy_field_name.c_str());
571  }
572  const auto cpy_bit_number = (*cpy_fit)->getBitNumber();
573 
574  auto miit = dofs->lower_bound(
576  if (miit == dofs->end()) {
577  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
578  "problem field < %s > not found, (top tip: check spelling)",
579  field_name.c_str());
580  }
581  auto hi_miit = dofs->upper_bound(
583  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
584  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
585  "fields have to have same space");
586  }
587  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
588  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
589  "fields have to have same rank");
590  }
591  switch (scatter_mode) {
592  case SCATTER_REVERSE: {
593  Vec V_glob;
594  VecScatter ctx;
595  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
596  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
597  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
598  int size;
599  CHKERR VecGetSize(V_glob, &size);
600  if (size != nb_dofs)
601  SETERRQ(
602  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
603  "data inconsistency: nb. of dofs and declared nb. dofs in database");
604  PetscScalar *array;
605  CHKERR VecGetArray(V_glob, &array);
606 
607  if (mode == INSERT_VALUES)
608  CHKERR SetOtherGlobalGhostVector<INSERT_VALUES>()(
609  dofs_ptr, array, miit, hi_miit, cpy_bit_number);
610  else if (mode == ADD_VALUES)
611  CHKERR SetOtherGlobalGhostVector<ADD_VALUES>()(dofs_ptr, array, miit,
612  hi_miit, cpy_bit_number);
613  else
614  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
615  "Operation mode not implemented");
616 
617  CHKERR VecRestoreArray(V_glob, &array);
618  CHKERR VecDestroy(&V_glob);
619  CHKERR VecScatterDestroy(&ctx);
620  } break;
621  case SCATTER_FORWARD: {
622  for (; miit != hi_miit; miit++) {
623  const auto uid = DofEntity::getUniqueIdCalculate(
624  (*miit)->getEntDofIdx(), FieldEntity::getLocalUniqueIdCalculate(
625  cpy_bit_number, (*miit)->getEnt()));
626  auto diiiit = dofs_ptr->get<Unique_mi_tag>().find(uid);
627  if (diiiit == dofs_ptr->get<Unique_mi_tag>().end()) {
628  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
629  "No data to fill the vector (top tip: you want scatter forward "
630  "or scatter reverse?)");
631  }
632  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
633  (*diiiit)->getFieldData(), mode);
634  }
635  CHKERR VecAssemblyBegin(V);
636  CHKERR VecAssemblyEnd(V);
637  } break;
638  default:
639  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
640  }
642 }
@ MOFEM_NOT_FOUND
Definition: definitions.h:46
NumeredDofEntity_multiIndex::index< Unique_mi_tag >::type NumeredDofEntityByUId
Numbered DoF multi-index by UId.
virtual const FieldEntity_multiIndex * get_field_ents() const =0
Get the field ents object.
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
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)

◆ setOtherGlobalGhostVector() [2/2]

MoFEMErrorCode MoFEM::VecManager::setOtherGlobalGhostVector ( const std::string  name,
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)

Deprecated:
VecManager

collective - need tu be run on all processors in communicator

Parameters
nameproblem name
field_namefield name used for indexing petsc vectors used in the problem
cpy_fieldfield name where data from vector are stored
RowColDatafor row or column
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities form V vector.

Definition at line 644 of file VecManager.cpp.

647  {
648  const MoFEM::Interface &m_field = cOre;
649  const Problem *problem_ptr;
651  CHKERR m_field.get_problem(name, &problem_ptr);
652  CHKERR setOtherGlobalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
653  V, mode, scatter_mode);
655 }
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:541

◆ setOtherLocalGhostVector() [1/2]

MoFEMErrorCode MoFEM::VecManager::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.

Parameters
pointerto problem multi_index
field_namefield name used for indexing petsc vectors used in the problem
cpy_fieldfield name where data from vector are stored
RowColDatafor row or column
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities form V vector.

Definition at line 410 of file VecManager.cpp.

413  {
414  const MoFEM::Interface &m_field = cOre;
415  auto fields_ptr = m_field.get_fields();
416  auto dofs_ptr = m_field.get_dofs();
419  DofsByUId *dofs;
420 
421  switch (rc) {
422  case ROW:
423  dofs = const_cast<DofsByUId *>(
424  &problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>());
425  break;
426  case COL:
427  dofs = const_cast<DofsByUId *>(
428  &problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>());
429  break;
430  default:
431  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
432  }
433 
434  auto cpy_fit = fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
435  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end())
436  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
437  "cpy field < %s > not found, (top tip: check spelling)",
438  cpy_field_name.c_str());
439  const auto cpy_bit_number = (*cpy_fit)->getBitNumber();
440 
441  auto miit = dofs->lower_bound(
443  if (miit == dofs->end()) {
444  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
445  "cpy field < %s > not found, (top tip: check spelling)",
446  field_name.c_str());
447  }
448  auto hi_miit = dofs->upper_bound(
450 
451  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
452  SETERRQ4(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
453  "fields have to have same space (%s) %s != (%s) %s",
454  (*miit)->getName().c_str(), FieldSpaceNames[(*miit)->getSpace()],
455  cpy_field_name.c_str(), FieldSpaceNames[(*cpy_fit)->getSpace()]);
456  }
457  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
458  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
459  "fields have to have same rank");
460  }
461 
462  switch (scatter_mode) {
463  case SCATTER_REVERSE: {
464 
465  PetscScalar *array;
466  CHKERR VecGetArray(V, &array);
467  if (mode == INSERT_VALUES)
468  CHKERR SetOtherLocalGhostVector<INSERT_VALUES>()(
469  dofs_ptr, array, miit, hi_miit, cpy_bit_number, cpy_field_name);
470  else if (mode == ADD_VALUES)
471  CHKERR SetOtherLocalGhostVector<ADD_VALUES>()(
472  dofs_ptr, array, miit, hi_miit, cpy_bit_number, cpy_field_name);
473  else
474  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
475  "Operation mode not implemented");
476  CHKERR VecRestoreArray(V, &array);
477  } break;
478  case SCATTER_FORWARD: {
479  for (; miit != hi_miit; miit++) {
480  if (!miit->get()->getHasLocalIndex())
481  continue;
482  const auto uid = DofEntity::getUniqueIdCalculate(
483  (*miit)->getEntDofIdx(), FieldEntity::getLocalUniqueIdCalculate(
484  cpy_bit_number, (*miit)->getEnt()));
485  auto diiiit = dofs_ptr->get<Unique_mi_tag>().find(uid);
486  if (diiiit == dofs_ptr->get<Unique_mi_tag>().end()) {
487  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
488  "no data to fill the vector (top tip: you want scatter forward "
489  "or scatter reverse?)");
490  }
491  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
492  (*diiiit)->getFieldData(), mode);
493  }
494  CHKERR VecAssemblyBegin(V);
495  CHKERR VecAssemblyEnd(V);
496  } break;
497  default:
498  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
499  }
501 }
static const char *const FieldSpaceNames[]
Definition: definitions.h:105

◆ setOtherLocalGhostVector() [2/2]

MoFEMErrorCode MoFEM::VecManager::setOtherLocalGhostVector ( const std::string  name,
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.

Parameters
nameproblem name
field_namefield name used for indexing petsc vectors used in the problem
cpy_fieldfield name where data from vector are stored
RowColDatafor row or column
Vvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities form V vector.

Definition at line 503 of file VecManager.cpp.

506  {
507  const MoFEM::Interface &m_field = cOre;
508  const Problem *problem_ptr;
510  CHKERR m_field.get_problem(name, &problem_ptr);
511  CHKERR setOtherLocalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
512  V, mode, scatter_mode);
514 }
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:410

◆ vecCreateGhost()

MoFEMErrorCode MoFEM::VecManager::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 communicator

Parameters
nameproblem name
RowColDataspecify what data is taken from Row, Col or Data
Vecthe vector where data is stored

Definition at line 74 of file VecManager.cpp.

75  {
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  CHKERR ::VecCreateGhost(PETSC_COMM_WORLD, nb_local_dofs, nb_dofs,
112  nb_ghost_dofs, &ghost_idx[0], V);
114 }

◆ vecCreateSeq()

MoFEMErrorCode MoFEM::VecManager::vecCreateSeq ( const std::string  name,
RowColData  rc,
Vec *  V 
) const

create local vector for problem

Parameters
nameproblem name
RowColDataspecify what data is taken from Row, Col or Data
Vecthe vector where data is stored

Definition at line 51 of file VecManager.cpp.

52  {
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 }

◆ vecScatterCreate() [1/2]

MoFEMErrorCode MoFEM::VecManager::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 one problem is scattered to another.

Parameters
xinvector
x_probleproblem name
x_fieldname
yinvector
y_problemproblem name
y_field_name
Return values
newctxscatter

Definition at line 125 of file VecManager.cpp.

128  {
129  const MoFEM::Interface &m_field = cOre;
131  std::vector<int> idx(0), idy(0);
132  CHKERR m_field.getInterface<ISManager>()
133  ->isCreateFromProblemFieldToOtherProblemField(
134  x_problem, x_field_name, x_rc, y_problem, y_field_name, y_rc, idx,
135  idy);
136  IS ix, iy;
137  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
138  PETSC_USE_POINTER, &ix);
139  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
140  PETSC_USE_POINTER, &iy);
141  CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
142  CHKERR ISDestroy(&ix);
143  CHKERR ISDestroy(&iy);
145 }
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ vecScatterCreate() [2/2]

MoFEMErrorCode MoFEM::VecManager::vecScatterCreate ( Vec  xin,
const std::string  x_problem,
RowColData  x_rc,
Vec  yin,
const std::string  y_problem,
RowColData  y_rc,
VecScatter *  newctx 
) const

create scatter for vectors from one to another problem (collective)

Parameters
xinvector
x_probleproblem name
yinvector
y_problemproblem name
Return values
newctxscatter

Definition at line 147 of file VecManager.cpp.

149  {
150  const MoFEM::Interface &m_field = cOre;
152  std::vector<int> idx(0), idy(0);
153  CHKERR m_field.getInterface<ISManager>()->isCreateFromProblemToOtherProblem(
154  x_problem, x_rc, y_problem, y_rc, idx, idy);
155  IS ix, iy;
156  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
157  PETSC_USE_POINTER, &ix);
158  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
159  PETSC_USE_POINTER, &iy);
160  if (dEbug) {
161  ISView(ix, PETSC_VIEWER_STDOUT_WORLD);
162  ISView(iy, PETSC_VIEWER_STDOUT_WORLD);
163  }
164  CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
165  CHKERR ISDestroy(&ix);
166  CHKERR ISDestroy(&iy);
168 }