v0.9.1
SeriesRecorder.cpp
Go to the documentation of this file.
1 /** \file SeriesRecorder.cpp
2  * \brief Record and save data series in time or load steps
3  */
4 
5 /* MoFEM is free software: you can redistribute it and/or modify it under
6  * the terms of the GNU Lesser General Public License as published by the
7  * Free Software Foundation, either version 3 of the License, or (at your
8  * option) any later version.
9  *
10  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13  * License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
17  */
18 
19 #define SeriesRecorderFunctionBegin \
20  MoFEMFunctionBegin; \
21  MOFEM_LOG_CHANNEL("WORLD"); \
22  MOFEM_LOG_CHANNEL("SYNC"); \
23  MOFEM_LOG_FUNCTION(); \
24  MOFEM_LOG_TAG("SYNC", "SeriesRecorder"); \
25  MOFEM_LOG_TAG("WORLD", "SeriesRecorder")
26 
27 namespace MoFEM {
28 
29 // const static int debug = 1;
30 
32  UnknownInterface **iface) const {
34  *iface = NULL;
35  if (uuid == IDD_MOFEMSeriesRecorder) {
36  *iface = const_cast<SeriesRecorder *>(this);
38  }
39  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
41 }
42 
44  : cOre(const_cast<Core &>(core)) {}
45 
47 
48  Interface &m_field = cOre;
49  moab::Interface &moab = m_field.get_moab();
51  const int def_val_len = 0;
52  CHKERR moab.tag_get_handle(
53  "_SeriesName", def_val_len, MB_TYPE_OPAQUE, th_SeriesName,
54  MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
56 }
57 
60  sEries.clear();
61  seriesSteps.clear();
63 }
64 
66  Interface &m_field = cOre;
67  moab::Interface &moab = m_field.get_moab();
69  Range meshsets;
70  const int def_val_len = 0;
71  CHKERR moab.tag_get_handle(
72  "_SeriesName", def_val_len, MB_TYPE_OPAQUE, th_SeriesName,
73  MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
74 
75  CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshsets, false);
76  for (Range::iterator mit = meshsets.begin(); mit != meshsets.end(); mit++) {
77  const void *tagName;
78  int tagNameSize;
79  rval = moab.tag_get_by_ptr(th_SeriesName, &*mit, 1, (const void **)&tagName,
80  &tagNameSize);
81  if (rval == MB_SUCCESS) {
82  std::pair<Series_multiIndex::iterator, bool> p =
83  sEries.insert(FieldSeries(moab, *mit));
84  if (verb > QUIET)
85  MOFEM_LOG("SYNC", LogManager::SeverityLevel::inform)
86  << "read series " << *p.first;
87  }
88  }
89  // //build series steps
90  for (Series_multiIndex::iterator sit = sEries.begin(); sit != sEries.end();
91  sit++) {
92  int nb_steps;
93  CHKERR sit->get_nb_steps(moab, nb_steps);
94  int ss = 0;
95  for (; ss < nb_steps; ss++) {
96  std::pair<SeriesStep_multiIndex::iterator, bool> p =
97  seriesSteps.insert(FieldSeriesStep(moab, &*sit, ss));
98  if (verb > QUIET)
99  MOFEM_LOG("SYNC", LogManager::SeverityLevel::inform)
100  << "add series step " << *p.first;
101  }
102  }
104 }
105 
107 SeriesRecorder::add_series_recorder(const std::string &series_name) {
108 
109  Interface &m_field = cOre;
110  moab::Interface &moab = m_field.get_moab();
112  EntityHandle meshset;
113  CHKERR moab.create_meshset(MESHSET_SET | MESHSET_TRACK_OWNER, meshset);
114  void const *tag_data[] = {series_name.c_str()};
115  int tag_sizes[1];
116  tag_sizes[0] = series_name.size();
117  CHKERR moab.tag_set_by_ptr(th_SeriesName, &meshset, 1, tag_data, tag_sizes);
118  std::pair<Series_multiIndex::iterator, bool> p =
119  sEries.insert(FieldSeries(moab, meshset));
120  if (!p.second) {
121  SETERRQ1(PETSC_COMM_SELF, 1, "series recorder <%s> is already there",
122  series_name.c_str());
123  }
125 }
126 
128 SeriesRecorder::delete_recorder_series(const std::string &series_name) {
129 
130  Interface &m_field = cOre;
131  moab::Interface &moab = m_field.get_moab();
133  SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator ssit, hi_ssit;
134  ssit = seriesSteps.get<SeriesName_mi_tag>().lower_bound(series_name);
135  hi_ssit = seriesSteps.get<SeriesName_mi_tag>().upper_bound(series_name);
136  seriesSteps.get<SeriesName_mi_tag>().erase(ssit, hi_ssit);
137  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit;
138  sit = sEries.get<SeriesName_mi_tag>().find(series_name);
139  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
140  SETERRQ1(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
141  "series recorder <%s> not exist and can be deleted",
142  series_name.c_str());
143  }
144  EntityHandle series_meshset = sit->getMeshset();
145  CHKERR moab.tag_delete(sit->th_SeriesTime);
146  CHKERR moab.tag_delete(sit->th_SeriesDataHandles);
147  CHKERR moab.tag_delete(sit->th_SeriesDataUIDs);
148  CHKERR moab.tag_delete(sit->th_SeriesData);
149  sEries.get<SeriesName_mi_tag>().erase(sit);
150  std::vector<EntityHandle> contained;
151  CHKERR moab.get_contained_meshsets(series_meshset, contained);
152  CHKERR moab.remove_entities(series_meshset, &contained[0], contained.size());
153  CHKERR moab.delete_entities(&contained[0], contained.size());
154  CHKERR moab.delete_entities(&series_meshset, 1);
156 }
157 
158 MoFEMErrorCode SeriesRecorder::record_problem(const std::string &serie_name,
159  const Problem *problemPtr,
160  RowColData rc) {
162 
163  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
164  sEries.get<SeriesName_mi_tag>().find(serie_name);
165  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
166  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
167  "series recorder <%s> not exist", serie_name.c_str());
168  }
169  switch (rc) {
170  case ROW:
171  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(
172  problemPtr->numeredDofsRows->begin(),
173  problemPtr->numeredDofsRows->end());
174  break;
175  case COL:
176  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(
177  problemPtr->numeredDofsCols->begin(),
178  problemPtr->numeredDofsCols->end());
179  break;
180  default:
181  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
182  }
184 }
185 
186 MoFEMErrorCode SeriesRecorder::record_problem(const std::string &serie_name,
187  const std::string &problem_name,
188  RowColData rc) {
189 
190  Interface &m_field = cOre;
192  const Problem *problem_ptr;
193  CHKERR m_field.get_problem(problem_name, &problem_ptr);
194  CHKERR record_problem(serie_name, problem_ptr, rc);
196 }
197 
198 MoFEMErrorCode SeriesRecorder::record_field(const std::string &serie_name,
199  const std::string &field_name,
200  const BitRefLevel &bit,
201  const BitRefLevel &mask) {
202 
203  Interface &m_field = cOre;
204  const DofEntity_multiIndex *dofs_ptr;
206  CHKERR m_field.get_dofs(&dofs_ptr);
207  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
208  sEries.get<SeriesName_mi_tag>().find(serie_name);
209  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
210  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
211  serie_name.c_str());
212  }
213  DofEntityByFieldName::iterator dit =
214  dofs_ptr->get<FieldName_mi_tag>().lower_bound(field_name);
215  if (dit == dofs_ptr->get<FieldName_mi_tag>().end()) {
216  SETERRQ1(PETSC_COMM_SELF, 1, "field <%s> not exist", field_name.c_str());
217  }
218  DofEntityByFieldName::iterator hi_dit =
219  dofs_ptr->get<FieldName_mi_tag>().upper_bound(field_name);
220  for (; dit != hi_dit; dit++) {
221  const BitRefLevel &dof_bit = (*dit)->getBitRefLevel();
222  if ((dof_bit & mask) != dof_bit)
223  continue;
224  if ((dof_bit & bit).any()) {
225  EntityHandle ent = (*dit)->getEnt();
226  ShortId uid = (*dit)->getNonNonuniqueShortId();
227  FieldData val = (*dit)->getFieldData();
228  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(ent, uid, val);
229  }
230  }
232 }
233 
234 MoFEMErrorCode SeriesRecorder::record_begin(const std::string &serie_name) {
235 
237  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
238  sEries.get<SeriesName_mi_tag>().find(serie_name);
239  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
240  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
241  serie_name.c_str());
242  }
243  CHKERR const_cast<FieldSeries *>(&*sit)->begin();
245 }
246 
247 MoFEMErrorCode SeriesRecorder::record_end(const std::string &serie_name,
248  double time) {
249 
251  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
252  sEries.get<SeriesName_mi_tag>().find(serie_name);
253  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
254  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
255  serie_name.c_str());
256  }
257  CHKERR const_cast<FieldSeries *>(&*sit)->end(time);
259 }
260 
262 SeriesRecorder::initialize_series_recorder(const std::string &serie_name) {
263 
264  Interface &m_field = cOre;
265  moab::Interface &moab = m_field.get_moab();
267  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
268  sEries.get<SeriesName_mi_tag>().find(serie_name);
269  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
270  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
271  serie_name.c_str());
272  }
273  CHKERR const_cast<FieldSeries *>(&*sit)->read(moab);
274  const_cast<FieldSeries *>(&*sit)->record_begin = false;
275  const_cast<FieldSeries *>(&*sit)->record_end = false;
277 }
278 
280 SeriesRecorder::finalize_series_recorder(const std::string &serie_name) {
281 
282  Interface &m_field = cOre;
283  moab::Interface &moab = m_field.get_moab();
285  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
286  sEries.get<SeriesName_mi_tag>().find(serie_name);
287  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
288  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
289  serie_name.c_str());
290  }
291  CHKERR sit->save(moab);
292  int nb_steps;
293  CHKERR sit->get_nb_steps(moab, nb_steps);
294  int ss = 0;
295  for (; ss < nb_steps; ss++) {
296  /*std::pair<SeriesStep_multiIndex::iterator,bool> p =*/seriesSteps.insert(
297  FieldSeriesStep(moab, &*sit, ss));
298  }
300 }
301 
304 
305  for (auto &sit : sEries.get<SeriesName_mi_tag>())
306  MOFEM_LOG("SYNC", LogManager::SeverityLevel::inform) << "series " << sit;
307 
308  for (auto &ssit : seriesSteps.get<SeriesName_mi_tag>())
309  MOFEM_LOG("SYNC", LogManager::SeverityLevel::inform)
310  << "series steps " << ssit;
311 
313 }
314 
315 bool SeriesRecorder::check_series(const std::string &name) const {
316  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
317  sEries.get<SeriesName_mi_tag>().find(name);
318  if (sit != sEries.get<SeriesName_mi_tag>().end())
319  return true;
320  return false;
321 }
322 
323 MoFEMErrorCode SeriesRecorder::load_series_data(const std::string &serie_name,
324  const int step_number) {
325 
326  Interface &m_field = cOre;
327  moab::Interface &moab = m_field.get_moab();
328  const DofEntity_multiIndex *dofs_ptr;
330  SeriesStep_multiIndex::index<
331  Composite_SeriesName_And_Step_mi_tag>::type::iterator sit;
333  boost::make_tuple(serie_name, step_number));
334  if (sit == seriesSteps.get<Composite_SeriesName_And_Step_mi_tag>().end()) {
335  SETERRQ2(PETSC_COMM_SELF, 1, "series <%s> and step %d not found",
336  serie_name.c_str(), step_number);
337  }
338  CHKERR m_field.get_dofs(&dofs_ptr);
339  CHKERR sit->get(moab, *(const_cast<DofEntity_multiIndex *>(dofs_ptr)));
341 }
342 
343 SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator
345  return seriesSteps.get<SeriesName_mi_tag>().lower_bound(name);
346 }
347 
348 SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator
350  return seriesSteps.get<SeriesName_mi_tag>().upper_bound(name);
351 }
352 
353 } // namespace MoFEM
MoFEMErrorCode initialiseDatabaseFromMesh(int verb=0)
virtual MoFEMErrorCode delete_recorder_series(const std::string &series_name)
virtual MoFEMErrorCode print_series_steps()
Structure for recording (time) series.
static const MOFEMuuid IDD_MOFEMSeriesRecorder
virtual MoFEMErrorCode finalize_series_recorder(const std::string &serie_name)
Deprecated interface functions.
MoFEM interface unique ID.
virtual moab::Interface & get_moab()=0
virtual MoFEMErrorCode load_series_data(const std::string &serie_name, const int step_number)
#define SeriesRecorderFunctionBegin
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:290
base class for all interface classes
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
virtual MoFEMErrorCode record_problem(const std::string &serie_name, const Problem *problemPtr, RowColData rc)
Core (interface) class.
Definition: Core.hpp:50
virtual MoFEMErrorCode record_end(const std::string &serie_name, double time=0)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
keeps basic data about problemThis is low level structure with information about problem,...
RowColData
RowColData.
Definition: definitions.h:192
virtual SeriesStep_multiIndex::index< SeriesName_mi_tag >::type::iterator get_series_steps_byName_end(const std::string &name)
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
virtual MoFEMErrorCode get_dofs(const DofEntity_multiIndex **dofs_ptr) const =0
Get dofs multi index.
MoFEMErrorCode getTags(int verb=-1)
get tags handlers used on meshsets
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:84
MoFEMErrorCode clearMap()
clear multi-index container
virtual SeriesStep_multiIndex::index< SeriesName_mi_tag >::type::iterator get_series_steps_byName_begin(const std::string &name)
virtual MoFEMErrorCode initialize_series_recorder(const std::string &serie_name)
virtual MoFEMErrorCode record_field(const std::string &serie_name, const std::string &field_name, const BitRefLevel &bit, const BitRefLevel &mask)
virtual bool check_series(const std::string &name) const
check if series is in database
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< DofEntity, UId, &DofEntity::globalUId > >, 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, const 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 > > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
Tag th_SeriesName
Recorded series name.
SeriesRecorder(const MoFEM::Core &core)
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsRows
store DOFs on rows for this problem
#define CHKERR
Inline error check.
Definition: definitions.h:602
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:50
virtual MoFEMErrorCode record_begin(const std::string &serie_name)
MultiIndex Tag for field name.
int ShortId
Unique Id in the field.
Definition: Types.hpp:42
Series_multiIndex sEries
recorded series
Structure for keeping time and step.
virtual MoFEMErrorCode add_series_recorder(const std::string &series_name)
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1788
SeriesStep_multiIndex seriesSteps
recorded series steps
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
virtual MPI_Comm & get_comm() const =0
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsCols
store DOFs on columns for this problem