--- /dev/null
+CMAKE_MINIMUM_REQUIRED(VERSION 2.8.9)
+INCLUDE(../setup_testsubproject.cmake)
+PROJECT(testsuite CXX)
+IF(DEAL_II_WITH_ADOLC)
+ DEAL_II_PICKUP_TESTS()
+ENDIF()
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 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:
+// - Taped doubles
+// - Gradient and hessian computations
+// Adapted from https://github.com/Homebrew/homebrew-science/blob/master/adol-c.rb
+
+#include "../tests.h"
+#include <adolc/adouble.h>
+#include <adolc/drivers/drivers.h>
+#include <adolc/taping.h>
+
+#include <math.h>
+
+int main(void)
+{
+ initlog();
+
+ 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;
+
+ // --- 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;
+
+ // --- 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;
+
+ // -- 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;
+}
--- /dev/null
+
+DEAL::Error (function 1): 0
+DEAL::Error (function 2): 0
+DEAL::Error (gradient): 0
+DEAL::Error (hessian): 0
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 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:
+// - Tapeless doubles
+// - Function and gradient computations
+// Adapted from https://github.com/Homebrew/homebrew-science/blob/master/adol-c.rb
+// See Adol-C manual v2.6 p65 figure 8
+
+#include "../tests.h"
+#include <adolc/adtl.h>
+#include <adolc/drivers/drivers.h>
+
+#include <math.h>
+
+int main(void)
+{
+ initlog();
+
+ const unsigned int n = 10;
+ adtl::setNumDir(n);
+
+ adtl::adouble *x = new adtl::adouble[n];
+ for (unsigned int i = 0; i < n; i++)
+ {
+ x[i] = (i + 1.0) / (2.0 + i);
+ x[i].setADValue(i,1);
+ }
+
+ adtl::adouble y = 1.0;
+ for (unsigned int i = 0; i < n; i++)
+ y *= x[i];
+
+ // --- Function ---
+
+ const double error_func = y.getValue() - 1.0 / (1.0 + n);
+
+ deallog << "Error (function): " << error_func << std::endl;
+
+ // --- Gradient ---
+
+ double err_grad = 0;
+ for (unsigned int i = 0; i < n; i++)
+ err_grad += std::abs(y.getADValue(i) - y.getValue() / x[i].getValue());
+
+ deallog << "Error (gradient): " << err_grad << std::endl;
+
+ // -- Cleanup ---
+
+ delete[] x;
+
+ return 0;
+}
--- /dev/null
+
+DEAL::Error (function): 0
+DEAL::Error (gradient): 0
--- /dev/null
+// ---------------------------------------------------------------------\r
+//\r
+// Copyright (C) 2016 - 2017 by the deal.II authors\r
+//\r
+// This file is part of the deal.II library.\r
+//\r
+// The deal.II library is free software; you can use it, redistribute\r
+// it, and/or modify it under the terms of the GNU Lesser General\r
+// Public License as published by the Free Software Foundation; either\r
+// version 2.1 of the License, or (at your option) any later version.\r
+// The full text of the license can be found in the file LICENSE at\r
+// the top level of the deal.II distribution.\r
+//\r
+// ---------------------------------------------------------------------\r
+\r
+// Test of basic functionality:\r
+// - Taped doubles\r
+// - Multiple dependent functions and Jacobian computations\r
+// - Update to using C++ data structures where possible\r
+// - Test that the tracing of active dependent/independent variables still\r
+// works with these structures\r
+// Adapted from https://github.com/Homebrew/homebrew-science/blob/master/adol-c.rb\r
+\r
+#include "../tests.h"\r
+#include <adolc/adouble.h>\r
+#include <adolc/drivers/drivers.h>\r
+#include <adolc/taping.h>\r
+\r
+#include <math.h>\r
+\r
+void test_reset_vector_values (const bool reset_values, const int tape_index)\r
+{\r
+ const unsigned int m = 5; // Dependents\r
+ const unsigned int n = 10; // Independents\r
+\r
+ std::vector<double> xp (n, 0.0);\r
+ std::vector<double> yp (m, 0.0);\r
+ std::vector<adouble> x (n, 1.0); // Dep. variable values initially set here\r
+ std::vector<adouble> y (m, 1.0);\r
+\r
+ if (reset_values == false)\r
+ for (unsigned int i = 0; i < n; i++)\r
+ xp[i] = (i + 1.0) / (2.0 + i);\r
+\r
+ trace_on(tape_index);\r
+ for (unsigned int i = 0; i < n; i++)\r
+ {\r
+ x[i] <<= xp[i];\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ y[j] *= (j+1)*x[i];\r
+ }\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ y[j] >>= yp[j];\r
+\r
+ trace_off();\r
+ // tapestats(1, tape_stats);\r
+\r
+ // --- Change values ---\r
+ if (reset_values == true)\r
+ for (unsigned int i = 0; i < n; i++)\r
+ xp[i] = (i + 1.0) / (2.0 + i);\r
+\r
+ // --- Functions ---\r
+ double *f = new double[m];\r
+ function(tape_index, m, n, xp.data(), f);\r
+\r
+ deallog << "Evaluation points:" << std::endl;\r
+ for (unsigned int i = 0; i < n; ++i)\r
+ deallog\r
+ << " x[" << i << "]: " << xp[i]\r
+ << std::endl;\r
+\r
+ deallog << "Function values:" << std::endl;\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ deallog\r
+ << " f[" << j << "]: " << f[j]\r
+ << " y[" << j << "]: " << yp[j]\r
+ << std::endl;\r
+\r
+ // --- Jacobian ---\r
+\r
+ double **J = new double*[m];\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ J[j] = new double[n];\r
+\r
+ jacobian(tape_index, m, n, xp.data(), J);\r
+\r
+ deallog << "Function jacobian J:" << std::endl;\r
+ for (unsigned int j = 0; j < m; j++)\r
+ {\r
+ for (unsigned int i = 0; i < n; i++)\r
+ deallog << J[j][i] << (i<n-1?",":"");\r
+\r
+ deallog << std::endl;\r
+ }\r
+\r
+ // -- Cleanup ---\r
+\r
+ delete[] f;\r
+ f = nullptr;\r
+\r
+ for (unsigned int j = 0; j < m; j++)\r
+ delete[] J[j];\r
+ delete[] J;\r
+ J = nullptr;\r
+}\r
+\r
+int main(void)\r
+{\r
+ initlog();\r
+\r
+ deallog << "Setting dependent variables a priori" << std::endl;\r
+ test_reset_vector_values(false, 1); // This works\r
+ deallog << std::endl;\r
+ deallog << "Setting dependent variables after function definition" << std::endl;\r
+ test_reset_vector_values(true, 2); // This doesn't\r
+}\r
--- /dev/null
+
+DEAL::Setting dependent variables a priori
+DEAL::Evaluation points:
+DEAL:: x[0]: 0.500000
+DEAL:: x[1]: 0.666667
+DEAL:: x[2]: 0.750000
+DEAL:: x[3]: 0.800000
+DEAL:: x[4]: 0.833333
+DEAL:: x[5]: 0.857143
+DEAL:: x[6]: 0.875000
+DEAL:: x[7]: 0.888889
+DEAL:: x[8]: 0.900000
+DEAL:: x[9]: 0.909091
+DEAL::Function values:
+DEAL:: f[0]: 0.0909091 y[0]: 0.0909091
+DEAL:: f[1]: 93.0909 y[1]: 93.0909
+DEAL:: f[2]: 5368.09 y[2]: 5368.09
+DEAL:: f[3]: 95325.1 y[3]: 95325.1
+DEAL:: f[4]: 887784. y[4]: 887784.
+DEAL::Function jacobian J:
+DEAL::0.181818,0.136364,0.121212,0.113636,0.109091,0.106061,0.103896,0.102273,0.101010,0.100000
+DEAL::186.182,139.636,124.121,116.364,111.709,108.606,106.390,104.727,103.434,102.400
+DEAL::10736.2,8052.14,7157.45,6710.11,6441.71,6262.77,6134.96,6039.10,5964.55,5904.90
+DEAL::190650.,142988.,127100.,119156.,114390.,111213.,108943.,107241.,105917.,104858.
+DEAL::1.77557e+06,1.33168e+06,1.18371e+06,1.10973e+06,1.06534e+06,1.03575e+06,1.01461e+06,998757.,986427.,976562.
+DEAL::
+DEAL::Setting dependent variables after function definition
+DEAL::Evaluation points:
+DEAL:: x[0]: 0.500000
+DEAL:: x[1]: 0.666667
+DEAL:: x[2]: 0.750000
+DEAL:: x[3]: 0.800000
+DEAL:: x[4]: 0.833333
+DEAL:: x[5]: 0.857143
+DEAL:: x[6]: 0.875000
+DEAL:: x[7]: 0.888889
+DEAL:: x[8]: 0.900000
+DEAL:: x[9]: 0.909091
+DEAL::Function values:
+DEAL:: f[0]: 0.0909091 y[0]: 0.00000
+DEAL:: f[1]: 93.0909 y[1]: 0.00000
+DEAL:: f[2]: 5368.09 y[2]: 0.00000
+DEAL:: f[3]: 95325.1 y[3]: 0.00000
+DEAL:: f[4]: 887784. y[4]: 0.00000
+DEAL::Function jacobian J:
+DEAL::0.181818,0.136364,0.121212,0.113636,0.109091,0.106061,0.103896,0.102273,0.101010,0.100000
+DEAL::186.182,139.636,124.121,116.364,111.709,108.606,106.390,104.727,103.434,102.400
+DEAL::10736.2,8052.14,7157.45,6710.11,6441.71,6262.77,6134.96,6039.10,5964.55,5904.90
+DEAL::190650.,142988.,127100.,119156.,114390.,111213.,108943.,107241.,105917.,104858.
+DEAL::1.77557e+06,1.33168e+06,1.18371e+06,1.10973e+06,1.06534e+06,1.03575e+06,1.01461e+06,998757.,986427.,976562.
--- /dev/null
+//\r
+// Copyright (C) 2016 - 2017 by the deal.II authors\r
+//\r
+// This file is part of the deal.II library.\r
+//\r
+// The deal.II library is free software; you can use it, redistribute\r
+// it, and/or modify it under the terms of the GNU Lesser General\r
+// Public License as published by the Free Software Foundation; either\r
+// version 2.1 of the License, or (at your option) any later version.\r
+// The full text of the license can be found in the file LICENSE at\r
+// the top level of the deal.II distribution.\r
+//\r
+// ---------------------------------------------------------------------\r
+\r
+// Test of basic functionality:\r
+// - Taped doubles\r
+// - Multiple dependent functions and Jacobian computations\r
+// Adapted from https://github.com/Homebrew/homebrew-science/blob/master/adol-c.rb\r
+\r
+#include "../tests.h"\r
+#include <adolc/adouble.h>\r
+#include <adolc/drivers/drivers.h>\r
+#include <adolc/taping.h>\r
+\r
+#include <math.h>\r
+\r
+int main(void)\r
+{\r
+ initlog();\r
+\r
+ const unsigned int m = 5; // Dependents\r
+ const unsigned int n = 10; // Independents\r
+ std::size_t tape_stats[STAT_SIZE];\r
+\r
+ double *xp = new double[n];\r
+ double *yp = new double[m];\r
+ adouble *x = new adouble[n];\r
+ adouble *y = new adouble[m];\r
+\r
+ for (unsigned int i = 0; i < n; i++)\r
+ xp[i] = (i + 1.0) / (2.0 + i);\r
+\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ y[j] = 1.0;\r
+\r
+ trace_on(1);\r
+ for (unsigned int i = 0; i < n; i++)\r
+ {\r
+ x[i] <<= xp[i];\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ y[j] *= (j+1)*x[i];\r
+ }\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ y[j] >>= yp[j];\r
+\r
+ delete[] x;\r
+ delete[] y;\r
+\r
+ trace_off();\r
+ tapestats(1, tape_stats);\r
+\r
+ // --- Functions ---\r
+ double *f = new double[m];\r
+ function(1, m, n, xp, f);\r
+\r
+ deallog << "Evaluation points:" << std::endl;\r
+ for (unsigned int i = 0; i < n; ++i)\r
+ deallog\r
+ << " x[" << i << "]: " << xp[i]\r
+ << std::endl;\r
+\r
+ deallog << "Function values:" << std::endl;\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ deallog\r
+ << " f[" << j << "]: " << f[j]\r
+ << " y[" << j << "]: " << yp[j]\r
+ << std::endl;\r
+\r
+ // --- Jacobian ---\r
+\r
+ double **J = new double*[m];\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ J[j] = new double[n];\r
+\r
+ jacobian(1, m, n, xp, J);\r
+\r
+ deallog << "Function jacobian J:" << std::endl;\r
+ for (unsigned int j = 0; j < m; j++)\r
+ {\r
+ for (unsigned int i = 0; i < n; i++)\r
+ deallog << J[j][i] << (i<n-1?",":"");\r
+\r
+ deallog << std::endl;\r
+ }\r
+\r
+ // -- Cleanup ---\r
+\r
+ delete[] xp;\r
+ delete[] yp;\r
+\r
+ delete[] f;\r
+ f = nullptr;\r
+\r
+ for (unsigned int j = 0; j < m; j++)\r
+ delete[] J[j];\r
+ delete[] J;\r
+ J = nullptr;\r
+\r
+ return 0;\r
+}\r
--- /dev/null
+
+DEAL::Evaluation points:
+DEAL:: x[0]: 0.500000
+DEAL:: x[1]: 0.666667
+DEAL:: x[2]: 0.750000
+DEAL:: x[3]: 0.800000
+DEAL:: x[4]: 0.833333
+DEAL:: x[5]: 0.857143
+DEAL:: x[6]: 0.875000
+DEAL:: x[7]: 0.888889
+DEAL:: x[8]: 0.900000
+DEAL:: x[9]: 0.909091
+DEAL::Function values:
+DEAL:: f[0]: 0.0909091 y[0]: 0.0909091
+DEAL:: f[1]: 93.0909 y[1]: 93.0909
+DEAL:: f[2]: 5368.09 y[2]: 5368.09
+DEAL:: f[3]: 95325.1 y[3]: 95325.1
+DEAL:: f[4]: 887784. y[4]: 887784.
+DEAL::Function jacobian J:
+DEAL::0.181818,0.136364,0.121212,0.113636,0.109091,0.106061,0.103896,0.102273,0.101010,0.100000
+DEAL::186.182,139.636,124.121,116.364,111.709,108.606,106.390,104.727,103.434,102.400
+DEAL::10736.2,8052.14,7157.45,6710.11,6441.71,6262.77,6134.96,6039.10,5964.55,5904.90
+DEAL::190650.,142988.,127100.,119156.,114390.,111213.,108943.,107241.,105917.,104858.
+DEAL::1.77557e+06,1.33168e+06,1.18371e+06,1.10973e+06,1.06534e+06,1.03575e+06,1.01461e+06,998757.,986427.,976562.
--- /dev/null
+// ---------------------------------------------------------------------\r
+//\r
+// Copyright (C) 2016 - 2017 by the deal.II authors\r
+//\r
+// This file is part of the deal.II library.\r
+//\r
+// The deal.II library is free software; you can use it, redistribute\r
+// it, and/or modify it under the terms of the GNU Lesser General\r
+// Public License as published by the Free Software Foundation; either\r
+// version 2.1 of the License, or (at your option) any later version.\r
+// The full text of the license can be found in the file LICENSE at\r
+// the top level of the deal.II distribution.\r
+//\r
+// ---------------------------------------------------------------------\r
+\r
+// Test of basic functionality:\r
+// - Tapeless doubles\r
+// - Multiple dependent functions and Jacobian computations\r
+// Adapted from https://github.com/Homebrew/homebrew-science/blob/master/adol-c.rb\r
+// See Adol-C manual v2.6 p65 figure 8\r
+\r
+#include "../tests.h"\r
+#include <adolc/adtl.h>\r
+#include <adolc/drivers/drivers.h>\r
+\r
+#include <math.h>\r
+\r
+int main(void)\r
+{\r
+ initlog();\r
+\r
+ const unsigned int m = 5; // Dependents\r
+ const unsigned int n = 10; // Independents\r
+ adtl::setNumDir(n);\r
+\r
+ adtl::adouble *x = new adtl::adouble[n];\r
+ for (unsigned int i = 0; i < n; i++)\r
+ {\r
+ x[i] = (i + 1.0) / (2.0 + i);\r
+ x[i].setADValue(i,1);\r
+ }\r
+\r
+ adtl::adouble *y = new adtl::adouble[m];\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ y[j] = 1.0;\r
+\r
+ for (unsigned int i = 0; i < n; i++)\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ y[j] *= (j+1)*x[i];\r
+\r
+ // --- Functions ---\r
+\r
+ deallog << "Evaluation points:" << std::endl;\r
+ for (unsigned int i = 0; i < n; ++i)\r
+ deallog\r
+ << " x[" << i << "]: " << x[i].getValue()\r
+ << std::endl;\r
+\r
+ deallog << "Function values:" << std::endl;\r
+ for (unsigned int j = 0; j < m; ++j)\r
+ deallog\r
+ << " y[" << j << "]: " << y[j].getValue()\r
+ << std::endl;\r
+\r
+ // --- Jacobian ---\r
+\r
+ deallog << "Function jacobian J:" << std::endl;\r
+ for (unsigned int j = 0; j < m; j++)\r
+ {\r
+ for (unsigned int i = 0; i < n; i++)\r
+ deallog << y[j].getADValue(i) << (i<n-1?",":"");\r
+\r
+ deallog << std::endl;\r
+ }\r
+\r
+ // -- Cleanup ---\r
+\r
+ delete[] x;\r
+ delete[] y;\r
+\r
+ return 0;\r
+}\r
--- /dev/null
+
+DEAL::Evaluation points:
+DEAL:: x[0]: 0.500000
+DEAL:: x[1]: 0.666667
+DEAL:: x[2]: 0.750000
+DEAL:: x[3]: 0.800000
+DEAL:: x[4]: 0.833333
+DEAL:: x[5]: 0.857143
+DEAL:: x[6]: 0.875000
+DEAL:: x[7]: 0.888889
+DEAL:: x[8]: 0.900000
+DEAL:: x[9]: 0.909091
+DEAL::Function values:
+DEAL:: y[0]: 0.0909091
+DEAL:: y[1]: 93.0909
+DEAL:: y[2]: 5368.09
+DEAL:: y[3]: 95325.1
+DEAL:: y[4]: 887784.
+DEAL::Function jacobian J:
+DEAL::0.181818,0.136364,0.121212,0.113636,0.109091,0.106061,0.103896,0.102273,0.101010,0.100000
+DEAL::186.182,139.636,124.121,116.364,111.709,108.606,106.390,104.727,103.434,102.400
+DEAL::10736.2,8052.14,7157.45,6710.11,6441.71,6262.77,6134.96,6039.10,5964.55,5904.90
+DEAL::190650.,142988.,127100.,119156.,114390.,111213.,108943.,107241.,105917.,104858.
+DEAL::1.77557e+06,1.33168e+06,1.18371e+06,1.10973e+06,1.06534e+06,1.03575e+06,1.01461e+06,998757.,986427.,976562.
projects / dealii.git / commitdiff