]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a test.
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 15 May 2023 22:09:38 +0000 (16:09 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 18 May 2023 20:41:06 +0000 (14:41 -0600)
tests/sundials/arkode_09.cc [new file with mode: 0644]
tests/sundials/arkode_09.output [new file with mode: 0644]
tests/sundials/arkode_09_in.prm [new file with mode: 0644]

diff --git a/tests/sundials/arkode_09.cc b/tests/sundials/arkode_09.cc
new file mode 100644 (file)
index 0000000..33a6458
--- /dev/null
@@ -0,0 +1,92 @@
+//-----------------------------------------------------------
+//
+//    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/arkode.h>
+
+#include "../tests.h"
+
+// Test explicit time stepper where we throw a recoverable error from
+// the function that computes the right hand side. The test case
+// solves the problem y'(t)=-y(t) with exact solution
+// y(t)=exp(-t). The right hand side function throws an exception if
+// evaluated with y<0, and we make sure that that happens by starting
+// with an overly large initial time step.
+int
+main()
+{
+  initlog();
+
+  using VectorType = Vector<double>;
+
+  ParameterHandler                             prm;
+  SUNDIALS::ARKode<VectorType>::AdditionalData data;
+  data.add_parameters(prm);
+
+  // Set to true to reset input file.
+  if (false)
+    {
+      std::ofstream ofile(SOURCE_DIR "/arkode_09_in.prm");
+      prm.print_parameters(ofile, ParameterHandler::ShortText);
+      ofile.close();
+    }
+
+  std::ifstream ifile(SOURCE_DIR "/arkode_09_in.prm");
+  prm.parse_input(ifile);
+
+  SUNDIALS::ARKode<VectorType> ode(data);
+
+  double kappa = 1.0;
+
+  ode.explicit_function =
+    [&](const double t, const VectorType &y, VectorType &ydot) {
+      deallog << "Evaluating right hand side at t=" << t << " with y=" << y[0]
+              << std::endl;
+
+      // Error out in a recoverable way if asked to evaluate at a
+      // point where y<0. This can happen with explicit methods if the
+      // time step is too large.
+      if (y[0] < 0)
+        {
+          deallog << "Reporting a recoverable error." << std::endl;
+          throw RecoverableUserCallbackError();
+        }
+
+      ydot[0] = -y[0];
+    };
+
+  ode.output_step =
+    [&](const double t, const VectorType &sol, const unsigned int step_number) {
+      deallog << "Monitor called at t=" << t << ' ' << sol[0] << std::endl;
+    };
+
+  Vector<double> y(1);
+  y[0] = 1;
+
+  try
+    {
+      ode.solve_ode(y);
+    }
+  catch (const std::exception &exc)
+    {
+      deallog << "ARKODE did not succeed with the following error message:"
+              << std::endl
+              << exc.what() << std::endl;
+    }
+}
diff --git a/tests/sundials/arkode_09.output b/tests/sundials/arkode_09.output
new file mode 100644 (file)
index 0000000..3df1467
--- /dev/null
@@ -0,0 +1,24 @@
+
+DEAL::Monitor called at t=0.00000 1.00000
+DEAL::Evaluating right hand side at t=0.00000 with y=1.00000
+DEAL::Evaluating right hand side at t=0.00000 with y=1.00000
+DEAL::Evaluating right hand side at t=0.600000 with y=0.400000
+DEAL::Evaluating right hand side at t=0.900000 with y=0.505000
+DEAL::Evaluating right hand side at t=1.80000 with y=-0.638000
+DEAL::Reporting a recoverable error.
+DEAL::ARKODE did not succeed with the following error message:
+DEAL::
+--------------------------------------------------------
+An error occurred in line <0> of file <> in function
+    
+The violated condition was: 
+    
+Additional information: 
+    A user call-back function encountered a recoverable error, but the
+    underlying library that called the call-back did not manage to recover
+    from the error and aborted its operation.
+    
+    See the glossary entry on user call-back functions for more
+    information.
+--------------------------------------------------------
+
diff --git a/tests/sundials/arkode_09_in.prm b/tests/sundials/arkode_09_in.prm
new file mode 100644 (file)
index 0000000..59d346a
--- /dev/null
@@ -0,0 +1,13 @@
+set Final time                        = 6
+set Initial time                      = 0
+set Time interval between each output = 2
+subsection Error control
+  set Absolute error tolerance = 0.001
+  set Relative error tolerance = 0.010
+end
+subsection Running parameters
+  set Initial step size                      = 3
+  set Maximum number of nonlinear iterations = 10
+  set Maximum order of ARK                   = 5
+  set Minimum step size                      = 0.01
+end

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.