]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add tests for recoverable and irrecoverable situations. 15263/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 24 May 2023 17:04:23 +0000 (11:04 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 24 May 2023 17:05:48 +0000 (11:05 -0600)
tests/sundials/ida_04.cc [new file with mode: 0644]
tests/sundials/ida_04.output [new file with mode: 0644]
tests/sundials/ida_04_in.prm [new file with mode: 0644]
tests/sundials/ida_05.cc [new file with mode: 0644]
tests/sundials/ida_05.output [new file with mode: 0644]

diff --git a/tests/sundials/ida_04.cc b/tests/sundials/ida_04.cc
new file mode 100644 (file)
index 0000000..8f91117
--- /dev/null
@@ -0,0 +1,164 @@
+//-----------------------------------------------------------
+//
+//    Copyright (C) 2017 - 2022 by the deal.II authors
+//
+//    This file is part of the deal.II library.
+//
+//    The deal.II library is free software; you can use it, redistribute
+//    it, and/or modify it under the terms of the GNU Lesser General
+//    Public License as published by the Free Software Foundation; either
+//    version 2.1 of the License, or (at your option) any later version.
+//    The full text of the license can be found in the file LICENSE.md at
+//    the top level directory of deal.II.
+//
+//-----------------------------------------------------------
+
+#include <deal.II/base/parameter_handler.h>
+
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+#include <deal.II/sundials/ida.h>
+
+#include "../tests.h"
+
+
+/**
+ * Like the _03 test, but throw an irrecoverable exception if the time
+ * step is longer than 0.1.
+ *
+ * Solve the exponential decay problem:
+ *
+ * y' = -k y
+ * y (0) = 1
+ *
+ * with k=2. That is
+ *
+ * F(y', y, t) = y' + k y = 0
+ * y_0  = 1
+ *
+ * The exact solution is
+ *
+ * y(t) = exp(-k*t)
+ *
+ * The Jacobian we then need to assemble is the following:
+ *
+ * J = dF/dy + alpha dF/dy'
+ *   = k + alpha 1
+ *   = (k+alpha)
+ */
+class HarmonicOscillator
+{
+public:
+  HarmonicOscillator(
+    double                                                        _kappa,
+    const typename SUNDIALS::IDA<Vector<double>>::AdditionalData &data)
+    : time_stepper(data)
+    , y(1)
+    , y_dot(1)
+    , kappa(_kappa)
+    , last_residual_eval_time(0)
+  {
+    using VectorType = Vector<double>;
+
+    time_stepper.reinit_vector = [&](VectorType &v) { v.reinit(1); };
+
+
+    time_stepper.residual = [&](const double      t,
+                                const VectorType &y,
+                                const VectorType &y_dot,
+                                VectorType &      res) {
+      deallog << "Evaluating residual at t=" << t << std::endl;
+
+      if (t > last_residual_eval_time + 0.1)
+        {
+          deallog << "Time step is too large: " << t - last_residual_eval_time
+                  << ". Failing irrecoverably." << std::endl;
+
+          throw ExcMessage("Irrecoverable exception.");
+        }
+
+      res[0]                  = y_dot[0] + kappa * y[0];
+      last_residual_eval_time = t;
+    };
+
+    time_stepper.setup_jacobian = [&](const double,
+                                      const VectorType &,
+                                      const VectorType &,
+                                      const double alpha) {
+      J = kappa + alpha;
+    };
+
+    // Used only in ver < 4.0.0
+    time_stepper.solve_jacobian_system =
+      [&](const VectorType &src, VectorType &dst) { dst[0] = src[0] / J; };
+
+    // Used in ver >= 4.0.0
+    time_stepper.solve_with_jacobian =
+      [&](const VectorType &src, VectorType &dst, const double) {
+        dst[0] = src[0] / J;
+      };
+
+    time_stepper.output_step = [&](const double       t,
+                                   const VectorType & sol,
+                                   const VectorType & sol_dot,
+                                   const unsigned int step_number) {
+      deallog << "Intermediate output:" << std::endl;
+      deallog << "  t =" << t << std::endl;
+      deallog << "  y =" << sol[0] << "  (exact: " << std::exp(-kappa * t)
+              << ')' << std::endl;
+      deallog << "  y'=" << sol_dot[0]
+              << "  (exact: " << -kappa * std::exp(-kappa * t) << ')'
+              << std::endl;
+    };
+  }
+
+  void
+  run()
+  {
+    y[0]     = 1;
+    y_dot[0] = -kappa;
+    try
+      {
+        time_stepper.solve_dae(y, y_dot);
+      }
+    catch (const std::exception &exc)
+      {
+        deallog << "Caught an irrecoverable exception in a callback:"
+                << std::endl
+                << exc.what() << std::endl;
+      }
+  }
+
+  SUNDIALS::IDA<Vector<double>> time_stepper;
+
+private:
+  Vector<double> y;
+  Vector<double> y_dot;
+  double         J;
+  double         kappa;
+  double         last_residual_eval_time;
+};
+
+
+int
+main()
+{
+  initlog();
+  deallog << std::setprecision(10);
+
+  SUNDIALS::IDA<Vector<double>>::AdditionalData data;
+  ParameterHandler                              prm;
+  data.add_parameters(prm);
+
+  // std::ofstream ofile(SOURCE_DIR "/ida_04.prm");
+  // prm.print_parameters(ofile, ParameterHandler::ShortText);
+  // ofile.close();
+
+  std::ifstream ifile(SOURCE_DIR "/ida_04_in.prm");
+  prm.parse_input(ifile);
+
+
+  HarmonicOscillator ode(2.0, data);
+  ode.run();
+}
diff --git a/tests/sundials/ida_04.output b/tests/sundials/ida_04.output
new file mode 100644 (file)
index 0000000..41b39a7
--- /dev/null
@@ -0,0 +1,18 @@
+
+DEAL::Intermediate output:
+DEAL::  t =0.000000000
+DEAL::  y =1.000000000  (exact: 1.000000000)
+DEAL::  y'=-2.000000000  (exact: -2.000000000)
+DEAL::Evaluating residual at t=0.2000000000
+DEAL::Time step is too large: 0.2000000000. Failing irrecoverably.
+DEAL::Caught an irrecoverable exception in a callback:
+DEAL::
+--------------------------------------------------------
+An error occurred in line <0> of file <> in function
+    
+The violated condition was: 
+    
+Additional information: 
+    Irrecoverable exception.
+--------------------------------------------------------
+
diff --git a/tests/sundials/ida_04_in.prm b/tests/sundials/ida_04_in.prm
new file mode 100644 (file)
index 0000000..957081d
--- /dev/null
@@ -0,0 +1,19 @@
+set Final time                        = 2
+set Initial time                      = 0
+set Time interval between each output = 0.2
+subsection Error control
+  set Absolute error tolerance                      = 1e-10
+  set Ignore algebraic terms for error computations = true
+  set Relative error tolerance                      = 1e-10
+end
+subsection Initial condition correction parameters
+  set Correction type at initial time        = none
+  set Correction type after restart          = none
+  set Maximum number of nonlinear iterations = 10
+end
+subsection Running parameters
+  set Initial step size                      = 0.2 # start with a large time step
+  set Maximum number of nonlinear iterations = 10
+  set Maximum order of BDF                   = 5
+  set Minimum step size                      = 1e-7
+end
diff --git a/tests/sundials/ida_05.cc b/tests/sundials/ida_05.cc
new file mode 100644 (file)
index 0000000..06ea411
--- /dev/null
@@ -0,0 +1,165 @@
+//-----------------------------------------------------------
+//
+//    Copyright (C) 2017 - 2022 by the deal.II authors
+//
+//    This file is part of the deal.II library.
+//
+//    The deal.II library is free software; you can use it, redistribute
+//    it, and/or modify it under the terms of the GNU Lesser General
+//    Public License as published by the Free Software Foundation; either
+//    version 2.1 of the License, or (at your option) any later version.
+//    The full text of the license can be found in the file LICENSE.md at
+//    the top level directory of deal.II.
+//
+//-----------------------------------------------------------
+
+#include <deal.II/base/parameter_handler.h>
+
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+#include <deal.II/sundials/ida.h>
+
+#include "../tests.h"
+
+
+/**
+ * Like the _03 test, but throw a recoverable exception if the time
+ * step is longer than 0.1. We will hopefully see that the solver
+ * reduces the time step and tries again.
+ *
+ * Solve the exponential decay problem:
+ *
+ * y' = -k y
+ * y (0) = 1
+ *
+ * with k=2. That is
+ *
+ * F(y', y, t) = y' + k y = 0
+ * y_0  = 1
+ *
+ * The exact solution is
+ *
+ * y(t) = exp(-k*t)
+ *
+ * The Jacobian we then need to assemble is the following:
+ *
+ * J = dF/dy + alpha dF/dy'
+ *   = k + alpha 1
+ *   = (k+alpha)
+ */
+class HarmonicOscillator
+{
+public:
+  HarmonicOscillator(
+    double                                                        _kappa,
+    const typename SUNDIALS::IDA<Vector<double>>::AdditionalData &data)
+    : time_stepper(data)
+    , y(1)
+    , y_dot(1)
+    , kappa(_kappa)
+    , last_residual_eval_time(0)
+  {
+    using VectorType = Vector<double>;
+
+    time_stepper.reinit_vector = [&](VectorType &v) { v.reinit(1); };
+
+
+    time_stepper.residual = [&](const double      t,
+                                const VectorType &y,
+                                const VectorType &y_dot,
+                                VectorType &      res) {
+      deallog << "Evaluating residual at t=" << t << std::endl;
+
+      if (t > last_residual_eval_time + 0.1)
+        {
+          deallog << "Time step is too large: " << t - last_residual_eval_time
+                  << ". Failing recoverably." << std::endl;
+
+          throw RecoverableUserCallbackError();
+        }
+
+      res[0]                  = y_dot[0] + kappa * y[0];
+      last_residual_eval_time = t;
+    };
+
+    time_stepper.setup_jacobian = [&](const double,
+                                      const VectorType &,
+                                      const VectorType &,
+                                      const double alpha) {
+      J = kappa + alpha;
+    };
+
+    // Used only in ver < 4.0.0
+    time_stepper.solve_jacobian_system =
+      [&](const VectorType &src, VectorType &dst) { dst[0] = src[0] / J; };
+
+    // Used in ver >= 4.0.0
+    time_stepper.solve_with_jacobian =
+      [&](const VectorType &src, VectorType &dst, const double) {
+        dst[0] = src[0] / J;
+      };
+
+    time_stepper.output_step = [&](const double       t,
+                                   const VectorType & sol,
+                                   const VectorType & sol_dot,
+                                   const unsigned int step_number) {
+      deallog << "Intermediate output:" << std::endl;
+      deallog << "  t =" << t << std::endl;
+      deallog << "  y =" << sol[0] << "  (exact: " << std::exp(-kappa * t)
+              << ')' << std::endl;
+      deallog << "  y'=" << sol_dot[0]
+              << "  (exact: " << -kappa * std::exp(-kappa * t) << ')'
+              << std::endl;
+    };
+  }
+
+  void
+  run()
+  {
+    y[0]     = 1;
+    y_dot[0] = -kappa;
+    try
+      {
+        time_stepper.solve_dae(y, y_dot);
+      }
+    catch (const std::exception &exc)
+      {
+        deallog << "Caught an irrecoverable exception in a callback:"
+                << std::endl
+                << exc.what() << std::endl;
+      }
+  }
+
+  SUNDIALS::IDA<Vector<double>> time_stepper;
+
+private:
+  Vector<double> y;
+  Vector<double> y_dot;
+  double         J;
+  double         kappa;
+  double         last_residual_eval_time;
+};
+
+
+int
+main()
+{
+  initlog();
+  deallog << std::setprecision(10);
+
+  SUNDIALS::IDA<Vector<double>>::AdditionalData data;
+  ParameterHandler                              prm;
+  data.add_parameters(prm);
+
+  // std::ofstream ofile(SOURCE_DIR "/ida_04.prm");
+  // prm.print_parameters(ofile, ParameterHandler::ShortText);
+  // ofile.close();
+
+  std::ifstream ifile(SOURCE_DIR "/ida_04_in.prm");
+  prm.parse_input(ifile);
+
+
+  HarmonicOscillator ode(2.0, data);
+  ode.run();
+}
diff --git a/tests/sundials/ida_05.output b/tests/sundials/ida_05.output
new file mode 100644 (file)
index 0000000..2422eb3
--- /dev/null
@@ -0,0 +1,282 @@
+
+DEAL::Intermediate output:
+DEAL::  t =0.000000000
+DEAL::  y =1.000000000  (exact: 1.000000000)
+DEAL::  y'=-2.000000000  (exact: -2.000000000)
+DEAL::Evaluating residual at t=0.2000000000
+DEAL::Time step is too large: 0.2000000000. Failing recoverably.
+DEAL::Evaluating residual at t=0.05000000000
+DEAL::Evaluating residual at t=0.05000000000
+DEAL::Evaluating residual at t=0.01250000000
+DEAL::Evaluating residual at t=0.01250000000
+DEAL::Evaluating residual at t=0.003125000000
+DEAL::Evaluating residual at t=0.003125000000
+DEAL::Evaluating residual at t=0.0007812500000
+DEAL::Evaluating residual at t=0.0007812500000
+DEAL::Evaluating residual at t=0.0001953125000
+DEAL::Evaluating residual at t=0.0001953125000
+DEAL::Evaluating residual at t=4.882812500e-05
+DEAL::Evaluating residual at t=4.882812500e-05
+DEAL::Evaluating residual at t=1.220703125e-05
+DEAL::Evaluating residual at t=1.220703125e-05
+DEAL::Evaluating residual at t=3.051757813e-06
+DEAL::Evaluating residual at t=3.051757813e-06
+DEAL::Evaluating residual at t=9.155273438e-06
+DEAL::Evaluating residual at t=9.155273438e-06
+DEAL::Evaluating residual at t=1.525878906e-05
+DEAL::Evaluating residual at t=2.136230469e-05
+DEAL::Evaluating residual at t=3.356933594e-05
+DEAL::Evaluating residual at t=3.356933594e-05
+DEAL::Evaluating residual at t=4.577636719e-05
+DEAL::Evaluating residual at t=5.798339844e-05
+DEAL::Evaluating residual at t=7.019042969e-05
+DEAL::Evaluating residual at t=9.460449219e-05
+DEAL::Evaluating residual at t=0.0001434326172
+DEAL::Evaluating residual at t=0.0001434326172
+DEAL::Evaluating residual at t=0.0002410888672
+DEAL::Evaluating residual at t=0.0002410888672
+DEAL::Evaluating residual at t=0.0004364013672
+DEAL::Evaluating residual at t=0.0004364013672
+DEAL::Evaluating residual at t=0.0006317138672
+DEAL::Evaluating residual at t=0.0008270263672
+DEAL::Evaluating residual at t=0.001022338867
+DEAL::Evaluating residual at t=0.001412963867
+DEAL::Evaluating residual at t=0.001412963867
+DEAL::Evaluating residual at t=0.001803588867
+DEAL::Evaluating residual at t=0.002194213867
+DEAL::Evaluating residual at t=0.002975463867
+DEAL::Evaluating residual at t=0.002975463867
+DEAL::Evaluating residual at t=0.003756713867
+DEAL::Evaluating residual at t=0.004537963867
+DEAL::Evaluating residual at t=0.005319213867
+DEAL::Evaluating residual at t=0.006881713867
+DEAL::Evaluating residual at t=0.006881713867
+DEAL::Evaluating residual at t=0.008444213867
+DEAL::Evaluating residual at t=0.01000671387
+DEAL::Evaluating residual at t=0.01156921387
+DEAL::Evaluating residual at t=0.01313171387
+DEAL::Evaluating residual at t=0.01625671387
+DEAL::Evaluating residual at t=0.01625671387
+DEAL::Evaluating residual at t=0.01938171387
+DEAL::Evaluating residual at t=0.02250671387
+DEAL::Evaluating residual at t=0.02563171387
+DEAL::Evaluating residual at t=0.02875671387
+DEAL::Evaluating residual at t=0.03500671387
+DEAL::Evaluating residual at t=0.03500671387
+DEAL::Evaluating residual at t=0.04125671387
+DEAL::Evaluating residual at t=0.04750671387
+DEAL::Evaluating residual at t=0.05375671387
+DEAL::Evaluating residual at t=0.06000671387
+DEAL::Evaluating residual at t=0.06625671387
+DEAL::Evaluating residual at t=0.07250671387
+DEAL::Evaluating residual at t=0.07250671387
+DEAL::Evaluating residual at t=0.07875671387
+DEAL::Evaluating residual at t=0.08500671387
+DEAL::Evaluating residual at t=0.09125671387
+DEAL::Evaluating residual at t=0.09750671387
+DEAL::Evaluating residual at t=0.1100067139
+DEAL::Evaluating residual at t=0.1100067139
+DEAL::Evaluating residual at t=0.1212567139
+DEAL::Evaluating residual at t=0.1212567139
+DEAL::Evaluating residual at t=0.1325067139
+DEAL::Evaluating residual at t=0.1437567139
+DEAL::Evaluating residual at t=0.1550067139
+DEAL::Evaluating residual at t=0.1662567139
+DEAL::Evaluating residual at t=0.1775067139
+DEAL::Evaluating residual at t=0.1887567139
+DEAL::Evaluating residual at t=0.2000067139
+DEAL::Intermediate output:
+DEAL::  t =0.2000000000
+DEAL::  y =0.6703200462  (exact: 0.6703200460)
+DEAL::  y'=-1.340640091  (exact: -1.340640092)
+DEAL::Evaluating residual at t=0.2112567139
+DEAL::Evaluating residual at t=0.2225067139
+DEAL::Evaluating residual at t=0.2337567139
+DEAL::Evaluating residual at t=0.2450067139
+DEAL::Evaluating residual at t=0.2562567139
+DEAL::Evaluating residual at t=0.2675067139
+DEAL::Evaluating residual at t=0.2787567139
+DEAL::Evaluating residual at t=0.2900067139
+DEAL::Evaluating residual at t=0.3012567139
+DEAL::Evaluating residual at t=0.3125067139
+DEAL::Evaluating residual at t=0.3237567139
+DEAL::Evaluating residual at t=0.3350067139
+DEAL::Evaluating residual at t=0.3462567139
+DEAL::Evaluating residual at t=0.3575067139
+DEAL::Evaluating residual at t=0.3687567139
+DEAL::Evaluating residual at t=0.3800067139
+DEAL::Evaluating residual at t=0.3912567139
+DEAL::Evaluating residual at t=0.4025067139
+DEAL::Intermediate output:
+DEAL::  t =0.4000000000
+DEAL::  y =0.4493289645  (exact: 0.4493289641)
+DEAL::  y'=-0.8986579278  (exact: -0.8986579282)
+DEAL::Evaluating residual at t=0.4137567139
+DEAL::Evaluating residual at t=0.4250067139
+DEAL::Evaluating residual at t=0.4362567139
+DEAL::Evaluating residual at t=0.4475067139
+DEAL::Evaluating residual at t=0.4587567139
+DEAL::Evaluating residual at t=0.4700067139
+DEAL::Evaluating residual at t=0.4812567139
+DEAL::Evaluating residual at t=0.4925067139
+DEAL::Evaluating residual at t=0.5037567139
+DEAL::Evaluating residual at t=0.5150067139
+DEAL::Evaluating residual at t=0.5262567139
+DEAL::Evaluating residual at t=0.5375067139
+DEAL::Evaluating residual at t=0.5487567139
+DEAL::Evaluating residual at t=0.5600067139
+DEAL::Evaluating residual at t=0.5712567139
+DEAL::Evaluating residual at t=0.5825067139
+DEAL::Evaluating residual at t=0.5937567139
+DEAL::Evaluating residual at t=0.6050067139
+DEAL::Intermediate output:
+DEAL::  t =0.6000000000
+DEAL::  y =0.3011942124  (exact: 0.3011942119)
+DEAL::  y'=-0.6023884237  (exact: -0.6023884238)
+DEAL::Evaluating residual at t=0.6162567139
+DEAL::Evaluating residual at t=0.6275067139
+DEAL::Evaluating residual at t=0.6387567139
+DEAL::Evaluating residual at t=0.6500067139
+DEAL::Evaluating residual at t=0.6612567139
+DEAL::Evaluating residual at t=0.6725067139
+DEAL::Evaluating residual at t=0.6837567139
+DEAL::Evaluating residual at t=0.6950067139
+DEAL::Evaluating residual at t=0.7062567139
+DEAL::Evaluating residual at t=0.7175067139
+DEAL::Evaluating residual at t=0.7287567139
+DEAL::Evaluating residual at t=0.7400067139
+DEAL::Evaluating residual at t=0.7512567139
+DEAL::Evaluating residual at t=0.7625067139
+DEAL::Evaluating residual at t=0.7737567139
+DEAL::Evaluating residual at t=0.7850067139
+DEAL::Evaluating residual at t=0.7962567139
+DEAL::Evaluating residual at t=0.8075067139
+DEAL::Intermediate output:
+DEAL::  t =0.8000000000
+DEAL::  y =0.2018965185  (exact: 0.2018965180)
+DEAL::  y'=-0.4037930361  (exact: -0.4037930360)
+DEAL::Evaluating residual at t=0.8187567139
+DEAL::Evaluating residual at t=0.8300067139
+DEAL::Evaluating residual at t=0.8412567139
+DEAL::Evaluating residual at t=0.8525067139
+DEAL::Evaluating residual at t=0.8637567139
+DEAL::Evaluating residual at t=0.8750067139
+DEAL::Evaluating residual at t=0.8862567139
+DEAL::Evaluating residual at t=0.8975067139
+DEAL::Evaluating residual at t=0.9087567139
+DEAL::Evaluating residual at t=0.9200067139
+DEAL::Evaluating residual at t=0.9312567139
+DEAL::Evaluating residual at t=0.9425067139
+DEAL::Evaluating residual at t=0.9537567139
+DEAL::Evaluating residual at t=0.9650067139
+DEAL::Evaluating residual at t=0.9762567139
+DEAL::Evaluating residual at t=0.9875067139
+DEAL::Evaluating residual at t=0.9987567139
+DEAL::Evaluating residual at t=1.010006714
+DEAL::Intermediate output:
+DEAL::  t =1.000000000
+DEAL::  y =0.1353352836  (exact: 0.1353352832)
+DEAL::  y'=-0.2706705668  (exact: -0.2706705665)
+DEAL::Evaluating residual at t=1.021256714
+DEAL::Evaluating residual at t=1.032506714
+DEAL::Evaluating residual at t=1.043756714
+DEAL::Evaluating residual at t=1.055006714
+DEAL::Evaluating residual at t=1.066256714
+DEAL::Evaluating residual at t=1.077506714
+DEAL::Evaluating residual at t=1.088756714
+DEAL::Evaluating residual at t=1.100006714
+DEAL::Evaluating residual at t=1.111256714
+DEAL::Evaluating residual at t=1.122506714
+DEAL::Evaluating residual at t=1.133756714
+DEAL::Evaluating residual at t=1.145006714
+DEAL::Evaluating residual at t=1.156256714
+DEAL::Evaluating residual at t=1.167506714
+DEAL::Evaluating residual at t=1.178756714
+DEAL::Evaluating residual at t=1.190006714
+DEAL::Evaluating residual at t=1.201256714
+DEAL::Intermediate output:
+DEAL::  t =1.200000000
+DEAL::  y =0.09071795362  (exact: 0.09071795329)
+DEAL::  y'=-0.1814359070  (exact: -0.1814359066)
+DEAL::Evaluating residual at t=1.212506714
+DEAL::Evaluating residual at t=1.223756714
+DEAL::Evaluating residual at t=1.235006714
+DEAL::Evaluating residual at t=1.246256714
+DEAL::Evaluating residual at t=1.257506714
+DEAL::Evaluating residual at t=1.268756714
+DEAL::Evaluating residual at t=1.280006714
+DEAL::Evaluating residual at t=1.291256714
+DEAL::Evaluating residual at t=1.302506714
+DEAL::Evaluating residual at t=1.313756714
+DEAL::Evaluating residual at t=1.325006714
+DEAL::Evaluating residual at t=1.336256714
+DEAL::Evaluating residual at t=1.347506714
+DEAL::Evaluating residual at t=1.358756714
+DEAL::Evaluating residual at t=1.370006714
+DEAL::Evaluating residual at t=1.381256714
+DEAL::Evaluating residual at t=1.392506714
+DEAL::Evaluating residual at t=1.403756714
+DEAL::Intermediate output:
+DEAL::  t =1.400000000
+DEAL::  y =0.06081006289  (exact: 0.06081006263)
+DEAL::  y'=-0.1216201256  (exact: -0.1216201253)
+DEAL::Evaluating residual at t=1.415006714
+DEAL::Evaluating residual at t=1.426256714
+DEAL::Evaluating residual at t=1.437506714
+DEAL::Evaluating residual at t=1.448756714
+DEAL::Evaluating residual at t=1.460006714
+DEAL::Evaluating residual at t=1.471256714
+DEAL::Evaluating residual at t=1.482506714
+DEAL::Evaluating residual at t=1.493756714
+DEAL::Evaluating residual at t=1.505006714
+DEAL::Evaluating residual at t=1.516256714
+DEAL::Evaluating residual at t=1.527506714
+DEAL::Evaluating residual at t=1.538756714
+DEAL::Evaluating residual at t=1.550006714
+DEAL::Evaluating residual at t=1.561256714
+DEAL::Evaluating residual at t=1.572506714
+DEAL::Evaluating residual at t=1.583756714
+DEAL::Evaluating residual at t=1.595006714
+DEAL::Evaluating residual at t=1.606256714
+DEAL::Intermediate output:
+DEAL::  t =1.600000000
+DEAL::  y =0.04076220418  (exact: 0.04076220398)
+DEAL::  y'=-0.08152440821  (exact: -0.08152440796)
+DEAL::Evaluating residual at t=1.617506714
+DEAL::Evaluating residual at t=1.628756714
+DEAL::Evaluating residual at t=1.640006714
+DEAL::Evaluating residual at t=1.651256714
+DEAL::Evaluating residual at t=1.662506714
+DEAL::Evaluating residual at t=1.673756714
+DEAL::Evaluating residual at t=1.685006714
+DEAL::Evaluating residual at t=1.707506714
+DEAL::Evaluating residual at t=1.707506714
+DEAL::Evaluating residual at t=1.727756714
+DEAL::Evaluating residual at t=1.727756714
+DEAL::Evaluating residual at t=1.748006714
+DEAL::Evaluating residual at t=1.768256714
+DEAL::Evaluating residual at t=1.788506714
+DEAL::Evaluating residual at t=1.808756714
+DEAL::Intermediate output:
+DEAL::  t =1.800000000
+DEAL::  y =0.02732372275  (exact: 0.02732372245)
+DEAL::  y'=-0.05464744359  (exact: -0.05464744489)
+DEAL::Evaluating residual at t=1.829006714
+DEAL::Evaluating residual at t=1.849256714
+DEAL::Evaluating residual at t=1.869506714
+DEAL::Evaluating residual at t=1.889756714
+DEAL::Evaluating residual at t=1.910006714
+DEAL::Evaluating residual at t=1.930256714
+DEAL::Evaluating residual at t=1.950506714
+DEAL::Evaluating residual at t=1.970756714
+DEAL::Evaluating residual at t=1.991006714
+DEAL::Evaluating residual at t=2.000000000
+DEAL::Evaluating residual at t=2.000000000
+DEAL::Intermediate output:
+DEAL::  t =2.000000000
+DEAL::  y =0.01831563933  (exact: 0.01831563889)
+DEAL::  y'=-0.03663127916  (exact: -0.03663127778)
+DEAL::Intermediate output:
+DEAL::  t =2.000000000
+DEAL::  y =0.01831563933  (exact: 0.01831563889)
+DEAL::  y'=-0.03663127916  (exact: -0.03663127778)

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.