]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a quick test for Adol-C.
authorJean-Paul Pelteret <jppelteret@gmail.com>
Wed, 6 Sep 2017 10:09:50 +0000 (12:09 +0200)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Wed, 6 Sep 2017 10:09:50 +0000 (12:09 +0200)
tests/quick_tests/CMakeLists.txt
tests/quick_tests/adolc.cc [new file with mode: 0644]

index e6a3529a024e6c2bc91a39b1c5877f171661c9be..678d2d6d07ff5fc62fc5f313216b0620412e9c3a 100644 (file)
@@ -146,6 +146,11 @@ IF (DEAL_II_WITH_ARPACK AND DEAL_II_WITH_UMFPACK)
   make_quicktest("arpack" ${_mybuild} "")
 ENDIF()
 
+# Test Adol-C
+IF (DEAL_II_WITH_ADOLC)
+  make_quicktest("adolc" ${_mybuild} "")
+ENDIF()
+
 
 
 # A custom test target:
diff --git a/tests/quick_tests/adolc.cc b/tests/quick_tests/adolc.cc
new file mode 100644 (file)
index 0000000..bddd47e
--- /dev/null
@@ -0,0 +1,111 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// Test of basic functionality for Adol-C:
+//  - Taped doubles
+//  - Gradient and hessian computations
+// Adapted from https://github.com/Homebrew/homebrew-science/blob/master/adol-c.rb
+
+#include <deal.II/base/logstream.h>
+
+#include <adolc/adouble.h>
+#include <adolc/drivers/drivers.h>
+#include <adolc/taping.h>
+
+#include <math.h>
+
+using namespace dealii;
+
+int main(void)
+{
+  const double tol = 1e-12;
+  const unsigned int n = 10;
+  std::size_t tape_stats[STAT_SIZE];
+
+  double *xp = new double[n];
+  double  yp = 0.0;
+  adouble *x = new adouble[n];
+  adouble  y = 1.0;
+
+  for (unsigned int i = 0; i < n; i++)
+    xp[i] = (i + 1.0) / (2.0 + i);
+
+  trace_on(1);
+  for (unsigned int i = 0; i < n; i++)
+    {
+      x[i] <<= xp[i];
+      y *= x[i];
+    }
+  y >>= yp;
+  delete[] x;
+  trace_off();
+  tapestats(1, tape_stats);
+
+  // --- Function ---
+  double *f = new double;
+  function(1, 1, n, xp, f);
+
+  const double error_func_1 = yp - 1.0 / (1.0 + n);
+  const double error_func_2 = *f - 1.0 / (1.0 + n);
+
+  deallog << "Error (function 1): " << error_func_1 << std::endl;
+  deallog << "Error (function 2): " << error_func_2 << std::endl;
+  Assert(error_func_1 < tol, ExcMessage("Should be zero!"));
+  Assert(error_func_2 < tol, ExcMessage("Should be zero!"));
+
+  // --- Gradient ---
+
+  double *g = new double[n];
+  gradient(1, n, xp, g);
+
+  double err_grad = 0;
+  for (unsigned int i = 0; i < n; i++)
+    err_grad += std::abs(g[i] - yp / xp[i]);
+
+  deallog << "Error (gradient): " <<  err_grad << std::endl;
+  Assert(err_grad < tol, ExcMessage("Should be zero!"));
+
+  // --- Hessian ---
+
+  double **H = new double*[n];
+  for (unsigned int i = 0; i < n; ++i)
+    H[i] = new double[i+1];
+
+  hessian(1, n, xp, H);
+
+  double error_hess = 0;
+  for (unsigned int i = 0; i < n; i++)
+    for (unsigned int j = 0; j < n; j++)
+      if (i > j)
+        error_hess += std::abs(H[i][j] - g[i] / xp[j]);
+
+  deallog << "Error (hessian): " <<  error_hess << std::endl;
+  Assert(error_hess < tol, ExcMessage("Should be zero!"));
+
+  // -- Cleanup ---
+
+  delete f;
+  f = nullptr;
+
+  delete[] g;
+  g = nullptr;
+
+  for (unsigned int i = 0; i < n; i++)
+    delete[] H[i];
+  delete[] H;
+  H = nullptr;
+
+  return 0;
+}

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.