OGS 6.1.0-1721-g6382411ad
ProcessOutput.cpp
Go to the documentation of this file.
1 
10 #include "ProcessOutput.h"
11 
12 #include "BaseLib/BuildInfo.h"
13 #include "MathLib/LinAlg/LinAlg.h"
16 
17 #include "IntegrationPointWriter.h"
18 
19 
22 static void addOgsVersion(MeshLib::Mesh& mesh)
23 {
24  auto& ogs_version_field = *MeshLib::getOrCreateMeshProperty<char>(
25  mesh, "OGS_VERSION", MeshLib::MeshItemType::IntegrationPoint, 1);
26 
27  ogs_version_field.assign(BaseLib::BuildInfo::ogs_version.begin(),
29 }
30 
32  double const t,
33  GlobalVector const& x,
34  NumLib::LocalToGlobalIndexMap const& dof_table,
36  std::string const& output_name,
37  MeshLib::Mesh& mesh)
38 {
39  DBUG(" secondary variable %s", output_name.c_str());
40 
41  auto& nodal_values_mesh = *MeshLib::getOrCreateMeshProperty<double>(
42  mesh, output_name, MeshLib::MeshItemType::Node,
43  var.fcts.num_components);
44  if (nodal_values_mesh.size() !=
45  mesh.getNumberOfNodes() * var.fcts.num_components)
46  {
47  OGS_FATAL(
48  "Nodal property `%s' does not have the right number of "
49  "components. Expected: %d, actual: %d",
50  output_name.c_str(),
51  mesh.getNumberOfNodes() * var.fcts.num_components,
52  nodal_values_mesh.size());
53  }
54 
55  std::unique_ptr<GlobalVector> result_cache;
56  auto const& nodal_values =
57  var.fcts.eval_field(t, x, dof_table, result_cache);
58 #ifdef USE_PETSC
59  std::size_t const global_vector_size =
60  nodal_values.getLocalSize() + nodal_values.getGhostSize();
61 #else
62  std::size_t const global_vector_size = nodal_values.size();
63 #endif
64  if (nodal_values_mesh.size() != global_vector_size)
65  {
66  OGS_FATAL(
67  "Secondary variable `%s' did not evaluate to the right "
68  "number of components. Expected: %d, actual: %d.",
69  var.name.c_str(), nodal_values_mesh.size(), global_vector_size);
70  }
71 
72  // Copy result
73  nodal_values.copyValues(nodal_values_mesh);
74 }
75 
77  double const t,
78  GlobalVector const& x,
79  NumLib::LocalToGlobalIndexMap const& dof_table,
81  std::string const& output_name,
82  MeshLib::Mesh& mesh)
83 {
84  if (!var.fcts.eval_residuals)
85  return;
86 
87  DBUG(" secondary variable %s residual", output_name.c_str());
88  auto const& property_name_res = output_name + "_residual";
89 
90  auto& residuals_mesh = *MeshLib::getOrCreateMeshProperty<double>(
91  mesh, property_name_res, MeshLib::MeshItemType::Cell,
92  var.fcts.num_components);
93  if (residuals_mesh.size() !=
94  mesh.getNumberOfElements() * var.fcts.num_components)
95  {
96  OGS_FATAL(
97  "Cell property `%s' does not have the right number of components. "
98  "Expected: %d, actual: %d",
99  property_name_res.c_str(),
100  mesh.getNumberOfElements() * var.fcts.num_components,
101  residuals_mesh.size());
102  }
103 
104  std::unique_ptr<GlobalVector> result_cache;
105  auto const& residuals =
106  var.fcts.eval_residuals(t, x, dof_table, result_cache);
107 #ifdef USE_PETSC
108  std::size_t const global_vector_size =
109  residuals.getLocalSize() + residuals.getGhostSize();
110 #else
111  std::size_t const global_vector_size = residuals.size();
112 #endif
113  if (residuals_mesh.size() != global_vector_size)
114  {
115  OGS_FATAL(
116  "The residual of secondary variable `%s' did not evaluate to the "
117  "right number of components. Expected: %d, actual: %d.",
118  var.name.c_str(), residuals_mesh.size(), global_vector_size);
119  }
120 
121  // Copy result
122  residuals.copyValues(residuals_mesh);
123 }
124 
125 namespace ProcessLib
126 {
128  const double t,
129  GlobalVector const& x,
130  MeshLib::Mesh& mesh,
131  NumLib::LocalToGlobalIndexMap const& dof_table,
132  std::vector<std::reference_wrapper<ProcessVariable>> const&
133  process_variables,
134  SecondaryVariableCollection const& secondary_variables,
135  bool const output_secondary_variable,
136  std::vector<std::unique_ptr<IntegrationPointWriter>> const&
137  integration_point_writer,
138  ProcessOutput const& process_output)
139 {
140  DBUG("Process output data.");
141 
142  addOgsVersion(mesh);
143 
144  // Copy result
145 #ifdef USE_PETSC
146  // TODO It is also possible directly to copy the data for single process
147  // variable to a mesh property. It needs a vector of global indices and
148  // some PETSc magic to do so.
149  std::vector<double> x_copy(x.getLocalSize() + x.getGhostSize());
150 #else
151  std::vector<double> x_copy(x.size());
152 #endif
153  x.copyValues(x_copy);
154 
155  auto const& output_variables = process_output.output_variables;
156  std::set<std::string> already_output;
157 
158  int global_component_offset = 0;
159  int global_component_offset_next = 0;
160 
161  const auto number_of_dof_variables = dof_table.getNumberOfVariables();
162  // primary variables
163  for (int variable_id = 0;
164  variable_id < static_cast<int>(process_variables.size());
165  ++variable_id)
166  {
167  ProcessVariable& pv = process_variables[variable_id];
168  int const n_components = pv.getNumberOfComponents();
169  // If (number_of_dof_variables==1), the case is either the staggered
170  // scheme being applied or a single PDE being solved.
171  const int sub_meshset_id =
172  (number_of_dof_variables == 1) ? 0 : variable_id;
173 
174  if (number_of_dof_variables > 1)
175  {
176  global_component_offset = global_component_offset_next;
177  global_component_offset_next += n_components;
178  }
179 
180  if (output_variables.find(pv.getName()) == output_variables.cend())
181  {
182  continue;
183  }
184 
185  already_output.insert(pv.getName());
186 
187  DBUG(" process variable %s", pv.getName().c_str());
188 
189  auto const num_comp = pv.getNumberOfComponents();
190  auto& output_data = *MeshLib::getOrCreateMeshProperty<double>(
191  mesh, pv.getName(), MeshLib::MeshItemType::Node, num_comp);
192 
193 
194  for (int component_id = 0; component_id < num_comp; ++component_id)
195  {
196  auto const& mesh_subset =
197  dof_table.getMeshSubset(sub_meshset_id, component_id);
198  auto const mesh_id = mesh_subset.getMeshID();
199  for (auto const* node : mesh_subset.getNodes())
200  {
202  node->getID());
203 
204  auto const global_component_id =
205  global_component_offset + component_id;
206  auto const index = dof_table.getLocalIndex(
207  l, global_component_id, x.getRangeBegin(), x.getRangeEnd());
208 
209  output_data[node->getID() * n_components + component_id] =
210  x_copy[index];
211  }
212  }
213  }
214 
215  if (output_secondary_variable)
216  {
217  for (auto const& external_variable_name : secondary_variables)
218  {
219  auto const& name = external_variable_name.first;
220  if (!already_output.insert(name).second)
221  {
222  // no insertion took place, output already done
223  continue;
224  }
225 
227  t, x, dof_table, secondary_variables.get(name), name, mesh);
228 
229  if (process_output.output_residuals)
230  {
232  t, x, dof_table, secondary_variables.get(name), name, mesh);
233  }
234  }
235  }
236 
237  addIntegrationPointWriter(mesh, integration_point_writer);
238 }
239 
240 void makeOutput(std::string const& file_name, MeshLib::Mesh& mesh,
241  bool const compress_output, int const data_mode)
242 {
243  // Write output file
244  DBUG("Writing output to '%s'.", file_name.c_str());
245  MeshLib::IO::VtuInterface vtu_interface(&mesh, data_mode, compress_output);
246  vtu_interface.writeToFile(file_name);
247 }
248 
249 } // namespace ProcessLib
const unsigned num_components
Number of components of the variable.
void addIntegrationPointWriter(MeshLib::Mesh &mesh, std::vector< std::unique_ptr< IntegrationPointWriter >> const &integration_point_writer)
SecondaryVariableFunctions fcts
Functions used for computing the secondary variable.
bool const output_residuals
Tells if also to output extrapolation residuals.
Definition: ProcessOutput.h:26
std::set< std::string > output_variables
All variables that shall be output.
Definition: ProcessOutput.h:23
std::size_t getMeshID() const
return this mesh ID
Definition: MeshSubset.h:71
Implementation of the VtuInterface class.
Function const eval_field
Computes the value of the field at every node of the underlying mesh.
Stores information about a specific secondary variable.
void makeOutput(std::string const &file_name, MeshLib::Mesh &mesh, bool const compress_output, int const data_mode)
Reads and writes VtkXMLUnstructuredGrid-files (vtu) to and from OGS data structures. This class is currently not inherited from Writer because VTK will implement writing to a string from 6.2 onwards.
Definition: VtuInterface.h:36
std::string const name
Name of the variable; used, e.g., for output.
void processOutputData(const double t, GlobalVector const &x, MeshLib::Mesh &mesh, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< std::reference_wrapper< ProcessVariable >> const &process_variables, SecondaryVariableCollection const &secondary_variables, bool const output_secondary_variable, std::vector< std::unique_ptr< IntegrationPointWriter >> const &integration_point_writer, ProcessOutput const &process_output)
Prepare the output data, i.e. add the solution to vtu data structure.
static void addSecondaryVariableResiduals(double const t, GlobalVector const &x, NumLib::LocalToGlobalIndexMap const &dof_table, ProcessLib::SecondaryVariable const &var, std::string const &output_name, MeshLib::Mesh &mesh)
std::string const & getName() const
int getNumberOfComponents() const
Returns the number of components of the process variable.
BASELIB_EXPORT const std::string ogs_version
static void addSecondaryVariableNodes(double const t, GlobalVector const &x, NumLib::LocalToGlobalIndexMap const &dof_table, ProcessLib::SecondaryVariable const &var, std::string const &output_name, MeshLib::Mesh &mesh)
MeshLib::MeshSubset const & getMeshSubset(int const variable_id, int const component_id) const
Handles configuration of several secondary variables from the project file.
#define OGS_FATAL(fmt,...)
Definition: Error.h:71
static void addOgsVersion(MeshLib::Mesh &mesh)
GlobalIndexType getLocalIndex(MeshLib::Location const &l, std::size_t const comp_id, std::size_t const range_begin, std::size_t const range_end) const
Holds information about which variables to write to output files.
Definition: ProcessOutput.h:20