v0.8.4
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 vectorsManaging problems, build and partitioning.
 

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)collective - need to be run on all processors in communicator 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)User specify what name of field on one problem is scattered to another. 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 form 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)collective - need tu be run on all processors in communicator 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)collective - need tu be run on all processors in communicator 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)collective - need tu be run on all processors in communicator. 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)collective - need tu be run on all processors in communicator. More...
 

Create vectors

DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::VecCreateSeq (const std::string &name, RowColData rc, Vec *V) const
 create local vector for problem More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::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 More...
 

Create IS

DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::ISCreateProblemFieldAndRank (const std::string &problem, RowColData rc, const std::string &field, int min_coeff_idx, int max_coeff_idx, IS *is, int verb=-1) const
 create IS for given problem, field and rank range (collective) More...
 

Scatter vectors

DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::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, int verb=-1) const
 create scatter for vectors form one to another problem (collective)User specify what name of field on one problem is scattered to another. More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::VecScatterCreate (Vec xin, const std::string &x_problem, RowColData x_rc, Vec yin, const std::string &y_problem, RowColData y_rc, VecScatter *newctx, int verb=-1) const
 create scatter for vectors form one to another problem (collective) More...
 

Set vector and mesh values

DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_local_ghost_vector (const Problem *problem_ptr, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
 set values of vector from/to meshdatabase More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_local_ghost_vector (const std::string &name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
 set values of vector from/to meshdatabase More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_global_ghost_vector (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 More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_global_ghost_vector (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 More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_local_ghost_vector (const Problem *problem_ptr, const std::string &fiel_name, const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode, int verb=-1)
 Copy vector to field which is not part of the problem. More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_local_ghost_vector (const std::string &name, const std::string &field_name, const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode, int verb=-1)
 Copy vector to field which is not part of the problem. More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_global_ghost_vector (const Problem *problem_ptr, const std::string &field_name, const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode, int verb=-1)
 Copy vector to field which is not part of the problem (collective)collective - need tu be run on all processors in communicator. More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_global_ghost_vector (const std::string &name, const std::string &field_name, const std::string &cpy_field_name, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode, int verb=-1)
 Copy vector to field which is not part of the problem (collective)collective - need tu be run on all processors in communicator. More...
 

Detailed Description

Creating and scattering vectors on the mesh for given problem.

Function Documentation

◆ ISCreateProblemFieldAndRank()

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::ISCreateProblemFieldAndRank ( const std::string &  problem,
RowColData  rc,
const std::string &  field,
int  min_coeff_idx,
int  max_coeff_idx,
IS *  is,
int  verb = -1 
) const

create IS for given problem, field and rank range (collective)

Deprecated:
Use ISManager
Parameters
problemname
rcROW or COL
fieldname
min_coeff_idx
max_coeff_idx
Return values
isout value

Definition at line 141 of file DeprecatedCoreInterface.cpp.

143  {
144  return getInterface<ISManager>()->isCreateProblemFieldAndRank(
145  problem, rc, field, min_coeff_idx, max_coeff_idx, is);
146 }

◆ set_global_ghost_vector() [1/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_global_ghost_vector ( 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

Deprecated:
use VecManager
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 194 of file DeprecatedCoreInterface.cpp.

196  {
197  return getInterface<VecManager>()->setGlobalGhostVector(problem_ptr, rc, V,
198  mode, scatter_mode);
199 }

◆ set_global_ghost_vector() [2/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_global_ghost_vector ( 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

Deprecated:
use VecManager
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 200 of file DeprecatedCoreInterface.cpp.

202  {
203  return getInterface<VecManager>()->setGlobalGhostVector(name, rc, V, mode,
204  scatter_mode);
205 }

◆ set_local_ghost_vector() [1/2]

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

set values of vector from/to meshdatabase

Deprecated:
use VecManager
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 182 of file DeprecatedCoreInterface.cpp.

184  {
185  return getInterface<VecManager>()->setLocalGhostVector(problem_ptr, rc, V,
186  mode, scatter_mode);
187 }

◆ set_local_ghost_vector() [2/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_local_ghost_vector ( const std::string &  name,
RowColData  rc,
Vec  V,
InsertMode  mode,
ScatterMode  scatter_mode 
) const

set values of vector from/to meshdatabase

Deprecated:
use VecManager
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 188 of file DeprecatedCoreInterface.cpp.

190  {
191  return getInterface<VecManager>()->setLocalGhostVector(name, rc, V, mode,
192  scatter_mode);
193 }

◆ set_other_global_ghost_vector() [1/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_global_ghost_vector ( const Problem problem_ptr,
const std::string &  field_name,
const std::string &  cpy_field_name,
RowColData  rc,
Vec  V,
InsertMode  mode,
ScatterMode  scatter_mode,
int  verb = -1 
)

Copy vector to field which is not part of the problem (collective)collective - need tu be run on all processors in communicator.

Deprecated:
use VecManager
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 220 of file DeprecatedCoreInterface.cpp.

223  {
224  return getInterface<VecManager>()->setOtherGlobalGhostVector(
225  problem_ptr, field_name, cpy_field_name, rc, V, mode, scatter_mode);
226 }

◆ set_other_global_ghost_vector() [2/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_global_ghost_vector ( const std::string &  name,
const std::string &  field_name,
const std::string &  cpy_field_name,
RowColData  rc,
Vec  V,
InsertMode  mode,
ScatterMode  scatter_mode,
int  verb = -1 
)

Copy vector to field which is not part of the problem (collective)collective - need tu be run on all processors in communicator.

Deprecated:
use VecManager
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 227 of file DeprecatedCoreInterface.cpp.

230  {
231  return getInterface<VecManager>()->setOtherGlobalGhostVector(
232  name, field_name, cpy_field_name, rc, V, mode, scatter_mode);
233 }

◆ set_other_local_ghost_vector() [1/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_local_ghost_vector ( const Problem problem_ptr,
const std::string &  fiel_name,
const std::string &  cpy_field_name,
RowColData  rc,
Vec  V,
InsertMode  mode,
ScatterMode  scatter_mode,
int  verb = -1 
)

Copy vector to field which is not part of the problem.

Deprecated:
use VecManager
Parameters
pointerto poroblem 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 206 of file DeprecatedCoreInterface.cpp.

209  {
210  return getInterface<VecManager>()->setOtherLocalGhostVector(
211  problem_ptr, field_name, cpy_field_name, rc, V, mode, scatter_mode);
212 }

◆ set_other_local_ghost_vector() [2/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_other_local_ghost_vector ( const std::string &  name,
const std::string &  field_name,
const std::string &  cpy_field_name,
RowColData  rc,
Vec  V,
InsertMode  mode,
ScatterMode  scatter_mode,
int  verb = -1 
)

Copy vector to field which is not part of the problem.

Deprecated:
use VecManager
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 213 of file DeprecatedCoreInterface.cpp.

216  {
217  return getInterface<VecManager>()->setOtherLocalGhostVector(
218  name, field_name, cpy_field_name, rc, V, mode, scatter_mode);
219 }

◆ 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 250 of file VecManager.cpp.

252  {
254  typedef NumeredDofEntity_multiIndex::index<PetscGlobalIdx_mi_tag>::type
255  DofsByGlobalIdx;
256  DofsByGlobalIdx *dofs;
257  DofIdx nb_dofs;
258  switch (rc) {
259  case ROW:
260  nb_dofs = problem_ptr->getNbDofsRow();
261  dofs = const_cast<DofsByGlobalIdx *>(
262  &problem_ptr->numeredDofsRows->get<PetscGlobalIdx_mi_tag>());
263  break;
264  case COL:
265  nb_dofs = problem_ptr->getNbDofsCol();
266  dofs = const_cast<DofsByGlobalIdx *>(
267  &problem_ptr->numeredDofsCols->get<PetscGlobalIdx_mi_tag>());
268  break;
269  default:
270  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
271  }
272  DofsByGlobalIdx::iterator miit = dofs->lower_bound(0);
273  DofsByGlobalIdx::iterator hi_miit = dofs->upper_bound(nb_dofs);
274  switch (scatter_mode) {
275  case SCATTER_REVERSE: {
276  VecScatter ctx;
277  Vec V_glob;
278  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
279  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
280  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
281  int size;
282  CHKERR VecGetSize(V_glob, &size);
283  if (size != nb_dofs) {
284  SETERRQ(
285  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
286  "data inconsistency: nb. of dofs and declared nb. dofs in database");
287  }
288  PetscScalar *array;
289  CHKERR VecGetArray(V_glob, &array);
290  switch (mode) {
291  case INSERT_VALUES:
292  for (; miit != hi_miit; miit++)
293  (*miit)->getFieldData() = array[(*miit)->getPetscGlobalDofIdx()];
294  break;
295  case ADD_VALUES:
296  for (; miit != hi_miit; miit++)
297  (*miit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
298  break;
299  default:
300  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
301  }
302  CHKERR VecRestoreArray(V_glob, &array);
303  CHKERR VecScatterDestroy(&ctx);
304  CHKERR VecDestroy(&V_glob);
305  break;
306  }
307  default:
308  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
309  }
311 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
int DofIdx
Index of DOF.
Definition: Common.hpp:126

◆ 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 314 of file VecManager.cpp.

316  {
317  const MoFEM::Interface &m_field = cOre;
318  const Problem *problem_ptr;
320  CHKERR m_field.get_problem(name, &problem_ptr);
321  CHKERR setGlobalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
323 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
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:250
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443

◆ 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 148 of file VecManager.cpp.

151  {
154  DofIdx nb_local_dofs, nb_ghost_dofs;
155  switch (rc) {
156  case ROW:
157  nb_local_dofs = problem_ptr->getNbLocalDofsRow();
158  nb_ghost_dofs = problem_ptr->getNbGhostDofsRow();
159  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
160  &problem_ptr->numeredDofsRows->get<PetscLocalIdx_mi_tag>());
161  break;
162  case COL:
163  nb_local_dofs = problem_ptr->getNbLocalDofsCol();
164  nb_ghost_dofs = problem_ptr->getNbGhostDofsCol();
165  dofs = const_cast<NumeredDofEntityByLocalIdx *>(
166  &problem_ptr->numeredDofsCols->get<PetscLocalIdx_mi_tag>());
167  break;
168  default:
169  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
170  }
171  Vec Vlocal;
172  CHKERR VecGhostGetLocalForm(V, &Vlocal);
173  int size;
174  CHKERR VecGetLocalSize(V, &size);
175  if (size != nb_local_dofs) {
176  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
177  "data inconsistency: check ghost vector, problem with nb. of "
178  "local nodes %d != %d",
179  size, nb_local_dofs);
180  }
181  CHKERR VecGetLocalSize(Vlocal, &size);
182  if (size != nb_local_dofs + nb_ghost_dofs) {
183  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
184  "data inconsistency: check ghost vector, problem with nb. of "
185  "ghost nodes %d != ",
186  size, nb_local_dofs + nb_ghost_dofs);
187  }
188  NumeredDofEntityByLocalIdx::iterator miit = dofs->lower_bound(0);
189  NumeredDofEntityByLocalIdx::iterator hi_miit =
190  dofs->upper_bound(nb_local_dofs + nb_ghost_dofs);
191  DofIdx ii = 0;
192  switch (scatter_mode) {
193  case SCATTER_FORWARD: {
194  PetscScalar *array;
195  CHKERR VecGetArray(Vlocal, &array);
196  switch (mode) {
197  case INSERT_VALUES:
198  for (; miit != hi_miit; ++miit, ++ii)
199  array[ii] = (*miit)->getFieldData();
200  break;
201  case ADD_VALUES:
202  for (; miit != hi_miit; ++miit, ++ii)
203  array[ii] += (*miit)->getFieldData();
204  break;
205  default:
206  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
207  }
208  CHKERR VecRestoreArray(Vlocal, &array);
209  }; break;
210  case SCATTER_REVERSE: {
211  const PetscScalar *array;
212  CHKERR VecGetArrayRead(Vlocal, &array);
213  switch (mode) {
214  case INSERT_VALUES:
215  for (; miit != hi_miit; ++miit, ++ii) {
216  // std::cerr << *miit << std::endl;
217  // std::cerr << array[ii] << std::endl;
218  (*miit)->getFieldData() = array[ii];
219  }
220  break;
221  case ADD_VALUES:
222  for (; miit != hi_miit; ++miit, ++ii)
223  (*miit)->getFieldData() += array[ii];
224  break;
225  default:
226  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
227  }
228  CHKERR VecRestoreArrayRead(Vlocal, &array);
229  }; break;
230  default:
231  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
232  }
233  CHKERR VecGhostRestoreLocalForm(V, &Vlocal);
235 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
NumeredDofEntity_multiIndex::index< PetscLocalIdx_mi_tag >::type NumeredDofEntityByLocalIdx
Numbered DoF multi-index by local index.
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
int DofIdx
Index of DOF.
Definition: Common.hpp:126

◆ 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 237 of file VecManager.cpp.

240  {
241  const MoFEM::Interface &m_field = cOre;
242  const Problem *problem_ptr;
244  CHKERR m_field.get_problem(name, &problem_ptr);
245  CHKERR setLocalGhostVector(problem_ptr, rc, V, mode, scatter_mode);
247 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
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:148
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443

◆ 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 451 of file VecManager.cpp.

454  {
455  const MoFEM::Interface &m_field = cOre;
456  const Field_multiIndex *fields_ptr;
457  const DofEntity_multiIndex *dofs_ptr;
458  const FieldEntity_multiIndex *field_ents;
460  CHKERR m_field.get_fields(&fields_ptr);
461  CHKERR m_field.get_dofs(&dofs_ptr);
462  CHKERR m_field.get_field_ents(&field_ents);
463  typedef NumeredDofEntityByFieldName DofsByName;
464  DofsByName *dofs;
465  DofIdx nb_dofs;
466  switch (rc) {
467  case ROW:
468  nb_dofs = problem_ptr->getNbDofsRow();
469  dofs = const_cast<DofsByName *>(
470  &problem_ptr->numeredDofsRows->get<FieldName_mi_tag>());
471  break;
472  case COL:
473  nb_dofs = problem_ptr->getNbDofsCol();
474  dofs = const_cast<DofsByName *>(
475  &problem_ptr->numeredDofsCols->get<FieldName_mi_tag>());
476  break;
477  default:
478  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
479  }
480  Field_multiIndex::index<FieldName_mi_tag>::type::iterator cpy_fit =
481  fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
482  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
483  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
484  "cpy field < %s > not found, (top tip: check spelling)",
485  cpy_field_name.c_str());
486  }
487  DofsByName::iterator miit = dofs->lower_bound(field_name);
488  if (miit == dofs->end()) {
489  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
490  "problem field < %s > not found, (top tip: check spelling)",
491  field_name.c_str());
492  }
493  DofsByName::iterator hi_miit = dofs->upper_bound(field_name);
494  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
495  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
496  "fields have to have same space");
497  }
498  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
499  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
500  "fields have to have same rank");
501  }
502  switch (scatter_mode) {
503  case SCATTER_REVERSE: {
504  Vec V_glob;
505  VecScatter ctx;
506  CHKERR VecScatterCreateToAll(V, &ctx, &V_glob);
507  CHKERR VecScatterBegin(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
508  CHKERR VecScatterEnd(ctx, V, V_glob, INSERT_VALUES, SCATTER_FORWARD);
509  int size;
510  CHKERR VecGetSize(V_glob, &size);
511  if (size != nb_dofs)
512  SETERRQ(
513  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
514  "data inconsistency: nb. of dofs and declared nb. dofs in database");
515  PetscScalar *array;
516  VecGetArray(V_glob, &array);
517  bool alpha = true;
518  switch (mode) {
519  case INSERT_VALUES:
520  break;
521  case ADD_VALUES:
522  alpha = false;
523  break;
524  default:
525  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
526  }
527  for (; miit != hi_miit; miit++) {
528  if ((*miit)->getPetscGlobalDofIdx() >= size) {
529  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
530  "data inconsistency: nb. of dofs and declared nb. dofs in "
531  "database");
532  }
533  DofEntity_multiIndex::index<
534  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator diiiit;
535  diiiit =
536  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
537  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
538  (*miit)->getEntDofIdx()));
539  if (diiiit ==
540  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
541  EntityHandle ent = (*miit)->getEnt();
542  rval = const_cast<moab::Interface &>(m_field.get_moab())
543  .add_entities((*cpy_fit)->getMeshset(), &ent, 1);
545  // create field moabent
546  ApproximationOrder order = (*miit)->getMaxOrder();
547  std::pair<FieldEntity_multiIndex::iterator, bool> p_e_miit;
548  try {
549  boost::shared_ptr<FieldEntity> moabent(
550  new FieldEntity(*cpy_fit, (*miit)->getRefEntityPtr()));
551  p_e_miit =
552  const_cast<FieldEntity_multiIndex *>(field_ents)->insert(moabent);
553  } catch (const std::exception &ex) {
554  std::ostringstream ss;
555  ss << "throw in method: " << ex.what() << " at line " << __LINE__
556  << " in file " << __FILE__ << std::endl;
557  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, ss.str().c_str());
558  }
559  if ((*p_e_miit.first)->getMaxOrder() < order) {
560  bool success =
561  const_cast<FieldEntity_multiIndex *>(field_ents)
562  ->modify(p_e_miit.first, FieldEntity_change_order(order));
563  if (!success)
564  SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
565  "modification unsuccessful");
566  }
567  // create field moabdof
568  DofEntityByNameAndEnt::iterator hi_diit, diit;
569  diit = dofs_ptr->get<Composite_Name_And_Ent_mi_tag>().lower_bound(
570  boost::make_tuple(field_name, (*miit)->getEnt()));
571  hi_diit = dofs_ptr->get<Composite_Name_And_Ent_mi_tag>().upper_bound(
572  boost::make_tuple(field_name, (*miit)->getEnt()));
573  for (; diit != hi_diit; diit++) {
574  boost::shared_ptr<DofEntity> mdof =
575  boost::shared_ptr<DofEntity>(new DofEntity(
576  *(p_e_miit.first), (*diit)->getDofOrder(),
577  (*diit)->getDofCoeffIdx(), (*diit)->getEntDofIdx()));
578  std::pair<DofEntity_multiIndex::iterator, bool> cpy_p_diit;
579  cpy_p_diit =
580  const_cast<DofEntity_multiIndex *>(dofs_ptr)->insert(mdof);
581  if (cpy_p_diit.second) {
582  bool success = const_cast<DofEntity_multiIndex *>(dofs_ptr)->modify(
583  cpy_p_diit.first, DofEntity_active_change(true));
584  if (!success)
585  SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
586  "modification unsuccessful");
587  }
588  }
589  diiiit =
590  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
591  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
592  (*miit)->getEntDofIdx()));
593  if (diiiit ==
594  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end())
595  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
596  "data inconsistency");
597  }
598  if (alpha)
599  (*diiiit)->getFieldData() = 0;
600  (*diiiit)->getFieldData() += array[(*miit)->getPetscGlobalDofIdx()];
601  if (dEbug) {
602  std::ostringstream ss;
603  ss << *(*diiiit) << "set " << array[(*miit)->getPetscGlobalDofIdx()]
604  << std::endl;
605  PetscPrintf(PETSC_COMM_WORLD, ss.str().c_str());
606  }
607  }
608  CHKERR VecRestoreArray(V_glob, &array);
609  CHKERR VecDestroy(&V_glob);
610  CHKERR VecScatterDestroy(&ctx);
611  } break;
612  case SCATTER_FORWARD: {
613  for (; miit != hi_miit; miit++) {
614  DofEntity_multiIndex::index<
615  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator diiiit;
616  diiiit =
617  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
618  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
619  (*miit)->getEntDofIdx()));
620  if (diiiit ==
621  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
622  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
623  "no data to fill the vector (top tip: you want scatter forward "
624  "or scatter reverse?)");
625  }
626  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
627  (*diiiit)->getFieldData(), mode);
628  }
629  CHKERR VecAssemblyBegin(V);
630  CHKERR VecAssemblyEnd(V);
631  } break;
632  default:
633  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
634  }
636 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:536
NumeredDofEntity_multiIndex::index< FieldName_mi_tag >::type NumeredDofEntityByFieldName
Numbered DoF multi-index by field name.
virtual moab::Interface & get_moab()=0
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.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
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.
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
virtual MoFEMErrorCode get_dofs(const DofEntity_multiIndex **dofs_ptr) const =0
Get dofs multi index.
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< DofEntity, UId, &DofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Composite_Ent_and_ShortId_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt >, const_mem_fun< DofEntity, ShortId, &DofEntity::getNonNonuniqueShortId > > >, 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, 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 > > >, ordered_non_unique< tag< Composite_Name_Ent_Order_And_CoeffIdx_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, ApproximationOrder, &DofEntity::getDofOrder >, const_mem_fun< DofEntity, FieldCoefficientsNumber, &DofEntity::getDofCoeffIdx > > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
#define CHKERR
Inline error check.
Definition: definitions.h:617
virtual MoFEMErrorCode get_field_ents(const FieldEntity_multiIndex **field_ents) const =0
Get field multi index.
int ApproximationOrder
Approximation on the entity.
Definition: Common.hpp:131
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
virtual MoFEMErrorCode get_fields(const Field_multiIndex **fields_ptr) const =0
Get fields multi-index from database.
int DofIdx
Index of DOF.
Definition: Common.hpp:126

◆ 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)collective - need tu be run on all processors in communicator.

Deprecated:
VecManager
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 638 of file VecManager.cpp.

641  {
642  const MoFEM::Interface &m_field = cOre;
643  const Problem *problem_ptr;
645  CHKERR m_field.get_problem(name, &problem_ptr);
646  CHKERR setOtherGlobalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
647  V, mode, scatter_mode);
649 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
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:451

◆ 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 325 of file VecManager.cpp.

328  {
329  const MoFEM::Interface &m_field = cOre;
330  const Field_multiIndex *fields_ptr;
331  const DofEntity_multiIndex *dofs_ptr;
333  CHKERR m_field.get_fields(&fields_ptr);
334  CHKERR m_field.get_dofs(&dofs_ptr);
335  typedef NumeredDofEntity_multiIndex::index<FieldName_mi_tag>::type DofsByName;
336  DofsByName *dofs;
337  switch (rc) {
338  case ROW:
339  dofs = const_cast<DofsByName *>(
340  &problem_ptr->numeredDofsRows->get<FieldName_mi_tag>());
341  break;
342  case COL:
343  dofs = const_cast<DofsByName *>(
344  &problem_ptr->numeredDofsCols->get<FieldName_mi_tag>());
345  break;
346  default:
347  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
348  }
349  Field_multiIndex::index<FieldName_mi_tag>::type::iterator cpy_fit =
350  fields_ptr->get<FieldName_mi_tag>().find(cpy_field_name);
351  if (cpy_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
352  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
353  "cpy field < %s > not found, (top tip: check spelling)",
354  cpy_field_name.c_str());
355  }
356  DofsByName::iterator miit = dofs->lower_bound(field_name);
357  if (miit == dofs->end()) {
358  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
359  "cpy field < %s > not found, (top tip: check spelling)",
360  field_name.c_str());
361  }
362  DofsByName::iterator hi_miit = dofs->upper_bound(field_name);
363  if ((*miit)->getSpace() != (*cpy_fit)->getSpace()) {
364  SETERRQ4(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
365  "fields have to have same space (%s) %s != (%s) %s",
366  (*miit)->getName().c_str(), FieldSpaceNames[(*miit)->getSpace()],
367  cpy_field_name.c_str(), FieldSpaceNames[(*cpy_fit)->getSpace()]);
368  }
369  if ((*miit)->getNbOfCoeffs() != (*cpy_fit)->getNbOfCoeffs()) {
370  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
371  "fields have to have same rank");
372  }
373  switch (scatter_mode) {
374  case SCATTER_REVERSE: {
375  bool alpha = true;
376  switch (mode) {
377  case INSERT_VALUES:
378  break;
379  case ADD_VALUES:
380  alpha = false;
381  break;
382  default:
383  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
384  }
385  PetscScalar *array;
386  VecGetArray(V, &array);
387  for (; miit != hi_miit; miit++) {
388  if (!miit->get()->getHasLocalIndex())
389  continue;
390  DofEntity_multiIndex::index<
391  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator diiiit;
392  diiiit =
393  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
394  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
395  (*miit)->getEntDofIdx()));
396  if (diiiit ==
397  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
398  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
399  "equivalent dof does not exist, use "
400  "set_other_global_ghost_vector to create dofs entries");
401  }
402  if (alpha) {
403  (*diiiit)->getFieldData() = array[(*miit)->getPetscLocalDofIdx()];
404  } else {
405  (*diiiit)->getFieldData() += array[(*miit)->getPetscLocalDofIdx()];
406  }
407  }
408  CHKERR VecRestoreArray(V, &array);
409  } break;
410  case SCATTER_FORWARD: {
411  for (; miit != hi_miit; miit++) {
412  if (!miit->get()->getHasLocalIndex())
413  continue;
414  DofEntity_multiIndex::index<
415  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator diiiit;
416  diiiit =
417  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().find(
418  boost::make_tuple(cpy_field_name, (*miit)->getEnt(),
419  (*miit)->getEntDofIdx()));
420  if (diiiit ==
421  dofs_ptr->get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().end()) {
422  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
423  "no data to fill the vector (top tip: you want scatter forward "
424  "or scatter reverse?)");
425  }
426  CHKERR VecSetValue(V, (*miit)->getPetscGlobalDofIdx(),
427  (*diiiit)->getFieldData(), mode);
428  }
429  CHKERR VecAssemblyBegin(V);
430  CHKERR VecAssemblyEnd(V);
431  } break;
432  default:
433  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "not implemented");
434  }
436 }
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.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
virtual MoFEMErrorCode get_dofs(const DofEntity_multiIndex **dofs_ptr) const =0
Get dofs multi index.
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< DofEntity, UId, &DofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Composite_Ent_and_ShortId_mi_tag >, composite_key< DofEntity, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt >, const_mem_fun< DofEntity, ShortId, &DofEntity::getNonNonuniqueShortId > > >, 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, 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 > > >, ordered_non_unique< tag< Composite_Name_Ent_Order_And_CoeffIdx_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, ApproximationOrder, &DofEntity::getDofOrder >, const_mem_fun< DofEntity, FieldCoefficientsNumber, &DofEntity::getDofCoeffIdx > > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
static const char *const FieldSpaceNames[]
Definition: definitions.h:184
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
virtual MoFEMErrorCode get_fields(const Field_multiIndex **fields_ptr) const =0
Get fields multi-index from database.

◆ 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 438 of file VecManager.cpp.

441  {
442  const MoFEM::Interface &m_field = cOre;
443  const Problem *problem_ptr;
445  CHKERR m_field.get_problem(name, &problem_ptr);
446  CHKERR setOtherLocalGhostVector(problem_ptr, field_name, cpy_field_name, rc,
447  V, mode, scatter_mode);
449 }
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:325
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443

◆ 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 59 of file VecManager.cpp.

60  {
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 = 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 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
NumeredDofEntity_multiIndex::index< PetscLocalIdx_mi_tag >::type NumeredDofEntityByLocalIdx
Numbered DoF multi-index by local index.
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
int DofIdx
Index of DOF.
Definition: Common.hpp:126

◆ VecCreateGhost()

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::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

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

Definition at line 130 of file DeprecatedCoreInterface.cpp.

132  {
133  return getInterface<VecManager>()->vecCreateGhost(name, rc, V);
134 }

◆ 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 36 of file VecManager.cpp.

37  {
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 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
int DofIdx
Index of DOF.
Definition: Common.hpp:126

◆ VecCreateSeq()

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::VecCreateSeq ( const std::string &  name,
RowColData  rc,
Vec *  V 
) const

create local vector for problem

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

Definition at line 125 of file DeprecatedCoreInterface.cpp.

127  {
128  return getInterface<VecManager>()->vecCreateSeq(name, rc, V);
129 }

◆ 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 102 of file VecManager.cpp.

106  {
107  const MoFEM::Interface &m_field = cOre;
109  std::vector<int> idx(0), idy(0);
110  CHKERR m_field.getInterface<ISManager>()
111  ->isCreateFromProblemFieldToOtherProblemField(
112  x_problem, x_field_name, x_rc, y_problem, y_field_name, y_rc, idx,
113  idy);
114  IS ix, iy;
115  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
116  PETSC_USE_POINTER, &ix);
117  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
118  PETSC_USE_POINTER, &iy);
119  CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
120  CHKERR ISDestroy(&ix);
121  CHKERR ISDestroy(&iy);
123 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443

◆ 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 form one to another problem (collective)

Parameters
xinvector
x_probleproblem name
yinvector
y_problemproblem name
Return values
newctxscatter

Definition at line 125 of file VecManager.cpp.

127  {
128  const MoFEM::Interface &m_field = cOre;
130  std::vector<int> idx(0), idy(0);
131  CHKERR m_field.getInterface<ISManager>()->isCreateFromProblemToOtherProblem(
132  x_problem, x_rc, y_problem, y_rc, idx, idy);
133  IS ix, iy;
134  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idx.size(), &idx[0],
135  PETSC_USE_POINTER, &ix);
136  CHKERR ISCreateGeneral(PETSC_COMM_WORLD, idy.size(), &idy[0],
137  PETSC_USE_POINTER, &iy);
138  if (dEbug) {
139  ISView(ix, PETSC_VIEWER_STDOUT_WORLD);
140  ISView(iy, PETSC_VIEWER_STDOUT_WORLD);
141  }
142  CHKERR ::VecScatterCreate(xin, ix, yin, iy, newctx);
143  CHKERR ISDestroy(&ix);
144  CHKERR ISDestroy(&iy);
146 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
const MoFEM::Interface & cOre
Definition: VecManager.hpp:39
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443

◆ VecScatterCreate() [1/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::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,
int  verb = -1 
) const

create scatter for vectors form one to another problem (collective)User specify what name of field on one problem is scattered to another.

Deprecated:
use VecManager
Parameters
xinvector
x_probleproblem name
x_fieldname
yinvector
y_problemproblem name
y_field_name
Return values
newctxscatter

Definition at line 166 of file DeprecatedCoreInterface.cpp.

170  {
171  return getInterface<VecManager>()->vecScatterCreate(
172  xin, x_problem, x_field_name, x_rc, yin, y_problem, y_field_name, y_rc,
173  newctx);
174 }

◆ VecScatterCreate() [2/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::VecScatterCreate ( Vec  xin,
const std::string &  x_problem,
RowColData  x_rc,
Vec  yin,
const std::string &  y_problem,
RowColData  y_rc,
VecScatter *  newctx,
int  verb = -1 
) const

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

Deprecated:
use VecManager
Parameters
xinvector
x_probleproblem name
yinvector
y_problemproblem name
Return values
newctxscatter

Definition at line 175 of file DeprecatedCoreInterface.cpp.

178  {
179  return getInterface<VecManager>()->vecScatterCreate(xin, x_problem, x_rc, yin,
180  y_problem, y_rc, newctx);
181 }