OGS 6.2.1-97-g73d1aeda3
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 (const auto& element_ip_values : ip_values)
37  {
38  std::copy(element_ip_values.begin(), element_ip_values.end(),
39  std::back_inserter(field_data));
40  }
41 
42  return {writer.name(), writer.numberOfComponents(),
43  writer.integrationOrder()};
44 }
45 
49  MeshLib::Mesh& mesh,
50  std::vector<ProcessLib::IntegrationPointMetaData> const& meta_data)
51 {
52  json json_meta_data;
53  json_meta_data["integration_point_arrays"] = json::array();
54 
55  for (auto const& md : meta_data)
56  {
57  json_meta_data["integration_point_arrays"].push_back(
58  {{"name", md.name},
59  {"number_of_components", md.n_components},
60  {"integration_order", md.integration_order}});
61  }
62 
63  // Store the field data.
64  std::string const json_string = json_meta_data.dump();
65  auto& dictionary = *MeshLib::getOrCreateMeshProperty<char>(
66  mesh, "IntegrationPointMetaData",
68  dictionary.clear();
69  std::copy(json_string.begin(), json_string.end(),
70  std::back_inserter(dictionary));
71 }
72 
76  json const& meta_data, std::string const& name)
77 {
78  for (auto const& md : meta_data["integration_point_arrays"])
79  {
80  if (md["name"] == name)
81  {
82  return {name, md["number_of_components"], md["integration_order"]};
83  }
84  }
85  OGS_FATAL("No integration point meta data with name '%s' found.",
86  name.c_str());
87 }
88 
89 namespace ProcessLib
90 {
92  MeshLib::Mesh& mesh,
93  std::vector<std::unique_ptr<IntegrationPointWriter>> const&
94  integration_point_writer)
95 {
96  std::vector<IntegrationPointMetaData> meta_data;
97  for (auto const& ip_writer : integration_point_writer)
98  {
99  meta_data.push_back(addIntegrationPointData(mesh, *ip_writer));
100  }
101  if (!meta_data.empty())
102  {
103  addIntegrationPointMetaData(mesh, meta_data);
104  }
105 }
106 
108  std::string const& name)
109 {
110  if (!mesh.getProperties().existsPropertyVector<char>(
111  "IntegrationPointMetaData"))
112  {
113  OGS_FATAL(
114  "Integration point data '%s' is present in the vtk field "
115  "data but the required 'IntegrationPointMetaData' array "
116  "is not available.",
117  name.c_str());
118  }
119  auto const& mesh_property_ip_meta_data =
120  *mesh.getProperties().template getPropertyVector<char>(
121  "IntegrationPointMetaData");
122 
123  if (mesh_property_ip_meta_data.getMeshItemType() !=
125  {
126  OGS_FATAL("IntegrationPointMetaData array must be field data.");
127  }
128 
129  // Find the current integration point data entry and extract the
130  // meta data.
131  auto const ip_meta_data = extractIntegrationPointMetaData(
132  json::parse(mesh_property_ip_meta_data.begin(),
133  mesh_property_ip_meta_data.end()),
134  name);
135 
136  return ip_meta_data;
137 }
138 } // namespace ProcessLib
static ProcessLib::IntegrationPointMetaData extractIntegrationPointMetaData(json const &meta_data, std::string const &name)
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)
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:63
void copy(MatrixOrVector const &x, MatrixOrVector &y)
Copies x to y.
Definition: LinAlg.h:36
std::vector< std::vector< double > > values() const