OGS 6.1.0-1721-g6382411ad
IntegrationPointWriter.cpp
Go to the documentation of this file.
1 
12 #include <nlohmann/json.hpp>
13 
14 #include "MeshLib/Mesh.h"
15 
16 #include "IntegrationPointWriter.h"
17 
18 using nlohmann::json;
19 
25 {
26  auto const& ip_values = writer.values(/*t, x, dof_table*/);
27  assert(ip_values.size() == mesh.getNumberOfElements());
28 
29  // create field data and fill it with nodal values, and an offsets cell
30  // array indicating where the cell's integration point data starts.
31  auto& field_data = *MeshLib::getOrCreateMeshProperty<double>(
33  writer.numberOfComponents());
34  field_data.clear();
35 
36  for (std::size_t e = 0; e < ip_values.size(); ++e)
37  {
38  auto const& element_ip_values = ip_values[e];
39  std::copy(element_ip_values.begin(), element_ip_values.end(),
40  std::back_inserter(field_data));
41  }
42 
43  return {writer.name(), writer.numberOfComponents(),
44  writer.integrationOrder()};
45 }
46 
50  MeshLib::Mesh& mesh,
51  std::vector<ProcessLib::IntegrationPointMetaData> const& meta_data)
52 {
53  json json_meta_data;
54  json_meta_data["integration_point_arrays"] = json::array();
55 
56  for (auto const& md : meta_data)
57  {
58  json_meta_data["integration_point_arrays"].push_back(
59  {{"name", md.name},
60  {"number_of_components", md.n_components},
61  {"integration_order", md.integration_order}});
62  }
63 
64  // Store the field data.
65  std::string const json_string = json_meta_data.dump();
66  auto& dictionary = *MeshLib::getOrCreateMeshProperty<char>(
67  mesh, "IntegrationPointMetaData",
69  dictionary.clear();
70  std::copy(json_string.begin(), json_string.end(),
71  std::back_inserter(dictionary));
72 }
73 
77  json const& meta_data, std::string const& name)
78 {
79  for (auto const& md : meta_data["integration_point_arrays"])
80  {
81  if (md["name"] == name)
82  {
83  return {name, md["number_of_components"], md["integration_order"]};
84  }
85  }
86  OGS_FATAL("No integration point meta data with name '%s' found.",
87  name.c_str());
88 }
89 
90 namespace ProcessLib
91 {
93  MeshLib::Mesh& mesh,
94  std::vector<std::unique_ptr<IntegrationPointWriter>> const&
95  integration_point_writer)
96 {
97  std::vector<IntegrationPointMetaData> meta_data;
98  for (auto const& ip_writer : integration_point_writer)
99  {
100  meta_data.push_back(addIntegrationPointData(mesh, *ip_writer));
101  }
102  if (!meta_data.empty())
103  {
104  addIntegrationPointMetaData(mesh, meta_data);
105  }
106 }
107 
109  std::string const& name)
110 {
111  if (!mesh.getProperties().existsPropertyVector<char>(
112  "IntegrationPointMetaData"))
113  {
114  OGS_FATAL(
115  "Integration point data '%s' is present in the vtk field "
116  "data but the required 'IntegrationPointMetaData' array "
117  "is not available.",
118  name.c_str());
119  }
120  auto const& mesh_property_ip_meta_data =
121  *mesh.getProperties().template getPropertyVector<char>(
122  "IntegrationPointMetaData");
123 
124  if (mesh_property_ip_meta_data.getMeshItemType() !=
126  {
127  OGS_FATAL("IntegrationPointMetaData array must be field data.");
128  }
129 
130  // Find the current integration point data entry and extract the
131  // meta data.
132  auto const ip_meta_data = extractIntegrationPointMetaData(
133  json::parse(mesh_property_ip_meta_data.begin(),
134  mesh_property_ip_meta_data.end()),
135  name);
136 
137  return ip_meta_data;
138 }
139 } // namespace ProcessLib
static ProcessLib::IntegrationPointMetaData extractIntegrationPointMetaData(json const &meta_data, std::string const &name)
virtual std::string name() const =0
IntegrationPointMetaData getIntegrationPointMetaData(MeshLib::Mesh const &mesh, std::string const &name)
void addIntegrationPointWriter(MeshLib::Mesh &mesh, std::vector< std::unique_ptr< IntegrationPointWriter >> const &integration_point_writer)
virtual int integrationOrder() const =0
virtual int numberOfComponents() const =0
MeshLib::Properties & getProperties()
Definition: Mesh.h:134
Definition of the Mesh class.
static ProcessLib::IntegrationPointMetaData addIntegrationPointData(MeshLib::Mesh &mesh, ProcessLib::IntegrationPointWriter const &writer)
static void addIntegrationPointMetaData(MeshLib::Mesh &mesh, std::vector< ProcessLib::IntegrationPointMetaData > const &meta_data)
bool existsPropertyVector(std::string const &name) const
Definition: Properties.h:79
std::size_t getNumberOfElements() const
Get the number of elements.
Definition: Mesh.h:96
#define OGS_FATAL(fmt,...)
Definition: Error.h:71
void copy(MatrixOrVector const &x, MatrixOrVector &y)
Copies x to y.
Definition: LinAlg.h:36
virtual std::vector< std::vector< double > > values() const =0