v0.9.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 namespace MoFEM {
20 
21 // const static int debug = 1;
22 
24  UnknownInterface **iface) const {
26  *iface = NULL;
27  if (uuid == IDD_MOFEMSeriesRecorder) {
28  *iface = const_cast<SeriesRecorder *>(this);
30  }
31  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
33 }
34 
36  : cOre(const_cast<Core &>(core)) {}
37 
39 
40  Interface &m_field = cOre;
41  moab::Interface &moab = m_field.get_moab();
43  const int def_val_len = 0;
44  CHKERR moab.tag_get_handle(
45  "_SeriesName", def_val_len, MB_TYPE_OPAQUE, th_SeriesName,
46  MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
48 }
49 
52  sEries.clear();
53  seriesSteps.clear();
55 }
56 
58  Interface &m_field = cOre;
59  moab::Interface &moab = m_field.get_moab();
61  Range meshsets;
62  const int def_val_len = 0;
63  CHKERR moab.tag_get_handle(
64  "_SeriesName", def_val_len, MB_TYPE_OPAQUE, th_SeriesName,
65  MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
66 
67  CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshsets, false);
68  for (Range::iterator mit = meshsets.begin(); mit != meshsets.end(); mit++) {
69  const void *tagName;
70  int tagNameSize;
71  rval = moab.tag_get_by_ptr(th_SeriesName, &*mit, 1, (const void **)&tagName,
72  &tagNameSize);
73  if (rval == MB_SUCCESS) {
74  std::pair<Series_multiIndex::iterator, bool> p =
75  sEries.insert(FieldSeries(moab, *mit));
76  if (verb > 0) {
77  std::ostringstream ss;
78  ss << "read series " << *p.first << std::endl;
79  PetscPrintf(m_field.get_comm(), ss.str().c_str());
80  }
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 > 0) {
93  std::ostringstream ss;
94  ss << "add series step " << *p.first << std::endl;
95  PetscPrintf(m_field.get_comm(), ss.str().c_str());
96  }
97  }
98  }
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_TRACK_OWNER, 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->numeredDofsRows->begin(),
169  problemPtr->numeredDofsRows->end());
170  break;
171  case COL:
172  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(
173  problemPtr->numeredDofsCols->begin(),
174  problemPtr->numeredDofsCols->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 
199  Interface &m_field = cOre;
200  const DofEntity_multiIndex *dofs_ptr;
202  CHKERR m_field.get_dofs(&dofs_ptr);
203  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
204  sEries.get<SeriesName_mi_tag>().find(serie_name);
205  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
206  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
207  serie_name.c_str());
208  }
209  DofEntityByFieldName::iterator dit =
210  dofs_ptr->get<FieldName_mi_tag>().lower_bound(field_name);
211  if (dit == dofs_ptr->get<FieldName_mi_tag>().end()) {
212  SETERRQ1(PETSC_COMM_SELF, 1, "field <%s> not exist", field_name.c_str());
213  }
214  DofEntityByFieldName::iterator hi_dit =
215  dofs_ptr->get<FieldName_mi_tag>().upper_bound(field_name);
216  for (; dit != hi_dit; dit++) {
217  const BitRefLevel &dof_bit = (*dit)->getBitRefLevel();
218  if ((dof_bit & mask) != dof_bit)
219  continue;
220  if ((dof_bit & bit).any()) {
221  EntityHandle ent = (*dit)->getEnt();
222  ShortId uid = (*dit)->getNonNonuniqueShortId();
223  FieldData val = (*dit)->getFieldData();
224  CHKERR const_cast<FieldSeries *>(&*sit)->push_dofs(ent, uid, val);
225  }
226  }
228 }
229 
230 MoFEMErrorCode SeriesRecorder::record_begin(const std::string &serie_name) {
231 
233  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
234  sEries.get<SeriesName_mi_tag>().find(serie_name);
235  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
236  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
237  serie_name.c_str());
238  }
239  CHKERR const_cast<FieldSeries *>(&*sit)->begin();
241 }
242 
243 MoFEMErrorCode SeriesRecorder::record_end(const std::string &serie_name,
244  double time) {
245 
247  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
248  sEries.get<SeriesName_mi_tag>().find(serie_name);
249  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
250  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
251  serie_name.c_str());
252  }
253  CHKERR const_cast<FieldSeries *>(&*sit)->end(time);
255 }
256 
258 SeriesRecorder::initialize_series_recorder(const std::string &serie_name) {
259 
260  Interface &m_field = cOre;
261  moab::Interface &moab = m_field.get_moab();
263  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
264  sEries.get<SeriesName_mi_tag>().find(serie_name);
265  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
266  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
267  serie_name.c_str());
268  }
269  CHKERR const_cast<FieldSeries *>(&*sit)->read(moab);
270  const_cast<FieldSeries *>(&*sit)->record_begin = false;
271  const_cast<FieldSeries *>(&*sit)->record_end = false;
273 }
274 
276 SeriesRecorder::finalize_series_recorder(const std::string &serie_name) {
277 
278  Interface &m_field = cOre;
279  moab::Interface &moab = m_field.get_moab();
281  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
282  sEries.get<SeriesName_mi_tag>().find(serie_name);
283  if (sit == sEries.get<SeriesName_mi_tag>().end()) {
284  SETERRQ1(PETSC_COMM_SELF, 1, "serie recorder <%s> not exist",
285  serie_name.c_str());
286  }
287  CHKERR sit->save(moab);
288  int nb_steps;
289  CHKERR sit->get_nb_steps(moab, nb_steps);
290  int ss = 0;
291  for (; ss < nb_steps; ss++) {
292  /*std::pair<SeriesStep_multiIndex::iterator,bool> p =*/seriesSteps.insert(
293  FieldSeriesStep(moab, &*sit, ss));
294  }
296 }
297 
299  //
300  Interface &m_field = cOre;
302  std::ostringstream ss;
303  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
304  sEries.get<SeriesName_mi_tag>().begin();
305  for (; sit != sEries.get<SeriesName_mi_tag>().end(); sit++) {
306  ss << "series " << *sit << std::endl;
307  }
308  SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator ssit =
309  seriesSteps.get<SeriesName_mi_tag>().begin();
310  for (; ssit != seriesSteps.get<SeriesName_mi_tag>().end(); ssit++) {
311  ss << "serises steps " << *ssit << std::endl;
312  }
313  PetscPrintf(m_field.get_comm(), ss.str().c_str());
315 }
316 
317 bool SeriesRecorder::check_series(const std::string &name) const {
318  Series_multiIndex::index<SeriesName_mi_tag>::type::iterator sit =
319  sEries.get<SeriesName_mi_tag>().find(name);
320  if (sit != sEries.get<SeriesName_mi_tag>().end())
321  return true;
322  return false;
323 }
324 
325 MoFEMErrorCode SeriesRecorder::load_series_data(const std::string &serie_name,
326  const int step_number) {
327 
328  Interface &m_field = cOre;
329  moab::Interface &moab = m_field.get_moab();
330  const DofEntity_multiIndex *dofs_ptr;
332  SeriesStep_multiIndex::index<
333  Composite_SeriesName_And_Step_mi_tag>::type::iterator sit;
335  boost::make_tuple(serie_name, step_number));
336  if (sit == seriesSteps.get<Composite_SeriesName_And_Step_mi_tag>().end()) {
337  SETERRQ2(PETSC_COMM_SELF, 1, "series <%s> and step %d not found",
338  serie_name.c_str(), step_number);
339  }
340  CHKERR m_field.get_dofs(&dofs_ptr);
341  CHKERR sit->get(moab, *(const_cast<DofEntity_multiIndex *>(dofs_ptr)));
343 }
344 
345 SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator
347  return seriesSteps.get<SeriesName_mi_tag>().lower_bound(name);
348 }
349 
350 SeriesStep_multiIndex::index<SeriesName_mi_tag>::type::iterator
352  return seriesSteps.get<SeriesName_mi_tag>().upper_bound(name);
353 }
354 
355 } // 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)
MoFEM interface unique ID.
virtual moab::Interface & get_moab()=0
virtual MoFEMErrorCode load_series_data(const std::string &serie_name, const int step_number)
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:477
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:508
keeps basic data about problemThis is low level structure with information about problem,...
RowColData
RowColData.
Definition: definitions.h:186
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:596
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
virtual MoFEMErrorCode record_begin(const std::string &serie_name)
MultiIndex Tag for field name.
Series_multiIndex sEries
recorded series
Structure for keeping time and step.
virtual MoFEMErrorCode add_series_recorder(const std::string &series_name)
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:407
virtual MPI_Comm & get_comm() const =0
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsCols
store DOFs on columns for this problem