add controlWriter

src/controlwriter.hh

@@ -0,0 +1,109 @@
+// SPDX-FileCopyrightText: 2021-2024 The Ikarus Developers mueller@ibb.uni-stuttgart.de
+// SPDX-License-Identifier: LGPL-3.0-or-later
+
+/**
+ * \file controlvtkwriter.hh
+ * \brief Observer implementation for writing vtk files when notified
+ */
+
+#pragma once
+#include <ikarus/utils/observer/observer.hh>
+#include <ikarus/utils/observer/observermessages.hh>
+
+#include <string>
+
+#include "dune/iga/io/igadatacollector.hh"
+#include "dune/vtk/writers/unstructuredgridwriter.hh"
+#include <dune/functions/gridfunctions/discreteglobalbasisfunction.hh>
+#include <dune/grid/io/file/vtk/subsamplingvtkwriter.hh>
+
+#pragma GCC diagnostic push
+#pragma GCC diagnostic ignored "-Wswitch-enum"
+
+namespace Ikarus {
+
+/**
+ * \brief ControlSubsamplingVertexVTKWriterIGA class for writing VTK files with subsampling based on control messages.
+ *
+ * \details It inherits from the IObserver class and is specifically designed for handling SOLUTION_CHANGED messages.
+ *
+ * \tparam B The type of the grid basis.
+ */
+template <typename B>
+class ControlSubsamplingVertexVTKWriterIGA : public IObserver<ControlMessages>
+{
+  using Basis                     = B;
+  static constexpr int components = Basis::LocalView::Tree::degree() == 0 ? 1 : Basis::LocalView::Tree::degree();
+
+public:
+  /**
+   * \brief Constructor for ControlSubsamplingVertexVTKWriterIGA.
+   *
+   * Initializes the VTK writer with the provided basis, solution, and refinement levels.
+   *
+   * \param basis The grid basis.
+   * \param sol The solution vector.
+   * \param refinementLevels The refinement levels for subsampling.
+   */
+  ControlSubsamplingVertexVTKWriterIGA(const Basis& basis, const Eigen::VectorXd& sol, int refinementLevels = 0)
+      : basis_{&basis},
+        dataCollector_(basis.gridView(), 0),
+        vtkWriter_(dataCollector_, Dune::Vtk::FormatTypes::ASCII),
+        solution_{&sol} {}
+
+  /**
+   * \brief Set field information for the VTK file.
+   *
+   * \param name The name of the field.
+   * \param type The type of the field.
+   * \param size The size of the field.
+   * \param prec The precision of the field.
+   * \return The field information.
+   */
+  auto setFieldInfo(std::string&& name, Dune::VTK::FieldInfo::Type type, std::size_t size,
+                    Dune::VTK::Precision prec = Dune::VTK::Precision::float32) {
+    fieldInfo_      = Dune::VTK::FieldInfo(std::move(name), type, size, prec);
+    isFieldInfoSet_ = true;
+  }
+
+  /**
+   * \brief Set the file name prefix for VTK files.
+   *
+   * \param p_name The file name prefix.
+   */
+  auto setFileNamePrefix(std::string&& name) { prefixString_ = std::move(name); }
+
+  /**
+   * \brief Implementation of the update method.
+   *
+   * This method is called upon receiving a SOLUTION_CHANGED control message.
+   * It writes VTK files with subsampling based on the provided field information.
+   *
+   * \param message The received control message.
+   */
+  void updateImpl(ControlMessages message) final {
+    assert(isFieldInfoSet_ && "You need to call setFieldInfo first!");
+    switch (message) {
+      case ControlMessages::SOLUTION_CHANGED: {
+        auto disp = Dune::Functions::makeDiscreteGlobalBasisFunction<Dune::FieldVector<double, components>>(*basis_,
+                                                                                                            *solution_);
+        vtkWriter_.addPointData(disp, fieldInfo_);
+        vtkWriter_.write(prefixString_ + std::to_string(step_++));
+      } break;
+      default:
+        break; // default: do nothing when notified
+    }
+  }
+
+private:
+  const Basis* basis_;
+  Dune::Vtk::DiscontinuousIgaDataCollector<typename Basis::GridView> dataCollector_;
+  Dune::Vtk::UnstructuredGridWriter<typename Basis::GridView, Dune::Vtk::DiscontinuousIgaDataCollector<typename Basis::GridView>> vtkWriter_;
+  const Eigen::VectorXd* solution_;
+  int step_{0};
+  Dune::VTK::FieldInfo fieldInfo_{"Default", Dune::VTK::FieldInfo::Type::scalar, 1};
+  std::string prefixString_{};
+  bool isFieldInfoSet_{false};
+};
+} // namespace Ikarus
+#pragma GCC diagnostic pop

src/iga_cpp.cpp

@@ -5,6 +5,7 @@
   #include "config.h"
 #endif
 
+#include "controlwriter.hh"
 #include "timer.h"
 
 #include <dune/functions/functionspacebases/subspacebasis.hh>
@@ -32,6 +33,7 @@
 #include <ikarus/utils/linearalgebrahelper.hh>
 #include <ikarus/utils/nonlinearoperator.hh>
 #include <ikarus/utils/observer/genericobserver.hh>
+#include <ikarus/utils/observer/nonlinearsolverlogger.hh>
 #include <ikarus/utils/observer/observermessages.hh>
 
 auto run_calculation(int degree, int refinement) {
@@ -112,7 +114,20 @@ auto run_calculation(int degree, int refinement) {
   settings.verbosity = 0;
   trustRegion->setup(settings);
 
-  auto informations = trustRegion->solve();
+  auto nonLinearSolverObserver = std::make_shared<Ikarus::NonLinearSolverLogger>();
+  auto controlWriter =
+      std::make_shared<Ikarus::ControlSubsamplingVertexVTKWriterIGA<std::remove_cvref_t<decltype(basis.flat())>>>(
+          basis.flat(), d, 0);
+  controlWriter->setFileNamePrefix("timeline");
+  controlWriter->setFieldInfo("displacement", Dune::VTK::FieldInfo::Type::vector, 3);
+
+  trustRegion->subscribeAll(nonLinearSolverObserver);
+
+  auto lc = Ikarus::LoadControl(trustRegion, 20, {0, 1});
+
+  // auto informations = trustRegion->solve();
+  lc.subscribeAll(controlWriter);
+  auto informations = lc.run();
 
   auto dispGlobalFunc = Dune::Functions::makeDiscreteGlobalBasisFunction<Dune::FieldVector<double, 3>>(basis.flat(), d);
 
@@ -145,7 +160,7 @@ auto run_calculation(int degree, int refinement) {
     ++i;
   }
 
-  return std::make_tuple(*std::ranges::max_element(nodalDisplacements), informations.iterations,
+  return std::make_tuple(*std::ranges::max_element(nodalDisplacements), informations.totalIterations,
                          sparseAssembler.reducedSize());
 }
 
@@ -156,8 +171,8 @@ int main(int argc, char* argv[]) {
   // If we are in testing mode (e.g. through GH Action, we only run one iteration)
   bool testing = argc > 1 && std::strcmp(argv[1], "testing") == 0;
 
-  auto degreeRange     = Dune::range(2, testing ? 3 : 5);
-  auto refinementRange = Dune::range(3, testing ? 4 : 7);
+  auto degreeRange     = Dune::range(3, testing ? 4 : 4);
+  auto refinementRange = Dune::range(4, testing ? 5 : 5);
 
   std::vector<std::tuple<int, int, double, int, int, Timer<>::Period>> results{};
   for (auto i : degreeRange) {