OGS
Adaptor.cpp
Go to the documentation of this file.
1
11#include "Adaptor.h"
12
13#include <vtkCPDataDescription.h>
14#include <vtkCPInputDataDescription.h>
15#include <vtkCPProcessor.h>
16#include <vtkCPPythonScriptPipeline.h>
17#include <vtkNew.h>
18#include <vtkUnstructuredGrid.h>
19
20#include <filesystem>
21
22#include "BaseLib/ConfigTree.h"
23#include "MeshLib/Mesh.h"
25
26namespace InSituLib
27{
28vtkCPProcessor* Processor = nullptr;
29
30void Initialize(BaseLib::ConfigTree const& scripts_config,
31 std::string const& path)
32{
33 if (Processor == nullptr)
34 {
35 Processor = vtkCPProcessor::New();
36 Processor->Initialize();
37 }
38 else
39 {
40 Processor->RemoveAllPipelines();
41 }
43 for (auto script_config : scripts_config.getConfigSubtreeList("script"))
44 {
45 auto scriptPath = std::filesystem::path(
47 script_config.getConfigParameter<std::string>("name"));
48 if (scriptPath.is_relative())
49 {
50 scriptPath = std::filesystem::path(path) / scriptPath;
51 }
52 if (!std::filesystem::exists(scriptPath))
53 {
54 ERR("In-situ script {:s} does not exist!", scriptPath.string());
55 }
56 INFO("Initializing in-situ script: {:s}", scriptPath.string());
57 vtkNew<vtkCPPythonScriptPipeline> pipeline;
58 pipeline->Initialize(scriptPath.c_str());
59 Processor->AddPipeline(pipeline.GetPointer());
60 }
61}
63{
64 if (Processor)
65 {
66 Processor->Delete();
67 Processor = nullptr;
68 }
69}
70void CoProcess(MeshLib::Mesh const& mesh, double const time,
71 unsigned int const timeStep, bool const lastTimeStep,
72 std::string output_directory)
73{
74 if (Processor == nullptr)
75 return;
76
77 vtkNew<vtkCPDataDescription> dataDescription;
78 dataDescription->AddInput("input");
79 dataDescription->SetTimeData(time, timeStep);
80 if (lastTimeStep == true)
81 {
82 // assume that we want to all the pipelines to execute if it
83 // is the last time step.
84 dataDescription->ForceOutputOn();
85 }
86 if (Processor->RequestDataDescription(dataDescription.GetPointer()) != 0)
87 {
88 INFO("Start InSitu process: timestep #{:d} (t={:g}, last={:d})",
89 timeStep, time, lastTimeStep);
90 vtkNew<MeshLib::VtkMappedMeshSource> vtkSource;
91 vtkSource->SetMesh(&mesh);
92 vtkSource->Update();
93 dataDescription->GetInputDescriptionByName("input")->SetGrid(
94 vtkSource->GetOutput());
95 auto const cwd = std::filesystem::current_path();
96 std::filesystem::current_path(std::move(output_directory));
97 Processor->CoProcess(dataDescription.GetPointer());
98 std::filesystem::current_path(cwd);
99 INFO("End InSitu process.");
100 }
101}
102} // namespace InSituLib
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:45
Definition of the Mesh class.
VtkMappedMeshSource is a source class to transform OGS meshes into complete vtkUnstructuredGrids....
Range< SubtreeIterator > getConfigSubtreeList(std::string const &root) const
void Initialize(BaseLib::ConfigTree const &scripts_config, std::string const &path)
Definition Adaptor.cpp:30
void Finalize()
Definition Adaptor.cpp:62
void CoProcess(MeshLib::Mesh const &mesh, double const time, unsigned int const timeStep, bool const lastTimeStep, std::string output_directory)
Definition Adaptor.cpp:70
vtkCPProcessor * Processor
Definition Adaptor.cpp:28