v0.10.0
Classes | Functions
Field Basic Algebra

Basic algebraic operation on fields. More...

Collaboration diagram for Field Basic Algebra:

Classes

struct  MoFEM::FieldBlas
 Basic algebra on fields. More...
 

Functions

MoFEMErrorCode MoFEM::FieldBlas::fieldLambda (TwoFieldFunction lambda, const std::string field_name_x, const std::string field_name_y, bool error_if_missing=false, bool creat_if_missing=false)
 filed lambdaDo calculation on two fields and save result to field fy More...
 
MoFEMErrorCode MoFEM::FieldBlas::fieldAxpy (const double alpha, const std::string field_name_x, const std::string field_name_y, bool error_if_missing=false, bool creat_if_missing=false)
 axpy fieldsfield_y = field_y + alpha*field_x More...
 
MoFEMErrorCode MoFEM::FieldBlas::fieldCopy (const double alpha, const std::string field_name_x, const std::string field_name_y, bool error_if_missing=false, bool creat_if_missing=false)
 copy and scale fieldsfield_y = alpha*field_x More...
 
MoFEMErrorCode MoFEM::FieldBlas::setVertexDofs (VertexCoordsFunction lambda, const std::string field_name, Range *verts=nullptr)
 Set DOFs on vertices using user functionExample: More...
 
MoFEMErrorCode MoFEM::FieldBlas::setField (const double val, const EntityType type, const std::string field_name)
 scale field More...
 
MoFEMErrorCode MoFEM::FieldBlas::setField (const double val, const EntityType type, const Range &ents, const std::string field_name)
 set fieldfield_y = val More...
 
MoFEMErrorCode MoFEM::FieldBlas::setField (const double val, const std::string field_name)
 set fieldfield_y = val More...
 
MoFEMErrorCode MoFEM::FieldBlas::fieldScale (const double alpha, const std::string field_name)
 set fieldfield_y = val More...
 

Field algebra

DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::field_axpy (const double alpha, const std::string &fiel_name_x, const std::string &field_name_y, bool error_if_missing=false, bool creat_if_missing=false)
 axpy fields More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::field_scale (const double alpha, const std::string &field_name)
 scale field More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_field (const double val, const EntityType type, const std::string &field_name)
 use FieldBlas More...
 
DEPRECATED MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_field (const double val, const EntityType type, const Range &ents, const std::string &field_name)
 set field More...
 

Detailed Description

Basic algebraic operation on fields.

Function Documentation

◆ field_axpy()

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::field_axpy ( const double  alpha,
const std::string &  fiel_name_x,
const std::string &  field_name_y,
bool  error_if_missing = false,
bool  creat_if_missing = false 
)

axpy fields

Deprecated:
use FieldBlas
Todo:
should be moved to independent interface, i.e. FieldAlgebra

field_y = field_y + alpha*field_x

Parameters
alpha
field_name_xname of field_x
field_name_yname of field_y
error_if_missingthrow error if entity/dof exist in field_x but not on field_y
create_if_missingcreat dof in field_y from fiedl_x if it is not database

Definition at line 70 of file DeprecatedCoreInterface.cpp.

73  {
74  return getInterface<FieldBlas>()->fieldAxpy(
75  alpha, field_name_x, field_name_y, error_if_missing, creat_if_missing);
76 }

◆ field_scale()

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::field_scale ( const double  alpha,
const std::string &  field_name 
)

scale field

Deprecated:
use FieldBlas
Todo:
should be moved to independent interface, i.e. FieldAlgebra
Parameters
alphais a scaling factor \field_name is a field name

Definition at line 89 of file DeprecatedCoreInterface.cpp.

90  {
91  return getInterface<FieldBlas>()->fieldScale(alpha, field_name);
92 }

◆ fieldAxpy()

MoFEMErrorCode MoFEM::FieldBlas::fieldAxpy ( const double  alpha,
const std::string  field_name_x,
const std::string  field_name_y,
bool  error_if_missing = false,
bool  creat_if_missing = false 
)

axpy fieldsfield_y = field_y + alpha*field_x

Parameters
alpha
field_name_xname of field_x
field_name_yname of field_y
error_if_missingthrow error if entity/dof exist in field_x but not on field_y
create_if_missingcreat dof in field_y from fiedl_x if it is not database
Examples
field_axpy_atom_test.cpp.

Definition at line 135 of file FieldBlas.cpp.

139  {
141  struct Axpy {
142  const double alpha;
143  Axpy(const double alpha) : alpha(alpha) {}
144  inline MoFEMErrorCode operator()(double &fy, const double fx) {
146  fy += alpha * fx;
148  }
149  };
150  CHKERR fieldLambda(Axpy(alpha), field_name_x, field_name_y, error_if_missing,
151  creat_if_missing);
153 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
MoFEMErrorCode fieldLambda(TwoFieldFunction lambda, const std::string field_name_x, const std::string field_name_y, bool error_if_missing=false, bool creat_if_missing=false)
filed lambdaDo calculation on two fields and save result to field fy
Definition: FieldBlas.cpp:38
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ fieldCopy()

MoFEMErrorCode MoFEM::FieldBlas::fieldCopy ( const double  alpha,
const std::string  field_name_x,
const std::string  field_name_y,
bool  error_if_missing = false,
bool  creat_if_missing = false 
)

copy and scale fieldsfield_y = alpha*field_x

Parameters
alpha
field_name_xname of field_x
field_name_yname of field_y
error_if_missingthrow error if entity/dof exist in field_x but not on field_y
create_if_missingcreat dof in field_y from fiedl_x if it is not database

Definition at line 155 of file FieldBlas.cpp.

159  {
161  struct Copy {
162  const double alpha;
163  Copy(const double alpha) : alpha(alpha) {}
164  inline MoFEMErrorCode operator()(double &fy, const double fx) {
166  fy = alpha * fx;
168  }
169  };
170  CHKERR fieldLambda(Copy(alpha), field_name_x, field_name_y, error_if_missing,
171  creat_if_missing);
173 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
MoFEMErrorCode fieldLambda(TwoFieldFunction lambda, const std::string field_name_x, const std::string field_name_y, bool error_if_missing=false, bool creat_if_missing=false)
filed lambdaDo calculation on two fields and save result to field fy
Definition: FieldBlas.cpp:38
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ fieldLambda()

MoFEMErrorCode MoFEM::FieldBlas::fieldLambda ( FieldBlas::TwoFieldFunction  lambda,
const std::string  field_name_x,
const std::string  field_name_y,
bool  error_if_missing = false,
bool  creat_if_missing = false 
)

filed lambdaDo calculation on two fields and save result to field fy

struct Axpy {
const double aLpha;
Axpy(const double alpha) : aLpha(alpha) {}
inline MoFEMErrorCode operator(double &fy, double fx) {
fy += Alpha * fx;
}
};
CHKERR m_fiel.getInterface<FieldBlas>()->fieldLambda(Axpy(aLpha),
field_name_x, field_name_y);
Parameters
functionf(double &x, double)
field_name_xname of field_x
field_name_yname of field_y
error_if_missingthrow error if entity/dof exist in field_x but not on field_y
create_if_missingcreat dof in field_y from field_x if it is not database

Definition at line 38 of file FieldBlas.cpp.

42  {
43  const MoFEM::Interface &m_field = cOre;
44  auto fields_ptr = m_field.get_fields();
45  auto field_ents = m_field.get_field_ents();
46  auto dofs_ptr = m_field.get_dofs();
48 
49  auto x_fit = fields_ptr->get<FieldName_mi_tag>().find(field_name_x);
50  if (x_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
51  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
52  "x field < %s > not found, (top tip: check spelling)",
53  field_name_x.c_str());
54  }
55  auto y_fit = fields_ptr->get<FieldName_mi_tag>().find(field_name_y);
56  if (y_fit == fields_ptr->get<FieldName_mi_tag>().end()) {
57  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
58  "y field < %s > not found, (top tip: check spelling)",
59  field_name_y.c_str());
60  }
61  if ((*x_fit)->getSpace() != (*y_fit)->getSpace()) {
62  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
63  "space for field < %s > and field <%s> are not compatible",
64  field_name_x.c_str(), field_name_y.c_str());
65  }
66  if ((*x_fit)->getNbOfCoeffs() != (*y_fit)->getNbOfCoeffs()) {
67  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
68  "rank for field < %s > and field <%s> are not compatible",
69  field_name_x.c_str(), field_name_y.c_str());
70  }
71 
72  typedef multi_index_container<
73  boost::shared_ptr<DofEntity>,
74  indexed_by<
75 
76  hashed_non_unique<
77  tag<Composite_Ent_Order_And_CoeffIdx_mi_tag>,
78  composite_key<
79 
80  DofEntity,
81  const_mem_fun<DofEntity, EntityHandle, &DofEntity::getEnt>,
82  const_mem_fun<DofEntity, ApproximationOrder,
84  const_mem_fun<DofEntity, FieldCoefficientsNumber,
86 
87  >>
88 
89  >>
90  DofEntity_multiIndex_composite_view;
91 
92  auto dof_lo_for_view = dofs_ptr->get<Unique_mi_tag>().lower_bound(
93  FieldEntity::getLoBitNumberUId((*y_fit)->getBitNumber()));
94  auto dof_hi_for_view = dofs_ptr->get<Unique_mi_tag>().upper_bound(
95  FieldEntity::getHiBitNumberUId((*y_fit)->getBitNumber()));
96 
97  DofEntity_multiIndex_composite_view dof_composite_view;
98  dof_composite_view.insert(dof_lo_for_view, dof_hi_for_view);
99 
100  auto x_eit = field_ents->get<Unique_mi_tag>().lower_bound(
101  FieldEntity::getLoBitNumberUId((*x_fit)->getBitNumber()));
102  auto x_eit_hi = field_ents->get<Unique_mi_tag>().upper_bound(
103  FieldEntity::getHiBitNumberUId((*x_fit)->getBitNumber()));
104  for (; x_eit != x_eit_hi; x_eit++) {
105  int nb_dofs_on_x_entity = (*x_eit)->getNbDofsOnEnt();
106  VectorAdaptor field_data = (*x_eit)->getEntFieldData();
107  for (int dd = 0; dd != nb_dofs_on_x_entity; ++dd) {
108  ApproximationOrder dof_order = (*x_eit)->getDofOrderMap()[dd];
109  FieldCoefficientsNumber dof_rank = dd % (*x_eit)->getNbOfCoeffs();
110  FieldData data = field_data[dd];
111  auto dit = dof_composite_view.find(
112  boost::make_tuple((*x_eit)->getEnt(), dof_order, dof_rank));
113  if (dit == dof_composite_view.end()) {
114  if (creat_if_missing) {
115  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
116  "not yet implemented");
117  } else {
118  if (error_if_missing) {
119  std::ostringstream ss;
120  ss << "dof on ent " << (*x_eit)->getEnt() << " order " << dof_order
121  << " rank " << dof_rank << " does not exist";
122  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
123  ss.str().c_str());
124  } else {
125  continue;
126  }
127  }
128  }
129  CHKERR lambda((*dit)->getFieldData(),data);
130  }
131  }
133 }
Deprecated interface functions.
FieldCoefficientsNumber getDofCoeffIdx() const
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:37
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.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
const MoFEM::Interface & cOre
Definition: FieldBlas.hpp:39
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:109
#define CHKERR
Inline error check.
Definition: definitions.h:604
double FieldData
Field data type.
Definition: Types.hpp:36
ApproximationOrder getDofOrder() const
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
int FieldCoefficientsNumber
Number of field coefficients.
Definition: Types.hpp:38
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
static UId getHiBitNumberUId(const FieldBitNumber bit_number)

◆ fieldScale()

MoFEMErrorCode MoFEM::FieldBlas::fieldScale ( const double  alpha,
const std::string  field_name 
)

set fieldfield_y = val

Parameters
val
entitytype
onenties
field_name
Examples
field_axpy_atom_test.cpp.

Definition at line 322 of file FieldBlas.cpp.

323  {
324  const MoFEM::Interface &m_field = cOre;
325  auto fields_ptr = m_field.get_fields();
326  auto dofs_ptr = m_field.get_dofs();
328 
329  auto fit = fields_ptr->get<FieldName_mi_tag>().find(field_name);
330  if (fit == fields_ptr->get<FieldName_mi_tag>().end())
331  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
332  " field < %s > not found, (top tip: check spelling)",
333  field_name.c_str());
334 
335  auto dit = dofs_ptr->get<Unique_mi_tag>().lower_bound(
336  FieldEntity::getLoBitNumberUId((*fit)->getBitNumber()));
337  auto hi_dit = dofs_ptr->get<Unique_mi_tag>().upper_bound(
338  FieldEntity::getHiBitNumberUId((*fit)->getBitNumber()));
339  for (; dit != hi_dit; dit++)
340  (*dit)->getFieldData() *= alpha;
341 
343 }
Deprecated interface functions.
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
const MoFEM::Interface & cOre
Definition: FieldBlas.hpp:39
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
static UId getHiBitNumberUId(const FieldBitNumber bit_number)

◆ set_field() [1/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_field ( const double  val,
const EntityType  type,
const std::string &  field_name 
)

use FieldBlas

set field

Todo:
should be moved to independent interface, i.e. FieldAlgebra

field_y = val

Parameters
val
entitytype
field_name

Definition at line 78 of file DeprecatedCoreInterface.cpp.

79  {
80  return getInterface<FieldBlas>()->setField(val, type, field_name);
81 }

◆ set_field() [2/2]

MoFEMErrorCode MoFEM::DeprecatedCoreInterface::set_field ( const double  val,
const EntityType  type,
const Range &  ents,
const std::string &  field_name 
)

set field

Deprecated:
use FieldBlas
Todo:
should be moved to independent interface, i.e. FieldAlgebra

field_y = val

Parameters
val
entitytype
onenties
field_name

Definition at line 83 of file DeprecatedCoreInterface.cpp.

85  {
86  return getInterface<FieldBlas>()->setField(val, type, ents, field_name);
87 }

◆ setField() [1/3]

MoFEMErrorCode MoFEM::FieldBlas::setField ( const double  val,
const EntityType  type,
const std::string  field_name 
)

scale field

Parameters
valis a set parameter \field_name is a field name
Examples
field_axpy_atom_test.cpp.

Definition at line 243 of file FieldBlas.cpp.

244  {
245  const MoFEM::Interface &m_field = cOre;
246  auto dofs_ptr = m_field.get_dofs();
248 
249  const auto bit_number = m_field.get_field_bit_number(field_name);
250  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
251  bit_number, get_id_for_min_type(type));
252  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
253  bit_number, get_id_for_max_type(type));
254 
255  auto dit = dofs_ptr->get<Unique_mi_tag>().lower_bound(lo_uid);
256  auto hi_dit = dofs_ptr->get<Unique_mi_tag>().upper_bound(hi_uid);
257 
258  for (; dit != hi_dit; dit++)
259  (*dit)->getFieldData() = val;
260 
262 }
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
Deprecated interface functions.
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
const MoFEM::Interface & cOre
Definition: FieldBlas.hpp:39
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number

◆ setField() [2/3]

MoFEMErrorCode MoFEM::FieldBlas::setField ( const double  val,
const EntityType  type,
const Range &  ents,
const std::string  field_name 
)

set fieldfield_y = val

Parameters
val
entitytype
field_name

Definition at line 264 of file FieldBlas.cpp.

266  {
267  const MoFEM::Interface &m_field = cOre;
268  auto dofs_ptr = m_field.get_dofs();
270 
271  const auto bit_number = m_field.get_field_bit_number(field_name);
272  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
273  bit_number, get_id_for_min_type(type));
274  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
275  bit_number, get_id_for_max_type(type));
276 
277  auto dit = dofs_ptr->get<Unique_mi_tag>().lower_bound(lo_uid);
278  auto hi_dit = dofs_ptr->get<Unique_mi_tag>().upper_bound(hi_uid);
279 
280  EntityHandle ent, last = 0;
281  bool cont = true;
282  for (; dit != hi_dit; dit++) {
283  ent = (*dit)->getEnt();
284  if (ent != last) {
285  if (ents.find(ent) == ents.end()) {
286  cont = true;
287  } else {
288  cont = false;
289  }
290  last = ent;
291  }
292  if (!cont)
293  (*dit)->getFieldData() = val;
294 
295  }
297 }
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
Deprecated interface functions.
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
const MoFEM::Interface & cOre
Definition: FieldBlas.hpp:39
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number

◆ setField() [3/3]

MoFEMErrorCode MoFEM::FieldBlas::setField ( const double  val,
const std::string  field_name 
)

set fieldfield_y = val

Parameters
val
field_name

Definition at line 299 of file FieldBlas.cpp.

300  {
301  const MoFEM::Interface &m_field = cOre;
302  auto fields_ptr = m_field.get_fields();
303  auto dofs_ptr = m_field.get_dofs();
305  auto fit = fields_ptr->get<FieldName_mi_tag>().find(field_name);
306  if (fit == fields_ptr->get<FieldName_mi_tag>().end()) {
307  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
308  " field < %s > not found, (top tip: check spelling)",
309  field_name.c_str());
310  }
311 
312  auto dit = dofs_ptr->get<Unique_mi_tag>().lower_bound(
313  FieldEntity::getLoBitNumberUId((*fit)->getBitNumber()));
314  auto hi_dit = dofs_ptr->get<Unique_mi_tag>().upper_bound(
315  FieldEntity::getHiBitNumberUId((*fit)->getBitNumber()));
316  for (; dit != hi_dit; dit++) {
317  (*dit)->getFieldData() = val;
318  }
320 }
Deprecated interface functions.
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
const MoFEM::Interface & cOre
Definition: FieldBlas.hpp:39
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
static UId getHiBitNumberUId(const FieldBitNumber bit_number)

◆ setVertexDofs()

MoFEMErrorCode MoFEM::FieldBlas::setVertexDofs ( FieldBlas::VertexCoordsFunction  lambda,
const std::string  field_name,
Range *  verts = nullptr 
)

Set DOFs on vertices using user functionExample:

auto do_something = [&](VectorAdaptor &field_data, double *x,
double *y, double *z) {
field_data[0] = (*x);
field_data[1] = (*y);
field_data[2] = (*z);
};
CHKERR m_field.getInterface<FieldBlas>()->setVertexDofs(set_distance,
"DISP");
Note
Function works both ways, using it coordinates can be set from field.
Parameters
lambdafunction evaluating field at points
field_nameis a field name
vertspointer to vertices if null all vertices in the field are evaluated)
Examples
field_blas_set_vertex_dofs.cpp, and lesson6_radiation.cpp.

Definition at line 175 of file FieldBlas.cpp.

177  {
178  const MoFEM::Interface &m_field = cOre;
180 
181  EntityHandle meshset = m_field.get_field_meshset(field_name);
182  Range verts;
183  CHKERR m_field.get_moab().get_entities_by_type(meshset, MBVERTEX, verts,
184  true);
185  if (sub_verts)
186  verts = intersect(*sub_verts, verts);
187 
188  struct LambdaMethod : EntityMethod {
189  LambdaMethod(MoFEM::Interface &m_field, Range &verts,
191  : EntityMethod(), mField(m_field), verts(verts), lambda(lambda),
192  count(0), total(0) {}
193  MoFEMErrorCode preProcess() {
194  vit = verts.begin();
195  return 0;
196  }
197  MoFEMErrorCode operator()() {
199  if (*vit != entPtr->getEnt())
200  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
201  "Inconsistent entity %ld != %ld", *vit, entPtr->getEnt());
202  if(!count)
203  CHKERR mField.get_moab().coords_iterate(vit, verts.end(), xCoord,
204  yCoord, zCoord, count);
205  CHKERR lambda(entPtr->getEntFieldData(), xCoord, yCoord, zCoord);
206  ++xCoord;
207  ++yCoord;
208  ++zCoord;
209  ++vit;
210  ++total;
211  --count;
213  }
214  MoFEMErrorCode postProcess() {
216  if(total != verts.size())
217  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
218  "Inconsistent number of vertices in field meshset and in the "
219  "field %d != %d",
220  total, verts.size());
222  }
223 
224  private:
225  MoFEM::Interface &mField;
226  Range &verts;
228  int count;
229  int total;
230  Range::iterator vit;
231  double *xCoord;
232  double *yCoord;
233  double *zCoord;
234  };
235 
236  LambdaMethod lambda_method(const_cast<MoFEM::Interface &>(m_field), verts,
237  lambda);
238  CHKERR const_cast<MoFEM::Interface &>(m_field).loop_entities(
239  field_name, lambda_method, &verts, QUIET);
241 }
Deprecated interface functions.
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
boost::function< MoFEMErrorCode(VectorAdaptor &&field_data, double *xcoord, double *ycoord, double *zcoord)> VertexCoordsFunction
Definition: FieldBlas.hpp:124
const MoFEM::Interface & cOre
Definition: FieldBlas.hpp:39
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
#define CHKERR
Inline error check.
Definition: definitions.h:604
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415