diff --git a/ProcessLib/Output/AddProcessDataToMesh.cpp b/ProcessLib/Output/AddProcessDataToMesh.cpp
index bff6347e64237eda711f3835911799bbd4456513..0e75327578596ef7eb8a251caea0e484070864e6 100644
--- a/ProcessLib/Output/AddProcessDataToMesh.cpp
+++ b/ProcessLib/Output/AddProcessDataToMesh.cpp
@@ -13,6 +13,7 @@
 #include "InfoLib/GitInfo.h"
 #include "MeshLib/Utils/IntegrationPointWriter.h"
 #include "MeshLib/Utils/getOrCreateMeshProperty.h"
+#include "NumLib/TimeStepping/Time.h"
 #include "ProcessLib/Output/SecondaryVariable.h"
 #include "ProcessLib/ProcessVariable.h"
 #ifdef USE_PETSC
@@ -33,7 +34,7 @@ static void addOgsVersion(MeshLib::Mesh& mesh)
 }
 
 static void addSecondaryVariableNodes(
-    double const t,
+    NumLib::Time const& t,
     std::vector<GlobalVector*> const& x,
     std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
     ProcessLib::SecondaryVariable const& var,
@@ -58,7 +59,7 @@ static void addSecondaryVariableNodes(
 
     std::unique_ptr<GlobalVector> result_cache;
     auto const& nodal_values =
-        var.fcts.eval_field(t, x, dof_tables, result_cache);
+        var.fcts.eval_field(t(), x, dof_tables, result_cache);
 
 #ifdef USE_PETSC
     std::size_t const global_vector_size =
@@ -79,7 +80,7 @@ static void addSecondaryVariableNodes(
 }
 
 static void addSecondaryVariableResiduals(
-    double const t,
+    NumLib::Time const& t,
     std::vector<GlobalVector*> const& x,
     std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
     ProcessLib::SecondaryVariable const& var,
@@ -110,7 +111,7 @@ static void addSecondaryVariableResiduals(
 
     std::unique_ptr<GlobalVector> result_cache;
     auto const& residuals =
-        var.fcts.eval_residuals(t, x, dof_table, result_cache);
+        var.fcts.eval_residuals(t(), x, dof_table, result_cache);
 #ifdef USE_PETSC
     std::size_t const global_vector_size =
         residuals.getLocalSize() + residuals.getGhostSize();
@@ -344,8 +345,9 @@ static std::set<std::string> addPrimaryVariablesToMesh(
 
 static void addSecondaryVariablesToMesh(
     ProcessLib::SecondaryVariableCollection const& secondary_variables,
-    std::set<std::string>& names_of_already_output_variables, const double t,
-    std::vector<GlobalVector*> const& xs, MeshLib::Mesh& mesh,
+    std::set<std::string>& names_of_already_output_variables,
+    NumLib::Time const& t, std::vector<GlobalVector*> const& xs,
+    MeshLib::Mesh& mesh,
     std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
     bool const output_residuals)
 {
@@ -371,7 +373,7 @@ static void addSecondaryVariablesToMesh(
 
 namespace ProcessLib
 {
-void addProcessDataToMesh(const double t,
+void addProcessDataToMesh(NumLib::Time const& t,
                           std::vector<GlobalVector*> const& xs,
                           int const process_id,
                           ProcessOutputData const& process_output_data,
diff --git a/ProcessLib/Output/AddProcessDataToMesh.h b/ProcessLib/Output/AddProcessDataToMesh.h
index bc56d8fe27811410d33770729ecf056ee730b2e4..2f42195668997e7b71041445bf6c401c14a1de74 100644
--- a/ProcessLib/Output/AddProcessDataToMesh.h
+++ b/ProcessLib/Output/AddProcessDataToMesh.h
@@ -18,7 +18,8 @@ namespace ProcessLib
 /// Adds output data to a mesh.
 ///
 /// \note The mesh is passed via \c process_output_data.
-void addProcessDataToMesh(const double t, std::vector<GlobalVector*> const& xs,
+void addProcessDataToMesh(NumLib::Time const& t,
+                          std::vector<GlobalVector*> const& xs,
                           int const process_id,
                           ProcessOutputData const& process_output_data,
                           bool const output_secondary_variables,
diff --git a/ProcessLib/Output/Output.cpp b/ProcessLib/Output/Output.cpp
index 799af5ef7c80fa3efccf7f7048e5a3baa795b02d..fdb7770ab221b39b3b07df3cfd8c4bd52a64a9d8 100644
--- a/ProcessLib/Output/Output.cpp
+++ b/ProcessLib/Output/Output.cpp
@@ -24,6 +24,7 @@
 #include "BaseLib/RunTime.h"
 #include "MeshLib/Mesh.h"
 #include "MeshLib/Utils/getOrCreateMeshProperty.h"
+#include "NumLib/TimeStepping/Time.h"
 #include "ProcessLib/Process.h"
 
 namespace ProcessLib
@@ -213,7 +214,7 @@ void Output::outputMeshes(
 
 MeshLib::Mesh const& Output::prepareSubmesh(
     std::string const& submesh_output_name, Process const& process,
-    const int process_id, double const t,
+    const int process_id, NumLib::Time const& t,
     std::vector<GlobalVector*> const& xs) const
 {
     auto& submesh = MeshLib::findMeshByName(_meshes.get(), submesh_output_name);
@@ -288,7 +289,7 @@ MeshLib::Mesh const& Output::prepareSubmesh(
 void Output::doOutputAlways(Process const& process,
                             const int process_id,
                             int const timestep,
-                            const double t,
+                            const NumLib::Time& t,
                             int const iteration,
                             std::vector<GlobalVector*> const& xs) const
 {
@@ -326,7 +327,7 @@ void Output::doOutputAlways(Process const& process,
         }
     }
 
-    outputMeshes(timestep, t, iteration, std::move(output_meshes));
+    outputMeshes(timestep, t(), iteration, std::move(output_meshes));
 
     INFO("[time] Output of timestep {:d} took {:g} s.", timestep,
          time_output.elapsed());
@@ -335,7 +336,7 @@ void Output::doOutputAlways(Process const& process,
 void Output::doOutput(Process const& process,
                       const int process_id,
                       int const timestep,
-                      const double t,
+                      const NumLib::Time& t,
                       int const iteration,
                       std::vector<GlobalVector*> const& xs) const
 {
@@ -354,7 +355,7 @@ void Output::doOutput(Process const& process,
 void Output::doOutputLastTimestep(Process const& process,
                                   const int process_id,
                                   int const timestep,
-                                  const double t,
+                                  const NumLib::Time& t,
                                   int const iteration,
                                   std::vector<GlobalVector*> const& xs) const
 {
@@ -370,7 +371,7 @@ void Output::doOutputLastTimestep(Process const& process,
 
 void Output::doOutputNonlinearIteration(
     Process const& process, const int process_id, int const timestep,
-    const double t, int const iteration,
+    const NumLib::Time& t, int const iteration,
     std::vector<GlobalVector*> const& xs) const
 {
     if (!_output_nonlinear_iteration_results)
@@ -394,7 +395,7 @@ void Output::doOutputNonlinearIteration(
     }
 
     std::string const output_file_name = _output_format->constructFilename(
-        process.getMesh().getName(), timestep, t, iteration);
+        process.getMesh().getName(), timestep, t(), iteration);
 
     std::string const output_file_path =
         BaseLib::joinPaths(_output_format->directory, output_file_name);
@@ -414,7 +415,7 @@ void Output::doOutputNonlinearIteration(
     INFO("[time] Output took {:g} s.", time_output.elapsed());
 }
 
-bool Output::isOutputStep(int const timestep, double const t) const
+bool Output::isOutputStep(int const timestep, NumLib::Time const& t) const
 {
     return _output_data_specification.isOutputStep(timestep, t);
 }
diff --git a/ProcessLib/Output/Output.h b/ProcessLib/Output/Output.h
index 72a0d21a5ebb5bab819101bacdeb7e166a9cf99a..f2e4b18ee27beb959978e5a9b48c45f45a3c3b7b 100644
--- a/ProcessLib/Output/Output.h
+++ b/ProcessLib/Output/Output.h
@@ -19,6 +19,11 @@
 #include "OutputDataSpecification.h"
 #include "OutputFormat.h"
 
+namespace NumLib
+{
+struct Time;
+}
+
 namespace ProcessLib
 {
 class Process;
@@ -55,14 +60,15 @@ public:
     //! Writes output for the given \c process if it should be written in
     //! the given \c timestep.
     void doOutput(Process const& process, const int process_id,
-                  int const timestep, const double t, int const iteration,
+                  int const timestep, const NumLib::Time& t,
+                  int const iteration,
                   std::vector<GlobalVector*> const& xs) const;
 
     //! Writes output for the given \c process if it has not been written yet.
     //! This method is intended for doing output after the last timestep in
     //! order to make sure that its results are written.
     void doOutputLastTimestep(Process const& process, const int process_id,
-                              int const timestep, const double t,
+                              int const timestep, const NumLib::Time& t,
                               int const iteration,
                               std::vector<GlobalVector*> const& xs) const;
 
@@ -70,18 +76,19 @@ public:
     //! This method will always write.
     //! It is intended to write output in error handling routines.
     void doOutputAlways(Process const& process, const int process_id,
-                        int const timestep, const double t, int const iteration,
+                        int const timestep, const NumLib::Time& t,
+                        int const iteration,
                         std::vector<GlobalVector*> const& xs) const;
 
     //! Writes output for the given \c process.
     //! To be used for debug output after an iteration of the nonlinear solver.
     void doOutputNonlinearIteration(Process const& process,
                                     const int process_id, int const timestep,
-                                    const double t, const int iteration,
+                                    const NumLib::Time& t, const int iteration,
                                     std::vector<GlobalVector*> const& xs) const;
 
     //! Tells if output will be written at the specified timestep/time.
-    bool isOutputStep(int const timestep, double const t) const;
+    bool isOutputStep(int const timestep, NumLib::Time const& t) const;
 
     std::vector<double> const& getFixedOutputTimes() const
     {
@@ -105,7 +112,7 @@ private:
 
     MeshLib::Mesh const& prepareSubmesh(
         std::string const& submesh_output_name, Process const& process,
-        const int process_id, double const t,
+        const int process_id, NumLib::Time const& t,
         std::vector<GlobalVector*> const& xs) const;
 
     std::unique_ptr<OutputFormat> _output_format;
diff --git a/ProcessLib/Output/OutputDataSpecification.cpp b/ProcessLib/Output/OutputDataSpecification.cpp
index 40dea119f21d7e8ab3a2012aad8598e9ca30a2f1..f1d29d4a54ba8f3b35f4fbeeb37c9d5ada17bf41 100644
--- a/ProcessLib/Output/OutputDataSpecification.cpp
+++ b/ProcessLib/Output/OutputDataSpecification.cpp
@@ -10,6 +10,7 @@
 
 #include "OutputDataSpecification.h"
 
+#include "NumLib/TimeStepping/Time.h"
 #include "NumLib/TimeStepping/TimeStep.h"
 
 namespace ProcessLib
@@ -47,16 +48,13 @@ OutputDataSpecification::OutputDataSpecification(
 }
 
 bool OutputDataSpecification::isOutputStep(int timestep,
-                                           double const time) const
+                                           NumLib::Time const& time) const
 {
-    auto isFixedOutputStep = [this](double const time) -> bool
+    auto isFixedOutputStep = [this](NumLib::Time const time) -> bool
     {
         return std::any_of(cbegin(fixed_output_times), cend(fixed_output_times),
                            [&](auto fixed_output_time)
-                           {
-                               return (std::abs(fixed_output_time - time) <
-                                       NumLib::TimeStep::minimalTimeStepSize);
-                           });
+                           { return NumLib::Time(fixed_output_time) == time; });
     };
 
     auto isPairRepeatsEachTimeStepOutput = [this](int timestep) -> bool
diff --git a/ProcessLib/Output/OutputDataSpecification.h b/ProcessLib/Output/OutputDataSpecification.h
index 8094ff22cadc3cfca79da035a8617a01f2b9025f..549e5562e6716eebb71d7b4af6dedba820eadecd 100644
--- a/ProcessLib/Output/OutputDataSpecification.h
+++ b/ProcessLib/Output/OutputDataSpecification.h
@@ -17,6 +17,11 @@
 
 #include "BaseLib/Error.h"
 
+namespace NumLib
+{
+struct Time;
+}
+
 namespace ProcessLib
 {
 struct PairRepeatEachSteps
@@ -52,7 +57,7 @@ struct OutputDataSpecification final
 
     //! Determines if there should be output at the given \c timestep or \c
     //! time.
-    bool isOutputStep(int timestep, double const time) const;
+    bool isOutputStep(int timestep, NumLib::Time const& time) const;
 };
 
 std::ostream& operator<<(std::ostream& os, OutputDataSpecification const& o);
diff --git a/ProcessLib/TimeLoop.cpp b/ProcessLib/TimeLoop.cpp
index 0bfe6f12bf9fcabc66d6c8e81955787931fed644..d66fe2343b2fd8393c340fc76578ab410d14671f 100644
--- a/ProcessLib/TimeLoop.cpp
+++ b/ProcessLib/TimeLoop.cpp
@@ -38,9 +38,10 @@ void updateDeactivatedSubdomains(
 }
 
 bool isOutputStep(std::vector<ProcessLib::Output> const& outputs,
-                  const int timestep, const double t, const double end_time)
+                  const int timestep, const NumLib::Time& t,
+                  const NumLib::Time& end_time)
 {
-    if (std::abs(end_time - t) < std::numeric_limits<double>::epsilon())
+    if (end_time == t)
     {
         // the last timestep is an output step
         return true;