]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Test all nonjacobian interfaces.
authorLuca Heltai <luca.heltai@sissa.it>
Wed, 27 Sep 2017 10:53:16 +0000 (12:53 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 27 Sep 2017 12:34:42 +0000 (14:34 +0200)
tests/sundials/harmonic_oscillator_02.cc [new file with mode: 0644]
tests/sundials/harmonic_oscillator_02.output [new file with mode: 0644]
tests/sundials/harmonic_oscillator_02.prm [new file with mode: 0644]
tests/sundials/harmonic_oscillator_03.cc [new file with mode: 0644]
tests/sundials/harmonic_oscillator_03.output [new file with mode: 0644]
tests/sundials/harmonic_oscillator_04.cc [new file with mode: 0644]
tests/sundials/harmonic_oscillator_04.output [new file with mode: 0644]
tests/sundials/harmonic_oscillator_04.prm [new file with mode: 0644]

diff --git a/tests/sundials/harmonic_oscillator_02.cc b/tests/sundials/harmonic_oscillator_02.cc
new file mode 100644 (file)
index 0000000..194cd2c
--- /dev/null
@@ -0,0 +1,105 @@
+//-----------------------------------------------------------
+//
+//    Copyright (C) 2017 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 at
+//    the top level of the deal.II distribution.
+//
+//-----------------------------------------------------------
+
+#include "../tests.h"
+#include <deal.II/sundials/arkode.h>
+#include <deal.II/base/parameter_handler.h>
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+// Test explicit time stepper. Only implements explicit_function.
+
+/**
+ * Solve the Harmonic oscillator problem.
+ *
+ * u'' = -k^2 u
+ * u (0) = 0
+ * u'(0) = k
+ *
+ * write in terms of a first order ode:
+ *
+ * y[0]' =       y[1]
+ * y[1]' = - k^2 y[0]
+ *
+ * That is
+ *
+ * y' = A y
+ *
+ * A = [ 0 , 1; -k^2, 0 ]
+ *
+ * y_0  = 0, k
+ *
+ * The exact solution is
+ *
+ * y[0](t) = sin(k t)
+ * y[1](t) = k cos(k t)
+ *
+ */
+int main (int argc, char **argv)
+{
+  std::ofstream out("output");
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, numbers::invalid_unsigned_int);
+
+  typedef Vector<double> VectorType;
+
+  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 "/harmonic_oscillator_02.prm");
+      prm.print_parameters(ofile, ParameterHandler::ShortText);
+      ofile.close();
+    }
+
+  std::ifstream ifile(SOURCE_DIR "/harmonic_oscillator_02.prm");
+  prm.parse_input(ifile);
+
+  SUNDIALS::ARKode<VectorType> ode(MPI_COMM_WORLD,data);
+
+  ode.reinit_vector = [&] (VectorType&v)
+  {
+    v.reinit(2);
+  };
+
+  double kappa = 1.0;
+
+  ode.explicit_function = [&] (double,
+                               const VectorType &y,
+                               VectorType &ydot) -> int
+  {
+    ydot[0] = y[1];
+    ydot[1] = -kappa*kappa*y[0];
+    return 0;
+  };
+
+  ode.output_step = [&](const double t,
+                        const VectorType &sol,
+                        const unsigned int step_number) -> int
+  {
+    out << t << " "
+    << sol[0] << " " << sol[1] << std::endl;
+    return 0;
+  };
+
+  Vector<double> y(2);
+  y[0] = 0;
+  y[1] = kappa;
+  ode.solve_ode(y);
+  return 0;
+}
diff --git a/tests/sundials/harmonic_oscillator_02.output b/tests/sundials/harmonic_oscillator_02.output
new file mode 100644 (file)
index 0000000..3571d90
--- /dev/null
@@ -0,0 +1,127 @@
+0 0 1
+0.05 0.0499791 0.99875
+0.1 0.0998334 0.995004
+0.15 0.149438 0.988771
+0.2 0.198669 0.980066
+0.25 0.247404 0.968911
+0.3 0.295516 0.955326
+0.35 0.34289 0.939353
+0.4 0.38941 0.92104
+0.45 0.43496 0.900433
+0.5 0.479424 0.877579
+0.55 0.522687 0.852524
+0.6 0.564637 0.825329
+0.65 0.605176 0.796071
+0.7 0.644206 0.764828
+0.75 0.681631 0.73168
+0.8 0.717355 0.696705
+0.85 0.751278 0.659982
+0.9 0.783303 0.621595
+0.95 0.813363 0.581652
+1 0.841396 0.540258
+1.05 0.867339 0.497522
+1.1 0.891131 0.453553
+1.15 0.912711 0.408458
+1.2 0.932015 0.362345
+1.25 0.948983 0.315321
+1.3 0.963554 0.267498
+1.35 0.975694 0.219003
+1.4 0.985391 0.169961
+1.45 0.992634 0.120495
+1.5 0.997413 0.07073
+1.55 0.999716 0.0207894
+1.6 0.999533 -0.0292026
+1.65 0.996854 -0.0791221
+1.7 0.991667 -0.128845
+1.75 0.983976 -0.178243
+1.8 0.973816 -0.227193
+1.85 0.961228 -0.275575
+1.9 0.946248 -0.323273
+1.95 0.928918 -0.370167
+2 0.909275 -0.41614
+2.05 0.88736 -0.461072
+2.1 0.86321 -0.504847
+2.15 0.836886 -0.547347
+2.2 0.808467 -0.588477
+2.25 0.778033 -0.62814
+2.3 0.745665 -0.666242
+2.35 0.711444 -0.702689
+2.4 0.675451 -0.737384
+2.45 0.637765 -0.770233
+2.5 0.598471 -0.801142
+2.55 0.557673 -0.830034
+2.6 0.515481 -0.856851
+2.65 0.472005 -0.881535
+2.7 0.427357 -0.904029
+2.75 0.381646 -0.924275
+2.8 0.334983 -0.942215
+2.85 0.287479 -0.957791
+2.9 0.239248 -0.970952
+2.95 0.190419 -0.981674
+3 0.141115 -0.989946
+3.05 0.091459 -0.995756
+3.1 0.0415761 -0.999091
+3.15 -0.00841 -0.999941
+3.2 -0.0583754 -0.998292
+3.25 -0.108196 -0.994133
+3.3 -0.157743 -0.987468
+3.35 -0.206893 -0.978328
+3.4 -0.255528 -0.96675
+3.45 -0.303527 -0.952769
+3.5 -0.350773 -0.93642
+3.55 -0.397144 -0.917738
+3.6 -0.442522 -0.89676
+3.65 -0.486788 -0.873522
+3.7 -0.529827 -0.848088
+3.75 -0.571542 -0.820533
+3.8 -0.611833 -0.790935
+3.85 -0.650604 -0.759371
+3.9 -0.687755 -0.725916
+3.95 -0.723188 -0.690649
+4 -0.756806 -0.653646
+4.05 -0.788517 -0.614997
+4.1 -0.818254 -0.574809
+4.15 -0.84595 -0.533188
+4.2 -0.871543 -0.490241
+4.25 -0.894968 -0.446074
+4.3 -0.916161 -0.400795
+4.35 -0.935058 -0.35451
+4.4 -0.951598 -0.307332
+4.45 -0.965751 -0.259384
+4.5 -0.977494 -0.210788
+4.55 -0.986805 -0.161669
+4.6 -0.993662 -0.112147
+4.65 -0.998044 -0.0623455
+4.7 -0.999928 -0.0123877
+4.75 -0.999294 0.0376024
+4.8 -0.996148 0.0874965
+4.85 -0.990514 0.137172
+4.9 -0.982413 0.186506
+4.95 -0.97187 0.235376
+5 -0.958908 0.28366
+5.05 -0.94355 0.331236
+5.1 -0.925819 0.37798
+5.15 -0.905761 0.423774
+5.2 -0.883439 0.468508
+5.25 -0.858915 0.512075
+5.3 -0.832253 0.554367
+5.35 -0.803517 0.595277
+5.4 -0.772769 0.634698
+5.45 -0.740077 0.672522
+5.5 -0.705531 0.708661
+5.55 -0.669225 0.743031
+5.6 -0.631253 0.775551
+5.65 -0.591709 0.80614
+5.7 -0.550686 0.834716
+5.75 -0.508281 0.861197
+5.8 -0.464599 0.885514
+5.85 -0.419756 0.907617
+5.9 -0.373868 0.927459
+5.95 -0.327048 0.944992
+6 -0.279413 0.960169
+6.05 -0.231078 0.972942
+6.1 -0.182163 0.983273
+6.15 -0.132791 0.991144
+6.2 -0.0830885 0.996543
+6.25 -0.0331784 0.999455
+6.28 -0.00318445 1
diff --git a/tests/sundials/harmonic_oscillator_02.prm b/tests/sundials/harmonic_oscillator_02.prm
new file mode 100644 (file)
index 0000000..d1da0f7
--- /dev/null
@@ -0,0 +1,13 @@
+set Final time                        = 6.28
+set Initial time                      = 0
+set Time interval between each output = 0.05
+subsection Error control
+  set Absolute error tolerance = 0.000001
+  set Relative error tolerance = 0.000010
+end
+subsection Running parameters
+  set Initial step size                      = 0.010000
+  set Maximum number of nonlinear iterations = 10
+  set Maximum order of ARK                   = 5
+  set Minimum step size                      = 0.000001
+end
diff --git a/tests/sundials/harmonic_oscillator_03.cc b/tests/sundials/harmonic_oscillator_03.cc
new file mode 100644 (file)
index 0000000..a32a28b
--- /dev/null
@@ -0,0 +1,98 @@
+//-----------------------------------------------------------
+//
+//    Copyright (C) 2015 by the deal2lkit authors
+//
+//    This file is part of the deal2lkit library.
+//
+//    The deal2lkit 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 at
+//    the top level of the deal2lkit distribution.
+//
+//-----------------------------------------------------------
+
+#include "../tests.h"
+#include <deal.II/sundials/arkode.h>
+#include <deal.II/base/parameter_handler.h>
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+// Test implicit time stepper, no jacobian. Only implements implicit_function.
+
+/**
+ * Solve the Harmonic oscillator problem.
+ *
+ * u'' = -k^2 u
+ * u (0) = 0
+ * u'(0) = k
+ *
+ * write in terms of a first order ode:
+ *
+ * y[0]' =       y[1]
+ * y[1]' = - k^2 y[0]
+ *
+ * That is
+ *
+ * y' = A y
+ *
+ * A = [ 0 , 1; -k^2, 0 ]
+ *
+ * y_0  = 0, k
+ *
+ * The exact solution is
+ *
+ * y[0](t) = sin(k t)
+ * y[1](t) = k cos(k t)
+ *
+ */
+int main (int argc, char **argv)
+{
+  std::ofstream out("output");
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, numbers::invalid_unsigned_int);
+
+  typedef Vector<double> VectorType;
+
+  ParameterHandler prm;
+  SUNDIALS::ARKode<VectorType>::AdditionalData data;
+  data.add_parameters(prm);
+
+  // Use the same parameters of test 2.
+  std::ifstream ifile(SOURCE_DIR "/harmonic_oscillator_02.prm");
+  prm.parse_input(ifile);
+
+  SUNDIALS::ARKode<VectorType> ode(MPI_COMM_WORLD,data);
+
+  ode.reinit_vector = [&] (VectorType&v)
+  {
+    v.reinit(2);
+  };
+
+  double kappa = 1.0;
+
+  ode.implicit_function = [&] (double,
+                               const VectorType &y,
+                               VectorType &ydot) -> int
+  {
+    ydot[0] = y[1];
+    ydot[1] = -kappa*kappa*y[0];
+    return 0;
+  };
+
+  ode.output_step = [&](const double t,
+                        const VectorType &sol,
+                        const unsigned int step_number) -> int
+  {
+    out << t << " "
+    << sol[0] << " " << sol[1] << std::endl;
+    return 0;
+  };
+
+  Vector<double> y(2);
+  y[0] = 0;
+  y[1] = kappa;
+  ode.solve_ode(y);
+  return 0;
+}
diff --git a/tests/sundials/harmonic_oscillator_03.output b/tests/sundials/harmonic_oscillator_03.output
new file mode 100644 (file)
index 0000000..194df9f
--- /dev/null
@@ -0,0 +1,127 @@
+0 0 1
+0.05 0.0499792 0.99875
+0.1 0.0998334 0.995004
+0.15 0.149437 0.988766
+0.2 0.198666 0.980053
+0.25 0.2474 0.968896
+0.3 0.295517 0.955325
+0.35 0.342897 0.93937
+0.4 0.389417 0.921059
+0.45 0.434946 0.900419
+0.5 0.479377 0.877514
+0.55 0.522611 0.852418
+0.6 0.564547 0.825205
+0.65 0.605086 0.795949
+0.7 0.644127 0.764724
+0.75 0.681571 0.731603
+0.8 0.717319 0.696661
+0.85 0.75127 0.659971
+0.9 0.783323 0.621608
+0.95 0.81334 0.581655
+1 0.841255 0.540225
+1.05 0.867034 0.497436
+1.1 0.890639 0.453406
+1.15 0.912037 0.408253
+1.2 0.931192 0.362094
+1.25 0.948069 0.315048
+1.3 0.962632 0.267233
+1.35 0.974846 0.218766
+1.4 0.984677 0.169766
+1.45 0.992087 0.12035
+1.5 0.997044 0.0706359
+1.55 0.99951 0.0207423
+1.6 0.999451 -0.0292132
+1.65 0.996832 -0.0791126
+1.7 0.991621 -0.128824
+1.75 0.983844 -0.178179
+1.8 0.973556 -0.227062
+1.85 0.960816 -0.275364
+1.9 0.945682 -0.322977
+1.95 0.928212 -0.369793
+2 0.908464 -0.415703
+2.05 0.886494 -0.460597
+2.1 0.862363 -0.504368
+2.15 0.836126 -0.546907
+2.2 0.807842 -0.588105
+2.25 0.77757 -0.627854
+2.3 0.745366 -0.666045
+2.35 0.711289 -0.702569
+2.4 0.675396 -0.737318
+2.45 0.637746 -0.770184
+2.5 0.598437 -0.801048
+2.55 0.557607 -0.829841
+2.6 0.515372 -0.856524
+2.65 0.47185 -0.88106
+2.7 0.427158 -0.903409
+2.75 0.381414 -0.923533
+2.8 0.334733 -0.941395
+2.85 0.287234 -0.956956
+2.9 0.239033 -0.970177
+2.95 0.190247 -0.981021
+3 0.140994 -0.989448
+3.05 0.0913914 -0.995421
+3.1 0.0415552 -0.998901
+3.15 -0.00839699 -0.99985
+3.2 -0.0583482 -0.998229
+3.25 -0.108147 -0.994018
+3.3 -0.157648 -0.987254
+3.35 -0.206738 -0.97799
+3.4 -0.255306 -0.966276
+3.45 -0.30324 -0.952166
+3.5 -0.350428 -0.93571
+3.55 -0.396758 -0.916961
+3.6 -0.442119 -0.89597
+3.65 -0.486399 -0.872789
+3.7 -0.529485 -0.847469
+3.75 -0.571267 -0.820064
+3.8 -0.611631 -0.790623
+3.85 -0.650467 -0.7592
+3.9 -0.687663 -0.725846
+3.95 -0.723107 -0.690613
+4 -0.756685 -0.653591
+4.05 -0.788325 -0.614912
+4.1 -0.817974 -0.574689
+4.15 -0.845582 -0.533032
+4.2 -0.871095 -0.490056
+4.25 -0.894462 -0.445871
+4.3 -0.915632 -0.40059
+4.35 -0.934552 -0.354325
+4.4 -0.95117 -0.307188
+4.45 -0.965435 -0.259293
+4.5 -0.977295 -0.21075
+4.55 -0.986699 -0.161672
+4.6 -0.993593 -0.112172
+4.65 -0.997932 -0.0623777
+4.7 -0.999719 -0.012438
+4.75 -0.998982 0.037526
+4.8 -0.995748 0.087394
+4.85 -0.990044 0.137045
+4.9 -0.9819 0.18636
+4.95 -0.971341 0.235217
+5 -0.958397 0.283495
+5.05 -0.943095 0.331076
+5.1 -0.925462 0.377837
+5.15 -0.905526 0.423659
+5.2 -0.883316 0.468421
+5.25 -0.858858 0.512003
+5.3 -0.832189 0.554277
+5.35 -0.803395 0.59512
+5.4 -0.772573 0.634451
+5.45 -0.739818 0.672191
+5.5 -0.705227 0.708262
+5.55 -0.668895 0.742585
+5.6 -0.630918 0.775083
+5.65 -0.591391 0.805678
+5.7 -0.55041 0.83429
+5.75 -0.508071 0.860841
+5.8 -0.464469 0.885255
+5.85 -0.419701 0.907451
+5.9 -0.373861 0.927352
+5.95 -0.327051 0.944886
+6 -0.279415 0.960036
+6.05 -0.23108 0.972786
+6.1 -0.18217 0.983113
+6.15 -0.132806 0.990994
+6.2 -0.0831124 0.996408
+6.25 -0.0332105 0.999333
+6.28 -0.00322054 0.999883
diff --git a/tests/sundials/harmonic_oscillator_04.cc b/tests/sundials/harmonic_oscillator_04.cc
new file mode 100644 (file)
index 0000000..729a7b3
--- /dev/null
@@ -0,0 +1,115 @@
+//-----------------------------------------------------------
+//
+//    Copyright (C) 2015 by the deal2lkit authors
+//
+//    This file is part of the deal2lkit library.
+//
+//    The deal2lkit 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 at
+//    the top level of the deal2lkit distribution.
+//
+//-----------------------------------------------------------
+
+#include "../tests.h"
+#include <deal.II/sundials/arkode.h>
+#include <deal.II/base/parameter_handler.h>
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+
+// Test implicit-explicit time stepper, no jacobian.
+// Brusselator benchmark
+
+/**
+ * This test problem is called "brusselator", and is a typical benchmark for
+ * ODE solvers. This problem has 3 dependent variables u, v and w, that depend
+ * on the independent variable t via the IVP system
+ *
+ * du/dt = a − (w + 1)u + v u^2
+ * dv/dt = w u − v u^2
+ * dw/dt = (b − w)/eps -w u
+ *
+ * We integrate over the interval 0 ≤ t ≤ 10, with the initial conditions
+ *
+ * u(0) = 3.9, v(0) = 1.1, w(0) = 2.8,
+ *
+ * and parameters
+ *
+ * a = 1.2, b = 2.5, and eps = 10−5
+ *
+ * The implicit part only contains the stiff part of the problem (the part with
+ * eps in right hand side of the third equation).
+ */
+int main (int argc, char **argv)
+{
+  std::ofstream out("output");
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, numbers::invalid_unsigned_int);
+
+  typedef Vector<double> VectorType;
+
+  ParameterHandler prm;
+  SUNDIALS::ARKode<VectorType>::AdditionalData data;
+  data.add_parameters(prm);
+
+  if (false)
+    {
+      std::ofstream ofile(SOURCE_DIR "/harmonic_oscillator_04.prm");
+      prm.print_parameters(ofile, ParameterHandler::ShortText);
+      ofile.close();
+    }
+
+  std::ifstream ifile(SOURCE_DIR "/harmonic_oscillator_04.prm");
+  prm.parse_input(ifile);
+
+  SUNDIALS::ARKode<VectorType> ode(MPI_COMM_WORLD,data);
+
+  ode.reinit_vector = [&] (VectorType&v)
+  {
+    // Three independent variables
+    v.reinit(3);
+  };
+
+  // Parameters
+  double u0 = 3.9, v0 = 1.1, w0 = 2.8, a = 1.2, b= 2.5, eps = 1e-5;
+
+  ode.implicit_function = [&] (double,
+                               const VectorType &y,
+                               VectorType &ydot) -> int
+  {
+    ydot[0] = 0;
+    ydot[1] = 0;
+    ydot[2] = (b-y[2])/eps;
+    return 0;
+  };
+
+
+  ode.explicit_function = [&] (double,
+                               const VectorType &y,
+                               VectorType &ydot) -> int
+  {
+    ydot[0] = a-(y[2]+1)*y[0]+y[1]*y[0]*y[0];
+    ydot[1] = y[2]*y[0]-y[1]*y[0]*y[0];
+    ydot[2] = -y[2]*y[0];
+    return 0;
+  };
+
+  ode.output_step = [&](const double t,
+                        const VectorType &sol,
+                        const unsigned int step_number) -> int
+  {
+    out << t << " "
+    << sol[0] << " " << sol[1] << " " << sol[2] << std::endl;
+    return 0;
+  };
+
+  Vector<double> y(3);
+  y[0] = u0;
+  y[1] = v0;
+  y[2] = w0;
+  ode.solve_ode(y);
+  return 0;
+}
diff --git a/tests/sundials/harmonic_oscillator_04.output b/tests/sundials/harmonic_oscillator_04.output
new file mode 100644 (file)
index 0000000..02827a7
--- /dev/null
@@ -0,0 +1,102 @@
+0 3.9 1.1 2.8
+0.1 3.99995 0.720432 2.49987
+0.2 3.78894 0.660967 2.49997
+0.3 3.52275 0.681501 2.49986
+0.4 3.26217 0.72296 2.49997
+0.5 3.01935 0.771864 2.4999
+0.6 2.79594 0.82467 2.49994
+0.7 2.59116 0.880241 2.49994
+0.8 2.4037 0.938097 2.49991
+0.9 2.23217 0.997959 2.49997
+1 2.07527 1.0596 2.49991
+1.1 1.93182 1.12281 2.49998
+1.2 1.80073 1.18737 2.49993
+1.3 1.68106 1.25304 2.49998
+1.4 1.57194 1.31959 2.49995
+1.5 1.47261 1.38677 2.49997
+1.6 1.38237 1.45434 2.49997
+1.7 1.3006 1.52203 2.49996
+1.8 1.2267 1.58963 2.49998
+1.9 1.16015 1.65689 2.49996
+2 1.10046 1.72361 2.49999
+2.1 1.04715 1.78959 2.49996
+2.2 0.999783 1.85466 2.49999
+2.3 0.957942 1.91866 2.49997
+2.4 0.921233 1.98145 2.49998
+2.5 0.889283 2.04291 2.49998
+2.6 0.861743 2.10294 2.49998
+2.7 0.838291 2.16142 2.49998
+2.8 0.818628 2.21827 2.49998
+2.9 0.802485 2.27338 2.49998
+3 0.78962 2.32667 2.49998
+3.1 0.779824 2.37802 2.49998
+3.2 0.772919 2.42731 2.49998
+3.3 0.768758 2.47441 2.49999
+3.4 0.767229 2.51916 2.49998
+3.5 0.768257 2.56138 2.5
+3.6 0.7718 2.60085 2.50001
+3.7 0.777859 2.63733 2.49998
+3.8 0.786474 2.67052 2.49993
+3.9 0.797733 2.70008 2.4999
+4 0.811773 2.72558 2.49999
+4.1 0.828789 2.74657 2.50008
+4.2 0.84904 2.76245 2.50006
+4.3 0.872862 2.77257 2.49987
+4.4 0.900673 2.77611 2.49997
+4.5 0.932987 2.77216 2.4999
+4.6 0.970423 2.75959 2.50008
+4.7 1.01371 2.73716 2.50002
+4.8 1.06366 2.70339 2.49987
+4.9 1.12115 2.65673 2.50006
+5 1.187 2.59554 2.50009
+5.1 1.26179 2.51839 2.49993
+5.2 1.34554 2.42435 2.49999
+5.3 1.43719 2.31362 2.49988
+5.4 1.53404 2.18823 2.49988
+5.5 1.6313 2.05269 2.49994
+5.6 1.72208 1.91414 2.50001
+5.7 1.79854 1.7815 2.50002
+5.8 1.85371 1.66352 2.49999
+5.9 1.88345 1.56671 2.49994
+6 1.88731 1.4941 2.49988
+6.1 1.86803 1.44543 2.49982
+6.2 1.83015 1.41827 2.50004
+6.3 1.77867 1.40922 2.50005
+6.4 1.71817 1.41482 2.49984
+6.5 1.65242 1.43201 2.49985
+6.6 1.58434 1.45824 2.49989
+6.7 1.51607 1.49149 2.50006
+6.8 1.44916 1.53017 2.49996
+6.9 1.38465 1.57301 2.49987
+7 1.32327 1.61902 2.50004
+7.1 1.26548 1.6674 2.5
+7.2 1.21156 1.7175 2.4999
+7.3 1.16164 1.76879 2.5
+7.4 1.11578 1.82081 2.50002
+7.5 1.07396 1.87319 2.49993
+7.6 1.03609 1.92558 2.49996
+7.7 1.00208 1.97772 2.50002
+7.8 0.971784 2.02935 2.49996
+7.9 0.945072 2.08025 2.49995
+8 0.921787 2.13022 2.50001
+8.1 0.901781 2.17907 2.49998
+8.2 0.884906 2.22664 2.49996
+8.3 0.871025 2.27275 2.49999
+8.4 0.860013 2.31723 2.49998
+8.5 0.851761 2.35992 2.49997
+8.6 0.84618 2.40062 2.49998
+8.7 0.843199 2.43916 2.49998
+8.8 0.842774 2.47531 2.49999
+8.9 0.844885 2.50883 2.50001
+9 0.849543 2.53948 2.50008
+9.1 0.856788 2.56694 2.49996
+9.2 0.866696 2.59088 2.4999
+9.3 0.879382 2.61091 2.50016
+9.4 0.895002 2.6266 2.49985
+9.5 0.913755 2.63743 2.49992
+9.6 0.935893 2.64284 2.50001
+9.7 0.961713 2.64217 2.49985
+9.8 0.991567 2.63469 2.50001
+9.9 1.02585 2.61958 2.49985
+10 1.06497 2.59596 2.49997
+10 1.06497 2.59596 2.49997
diff --git a/tests/sundials/harmonic_oscillator_04.prm b/tests/sundials/harmonic_oscillator_04.prm
new file mode 100644 (file)
index 0000000..d0ea5d3
--- /dev/null
@@ -0,0 +1,15 @@
+set Final time                        = 10
+set Initial time                      = 0
+set Time interval between each output = 0.1
+subsection Error control
+  set Absolute error tolerance = 0.000001
+  set Relative error tolerance = 0.000010
+end
+subsection Running parameters
+  set Implicit function is linear            = true
+  set Implicit function is time independent  = true
+  set Initial step size                      = 0.010000
+  set Maximum number of nonlinear iterations = 10
+  set Maximum order of ARK                   = 5
+  set Minimum step size                      = 0.000001
+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.