31 : cOre(const_cast<
MoFEM::
Core &>(core)), dEbug(false) {}
35 const std::string field_name, Range *ents_ptr) {
43 "field < %s > not found, (top tip: check spelling)",
47 auto wrap_lambda_on_enties =
48 [
lambda](boost::shared_ptr<FieldEntity> field_entity_ptr) {
51 auto field_data = field_entity_ptr->getEntFieldData();
52 for (
auto &
v : field_data)
64 const std::string field_name,
75 "field < %s > not found, (top tip: check spelling)",
80 const auto bit_number = (*fit)->getBitNumber();
89 auto eit = field_ents->get<
Unique_mi_tag>().lower_bound(lo_uid);
91 auto eit_hi = field_ents->get<
Unique_mi_tag>().upper_bound(hi_uid);
93 for (; eit != eit_hi; ++eit) {
101 for (
auto p = ents_ptr->const_pair_begin();
p != ents_ptr->const_pair_end();
103 CHKERR loop_ents(
p->first,
p->second);
108 CHKERR loop_ents(get_id_for_min_type<MBVERTEX>(),
109 get_id_for_max_type<MBENTITYSET>());
117 const std::string field_name_y,
bool error_if_missing, Range *ents_ptr) {
125 "x field < %s > not found, (top tip: check spelling)",
126 field_name_x.c_str());
131 "y field < %s > not found, (top tip: check spelling)",
132 field_name_y.c_str());
134 if ((*x_fit)->getSpace() != (*y_fit)->getSpace()) {
136 "space for field < %s > and field <%s> are not compatible for "
138 field_name_x.c_str(), field_name_y.c_str());
140 if ((*x_fit)->getNbOfCoeffs() != (*y_fit)->getNbOfCoeffs()) {
143 "rank for field < %s > and field <%s> are not compatible for fieldblas",
144 field_name_x.c_str(), field_name_y.c_str());
147 auto wrap_lambda_on_enties =
148 [
lambda](boost::shared_ptr<FieldEntity> y_field_entity_ptr,
149 const boost::shared_ptr<FieldEntity> x_field_entity_ptr) {
152 auto x_field_data = x_field_entity_ptr->getEntFieldData();
153 auto y_field_data = y_field_entity_ptr->getEntFieldData();
154 const auto size_x = x_field_data.size();
155 const auto size_y = y_field_data.size();
158 for (;
dd != std::min(size_x, size_y); ++
dd)
159 y_field_data[
dd] =
lambda(y_field_data[
dd], x_field_data[
dd]);
160 for (;
dd < size_y; ++
dd)
161 y_field_data[
dd] = 0;
167 field_name_y, error_if_missing, ents_ptr);
173 const std::string field_name_x,
174 const std::string field_name_y,
175 bool error_if_missing, Range *ents_ptr) {
185 "x field < %s > not found, (top tip: check spelling)",
186 field_name_x.c_str());
191 "y field < %s > not found, (top tip: check spelling)",
192 field_name_y.c_str());
196 const auto x_bit_number = (*x_fit)->getBitNumber();
197 const auto y_bit_number = (*y_fit)->getBitNumber();
198 const auto y_eit_hi = field_ents->get<
Unique_mi_tag>().upper_bound(
204 const auto x_lo_uid =
206 const auto x_hi_uid =
210 auto x_eit = field_ents->get<
Unique_mi_tag>().lower_bound(x_lo_uid);
212 auto x_eit_hi = field_ents->get<
Unique_mi_tag>().upper_bound(x_hi_uid);
214 for (; x_eit != x_eit_hi;) {
217 (*y_fit)->getBitNumber(), (*x_eit)->getEnt());
218 auto y_eit = field_ents->get<
Unique_mi_tag>().find(y_lo_uid);
220 if (y_eit == field_ents->end()) {
222 if (error_if_missing) {
224 "Missing entity in y_field.");
232 if (x_eit == x_eit_hi)
234 if (y_eit == y_eit_hi)
236 if ((*y_eit)->getEnt() != (*x_eit)->getEnt())
256 for (
auto p = ents_ptr->const_pair_begin();
p != ents_ptr->const_pair_end();
258 CHKERR loop_ents(
p->first,
p->second);
263 CHKERR loop_ents(get_id_for_min_type<MBVERTEX>(),
264 get_id_for_max_type<MBENTITYSET>());
271 const std::string field_name_x,
272 const std::string field_name_y,
273 bool error_if_missing,
274 bool create_if_missing) {
275 if (create_if_missing)
277 <<
"Option create_if_missing is set to true but have no effect";
283 const std::string field_name_x,
284 const std::string field_name_y,
285 bool error_if_missing, Range *ents_ptr) {
287 auto axpy = [
alpha](
const double fy,
const double fx) {
288 return fy +
alpha * fx;
296 const std::string field_name_x,
297 const std::string field_name_y,
298 bool error_if_missing,
299 bool create_if_missing) {
300 return fieldAxpy(
alpha, field_name_x, field_name_y, error_if_missing);
304 const std::string field_name_x,
305 const std::string field_name_y,
306 bool error_if_missing, Range *ents_ptr) {
308 auto copy = [
alpha](
const double fy,
const double fx) {
return alpha * fx; };
315 const std::string field_name_x,
316 const std::string field_name_y,
317 bool error_if_missing,
318 bool create_if_missing) {
319 return fieldCopy(
alpha, field_name_x, field_name_y, error_if_missing);
323 const std::string field_name,
334 const std::string field_name,
341 CHKERR m_field.
get_moab().get_entities_by_type(meshset, MBVERTEX, verts,
344 verts = intersect(*sub_verts, verts);
350 count(0), total(0) {}
357 if (*vit != entPtr->getEnt())
359 "Inconsistent entity %ld != %ld", *vit, entPtr->getEnt());
361 CHKERR mField.get_moab().coords_iterate(vit, verts.end(), xCoord,
362 yCoord, zCoord, count);
363 CHKERR lambda(entPtr->getEntFieldData(), xCoord, yCoord, zCoord);
374 if (total != verts.size())
376 "Inconsistent number of vertices in field meshset and in the "
378 total, verts.size());
397 field_name, lambda_method, &verts,
QUIET);
402 const std::string field_name) {
413 auto dit = dofs_ptr->get<
Unique_mi_tag>().lower_bound(lo_uid);
414 auto hi_dit = dofs_ptr->get<
Unique_mi_tag>().upper_bound(hi_uid);
416 for (; dit != hi_dit; dit++)
417 (*dit)->getFieldData() = val;
424 const std::string field_name) {
435 auto dit = dofs_ptr->get<
Unique_mi_tag>().lower_bound(lo_uid);
436 auto hi_dit = dofs_ptr->get<
Unique_mi_tag>().upper_bound(hi_uid);
440 for (; dit != hi_dit; dit++) {
441 ent = (*dit)->getEnt();
443 if (ents.find(ent) == ents.end()) {
451 (*dit)->getFieldData() = val;
457 const std::string field_name) {
465 " field < %s > not found, (top tip: check spelling)",
473 for (; dit != hi_dit; dit++) {
474 (*dit)->getFieldData() = val;
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
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.
#define MOFEM_LOG(channel, severity)
Log.
double v
phase velocity of light in medium (cm/ns)
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)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
Data structure to exchange data between mofem and User Loop Methods on entities.
DEPRECATED MoFEMErrorCode fieldLambda(TwoFieldFunctionOnValues lambda, const std::string field_name_x, const std::string field_name_y, bool error_if_missing, bool create_if_missing)
MoFEMErrorCode fieldAxpy(const double alpha, const std::string field_name_x, const std::string field_name_y, bool error_if_missing=false, Range *ents_ptr=nullptr)
axpy fields
boost::function< MoFEMErrorCode(boost::shared_ptr< FieldEntity > field_entity_ptr)> OneFieldFunctionOnEntities
MoFEMErrorCode fieldLambdaOnEntities(OneFieldFunctionOnEntities lambda, const std::string field_name, Range *ents_ptr=nullptr)
filed lambda
boost::function< double(const double y_field_value_reference, const double x_field_value)> TwoFieldFunctionOnValues
MoFEMErrorCode setVertexDofs(VertexCoordsFunction lambda, const std::string field_name, Range *verts=nullptr)
Set DOFs on vertices using user function.
boost::function< MoFEMErrorCode(VectorAdaptor &&field_data, double *xcoord, double *ycoord, double *zcoord)> VertexCoordsFunction
MoFEMErrorCode setField(const double val, const EntityType type, const std::string field_name)
boost::function< double(const double field_val)> OneFieldFunctionOnValues
function to set a field value
MoFEMErrorCode fieldScale(const double alpha, const std::string field_name, Range *ents_ptr=nullptr)
field scale
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode fieldLambdaOnValues(OneFieldFunctionOnValues lambda, const std::string field_name, Range *ents_ptr=nullptr)
filed lambda
FieldBlas(const MoFEM::Core &core)
const MoFEM::Interface & cOre
boost::function< MoFEMErrorCode(boost::shared_ptr< FieldEntity > y_field_entity_ptr, const boost::shared_ptr< FieldEntity > x_field_entity_ptr)> TwoFieldFunctionOnEntities
MoFEMErrorCode fieldCopy(const double alpha, const std::string field_name_x, const std::string field_name_y, bool error_if_missing=false, Range *ents_ptr=nullptr)
copy and scale fields
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
static UId getLoLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
static UId getHiLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
MultiIndex Tag for field name.
base class for all interface classes