v0.13.0
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 SeriesRecorder::query_interface(boost::typeindex::type_index type_index,
33  UnknownInterface **iface) const {
34  *iface = const_cast<SeriesRecorder *>(this);
35  return 0;
36 }
37 
39  : cOre(const_cast<Core &>(core)) {}
40 
42 
43  Interface &m_field = cOre;
44  moab::Interface &moab = m_field.get_moab();
46  const int def_val_len = 0;
47  CHKERR moab.tag_get_handle(
48  "_SeriesName", def_val_len, MB_TYPE_OPAQUE, th_SeriesName,
49  MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
51 }
52 
55  sEries.clear();
56  seriesSteps.clear();
58 }
59 
61  Interface &m_field = cOre;
62  moab::Interface &moab = m_field.get_moab();
64  Range meshsets;
65  const int def_val_len = 0;
66  CHKERR moab.tag_get_handle(
67  "_SeriesName", def_val_len, MB_TYPE_OPAQUE, th_SeriesName,
68  MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
69 
70  CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshsets, false);
71  for (Range::iterator mit = meshsets.begin(); mit != meshsets.end(); mit++) {
72  const void *tagName;
73  int tagNameSize;
74  rval = moab.tag_get_by_ptr(th_SeriesName, &*mit, 1, (const void **)&tagName,
75  &tagNameSize);
76  if (rval == MB_SUCCESS) {
77  std::pair<Series_multiIndex::iterator, bool> p =
78  sEries.insert(FieldSeries(moab, *mit));
79  if (verb > QUIET)
80  MOFEM_LOG("SYNC", Sev::inform) << "read series " << *p.first;
81  }
82  }
83  // //build series steps
84  for (Series_multiIndex::iterator sit = sEries.begin(); sit != sEries.end();
85  sit++) {
86  int nb_steps;
87  CHKERR sit->get_nb_steps(moab, nb_steps);
88  int ss = 0;
89  for (; ss < nb_steps; ss++) {
90  std::pair<SeriesStep_multiIndex::iterator, bool> p =
91  seriesSteps.insert(FieldSeriesStep(moab, &*sit, ss));
92  if (verb > QUIET)
93  MOFEM_LOG("SYNC", Sev::inform) << "add series step " << *p.first;
94  }
95  }
96  MOFEM_LOG_CHANNEL("WORLD");
97  MOFEM_LOG_CHANNEL("SYNC");
98  MOFEM_LOG_CHANNEL("SELF");
100 }
101 
103 SeriesRecorder::add_series_recorder(const std::string &series_name) {
104 
105  Interface &m_field = cOre;
106  moab::Interface &moab = m_field.get_moab();
108  EntityHandle meshset;
109  CHKERR moab.create_meshset(MESHSET_SET, meshset);
110  void const *tag_data[] = {series_name.c_str()};
111  int tag_sizes[1];
112  tag_sizes[0] = series_name.size();
113  CHKERR moab.tag_set_by_ptr(th_SeriesName, &meshset, 1, tag_data, tag_sizes);
114  std::pair<Series_multiIndex::iterator, bool> p =
115  sEries.insert(FieldSeries(moab, meshset));
116  if (!p.second) {
117  SETERRQ1(PETSC_COMM_SELF, 1, "series recorder <%s> is already there",
118  series_name.c_str());
119  }
121 }
122 
124 SeriesRecorder::delete_recorder_series(const std::string &series_name) {
125 
126  Interface &m_field = cOre;
127  moab::Interface &moab = m_field.get_moab();
129  SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator ssit, hi_ssit;
130  ssit = seriesSteps.get<SeriesName_mi_tag>().lower_bound(series_name);
131  hi_ssit = seriesSteps.get<SeriesName_mi_tag>().upper_bound(series_name);
132  seriesSteps.get<SeriesName_mi_tag>().erase(ssit, hi_ssit);
133  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit;
134  sit = sEries.get<SeriesName_mi_tag>().find(series_name);
135  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
136  SETERRQ1(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
137  "series recorder <%s> not exist and can be deleted",
138  series_name.c_str());
139  }
140  EntityHandle series_meshset = sit->getMeshset();
141  CHKERR moab.tag_delete(sit->th_SeriesTime);
142  CHKERR moab.tag_delete(sit->th_SeriesDataHandles);
143  CHKERR moab.tag_delete(sit->th_SeriesDataUIDs);
144  CHKERR moab.tag_delete(sit->th_SeriesData);
145  sEries.get<SeriesName_mi_tag>().erase(sit);
146  std::vector<EntityHandle> contained;
147  CHKERR moab.get_contained_meshsets(series_meshset, contained);
148  CHKERR moab.remove_entities(series_meshset, &contained[0], contained.size());
149  CHKERR moab.delete_entities(&contained[0], contained.size());
150  CHKERR moab.delete_entities(&series_meshset, 1);
152 }
153 
154 MoFEMErrorCode SeriesRecorder::record_problem(const std::string &serie_name,
155  const Problem *problemPtr,
156  RowColData rc) {
158 
159  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
160  sEries.get<SeriesName_mi_tag>().find(serie_name);
161  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
162  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
163  "series recorder <%s> not exist", serie_name.c_str());
164  }
165  switch (rc) {
166  case ROW:
167  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(
168  problemPtr->numeredRowDofsPtr->begin(),
169  problemPtr->numeredRowDofsPtr->end());
170  break;
171  case COL:
172  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(
173  problemPtr->numeredColDofsPtr->begin(),
174  problemPtr->numeredColDofsPtr->end());
175  break;
176  default:
177  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
178  }
180 }
181 
182 MoFEMErrorCode SeriesRecorder::record_problem(const std::string &serie_name,
183  const std::string &problem_name,
184  RowColData rc) {
185 
186  Interface &m_field = cOre;
188  const Problem *problem_ptr;
189  CHKERR m_field.get_problem(problem_name, &problem_ptr);
190  CHKERR record_problem(serie_name, problem_ptr, rc);
192 }
193 
194 MoFEMErrorCode SeriesRecorder::record_field(const std::string &serie_name,
195  const std::string &field_name,
196  const BitRefLevel &bit,
197  const BitRefLevel &mask) {
198  Interface &m_field = cOre;
199  auto dofs_ptr = m_field.get_dofs();
201  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
202  sEries.get<SeriesName_mi_tag>().find(serie_name);
203  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
204  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "serie recorder <%s> not exist",
205  serie_name.c_str());
206  }
207 
208  const auto bit_number = m_field.get_field_bit_number(field_name);
209 
210  auto dit = dofs_ptr->get<Unique_mi_tag>().lower_bound(
211  FieldEntity::getLoBitNumberUId(bit_number));
212  if (dit == dofs_ptr->get<Unique_mi_tag>().end())
213  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "field <%s> not exist", field_name.c_str());
214  auto hi_dit = dofs_ptr->get<Unique_mi_tag>().upper_bound(
215  FieldEntity::getHiBitNumberUId(bit_number));
216 
217  for (; dit != hi_dit; dit++) {
218  const BitRefLevel &dof_bit = (*dit)->getBitRefLevel();
219  if ((dof_bit & mask) != dof_bit)
220  continue;
221  if ((dof_bit & bit).any()) {
222  EntityHandle ent = (*dit)->getEnt();
223  ShortId uid = (*dit)->getNonNonuniqueShortId();
224  FieldData val = (*dit)->getFieldData();
225  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(ent, uid, val);
226  }
227  }
229 }
230 
231 MoFEMErrorCode SeriesRecorder::record_begin(const std::string &serie_name) {
232 
234  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
235  sEries.get<SeriesName_mi_tag>().find(serie_name);
236  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
237  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
238  serie_name.c_str());
239  }
240  CHKERR const_cast<FieldSeries *>(&*sit)->begin();
242 }
243 
244 MoFEMErrorCode SeriesRecorder::record_end(const std::string &serie_name,
245  double time) {
246 
248  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
249  sEries.get<SeriesName_mi_tag>().find(serie_name);
250  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
251  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
252  serie_name.c_str());
253  }
254  CHKERR const_cast<FieldSeries *>(&*sit)->end(time);
256 }
257 
259 SeriesRecorder::initialize_series_recorder(const std::string &serie_name) {
260 
261  Interface &m_field = cOre;
262  moab::Interface &moab = m_field.get_moab();
264  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
265  sEries.get<SeriesName_mi_tag>().find(serie_name);
266  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
267  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
268  serie_name.c_str());
269  }
270  CHKERR const_cast<FieldSeries *>(&*sit)->read(moab);
271  const_cast<FieldSeries *>(&*sit)->record_begin = false;
272  const_cast<FieldSeries *>(&*sit)->record_end = false;
274 }
275 
277 SeriesRecorder::finalize_series_recorder(const std::string &serie_name) {
278 
279  Interface &m_field = cOre;
280  moab::Interface &moab = m_field.get_moab();
282  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
283  sEries.get<SeriesName_mi_tag>().find(serie_name);
284  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
285  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
286  serie_name.c_str());
287  }
288  CHKERR sit->save(moab);
289  int nb_steps;
290  CHKERR sit->get_nb_steps(moab, nb_steps);
291  int ss = 0;
292  for (; ss < nb_steps; ss++) {
293  /*std::pair<SeriesStep_multiIndex::iterator,bool> p =*/seriesSteps.insert(
294  FieldSeriesStep(moab, &*sit, ss));
295  }
297 }
298 
301 
302  for (auto &sit : sEries.get<SeriesName_mi_tag>())
303  MOFEM_LOG("SYNC", Sev::inform) << "series " << sit;
304 
305  for (auto &ssit : seriesSteps.get<SeriesName_mi_tag>())
306  MOFEM_LOG("SYNC", Sev::inform) << "series steps " << ssit;
307 
308  MOFEM_LOG_CHANNEL("WORLD");
309  MOFEM_LOG_CHANNEL("SYNC");
310  MOFEM_LOG_CHANNEL("SELF");
312 }
313 
314 bool SeriesRecorder::check_series(const std::string &name) const {
315  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
316  sEries.get<SeriesName_mi_tag>().find(name);
317  if (sit != sEries.get<SeriesName_mi_tag>().end())
318  return true;
319  return false;
320 }
321 
322 MoFEMErrorCode SeriesRecorder::load_series_data(const std::string &serie_name,
323  const int step_number) {
324 
325  Interface &m_field = cOre;
326  moab::Interface &moab = m_field.get_moab();
327  auto dofs_ptr = m_field.get_dofs();
330  Composite_SeriesName_And_Step_mi_tag>::type::iterator sit;
332  boost::make_tuple(serie_name, step_number));
333  if (sit == seriesSteps.get<Composite_SeriesName_And_Step_mi_tag>().end()) {
334  SETERRQ2(PETSC_COMM_SELF, 1, "series <%s> and step %d not found",
335  serie_name.c_str(), step_number);
336  }
337  CHKERR sit->get(moab, *(const_cast<DofEntity_multiIndex *>(dofs_ptr)));
339 }
340 
341 SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator
343  return seriesSteps.get<SeriesName_mi_tag>().lower_bound(name);
344 }
345 
346 SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator
348  return seriesSteps.get<SeriesName_mi_tag>().upper_bound(name);
349 }
350 
351 } // namespace MoFEM
static Index< 'p', 3 > p
#define SeriesRecorderFunctionBegin
@ QUIET
Definition: definitions.h:221
RowColData
RowColData.
Definition: definitions.h:136
@ COL
Definition: definitions.h:136
@ ROW
Definition: definitions.h:136
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_NOT_FOUND
Definition: definitions.h:46
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< DofEntity, UId, &DofEntity::getLocalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:287
virtual MoFEMErrorCode load_series_data(const std::string &serie_name, const int step_number)
virtual MoFEMErrorCode print_series_steps()
virtual MoFEMErrorCode finalize_series_recorder(const std::string &serie_name)
virtual bool check_series(const std::string &name) const
check if series is in database
virtual MoFEMErrorCode record_problem(const std::string &serie_name, const Problem *problemPtr, RowColData rc)
virtual MoFEMErrorCode record_field(const std::string &serie_name, const std::string &field_name, const BitRefLevel &bit, const BitRefLevel &mask)
virtual MoFEMErrorCode initialize_series_recorder(const std::string &serie_name)
virtual MoFEMErrorCode add_series_recorder(const std::string &series_name)
virtual MoFEMErrorCode record_begin(const std::string &serie_name)
virtual MoFEMErrorCode record_end(const std::string &serie_name, double time=0)
virtual MoFEMErrorCode delete_recorder_series(const std::string &series_name)
auto bit
set bit
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:85
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
double FieldData
Field data type.
Definition: Types.hpp:36
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
int ShortId
Unique Id in the field.
Definition: Types.hpp:43
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
Core (interface) class.
Definition: Core.hpp:92
Deprecated interface functions.
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
Structure for recording (time) series.
Structure for keeping time and step.
keeps basic data about problem
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredColDofsPtr
store DOFs on columns for this problem
virtual SeriesStep_multiIndex::index< SeriesName_mi_tag >::type::iterator get_series_steps_byName_begin(const std::string &name)
Series_multiIndex sEries
recorded series
MoFEMErrorCode clearMap()
clear multi-index container
SeriesStep_multiIndex seriesSteps
recorded series steps
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
SeriesRecorder(const MoFEM::Core &core)
Tag th_SeriesName
Recorded series name.
MoFEMErrorCode getTags(int verb=-1)
get tags handlers used on meshsets
virtual SeriesStep_multiIndex::index< SeriesName_mi_tag >::type::iterator get_series_steps_byName_end(const std::string &name)
MoFEMErrorCode initialiseDatabaseFromMesh(int verb=0)
base class for all interface classes