commit b51b72a7b51e30501ac1cee2974d17a2e75716cf Author: Zhuo Yang Date: Wed Apr 24 15:22:10 2024 +0800 add test for autodiff about forwarding func diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..7194ea7 --- /dev/null +++ b/.gitignore @@ -0,0 +1,2 @@ +.cache +build diff --git a/.vscode/launch.json b/.vscode/launch.json new file mode 100644 index 0000000..b15f2ff --- /dev/null +++ b/.vscode/launch.json @@ -0,0 +1,21 @@ +{ + // Use IntelliSense to learn about possible attributes. + // Hover to view descriptions of existing attributes. + // For more information, visit: https://go.microsoft.com/fwlink/?linkid=830387 + "version": "0.2.0", + "configurations": [ + { + "type": "vgdb", + "request": "launch", + "name": "C/C++ Debug", + "program": "${command:cmake.launchTargetPath}", + "args": [], + // "cwd": "${command:cmake.getLaunchTargetDirectory}", + "cwd": "${workspaceFolder}", + "env": { + "name": "PATH", + "value": "${env:PATH}:${command:cmake.getLaunchTargetDirectory}" + } + } + ] +} \ No newline at end of file diff --git a/CMakeLists.txt b/CMakeLists.txt new file mode 100644 index 0000000..34c7e7c --- /dev/null +++ b/CMakeLists.txt @@ -0,0 +1,24 @@ +cmake_minimum_required(VERSION 3.12) +set(CMAKE_CXX_STANDARD 17) +set(CMAKE_CXX_STANDARD_REQUIRED ON) + +set(CMAKE_EXPORT_COMPILE_COMMANDS ON) + +set(CMAKE_BUILD_TYPE Debug) + +set(conda_env_path "/share/software/anaconda3/envs/unisolver") +if(DEFINED CMAKE_PREFIX_PATH) + list(APPEND CMAKE_PREFIX_PATH ${conda_path}) +else() + set(CMAKE_PREFIX_PATH ${conda_env_path}) +endif() +message(STATUS "cmake_prefix_path: ${CMAKE_PREFIX_PATH}") + +project(test_autodiff LANGUAGES CXX) + +find_package(autodiff) +find_package(Eigen3 3.4 REQUIRED) +include_directories(${EIGEN3_INCLUDE_DIRS}) + + +add_subdirectory(test) diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt new file mode 100644 index 0000000..21b2b04 --- /dev/null +++ b/test/CMakeLists.txt @@ -0,0 +1,7 @@ +file(GLOB V_GLOB LIST_DIRECTORIES true "*") +# list(FILTER V_GLOB EXCLUDE REGEX ".*ThirdParty") +foreach(item ${V_GLOB}) + if(IS_DIRECTORY ${item}) + add_subdirectory(${item}) + endif() +endforeach() \ No newline at end of file diff --git a/test/forward/10jacobian_matrix_of_vector_func_using_memory_maps.cpp b/test/forward/10jacobian_matrix_of_vector_func_using_memory_maps.cpp new file mode 100644 index 0000000..f507536 --- /dev/null +++ b/test/forward/10jacobian_matrix_of_vector_func_using_memory_maps.cpp @@ -0,0 +1,37 @@ +// NOTE: 使用内存映射的向量函数的雅可比矩阵 + +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The vector function for which the Jacobian is needed +VectorXreal f(const VectorXreal& x) +{ + return x * x.sum(); +} + +int main() +{ + using Eigen::Map; + using Eigen::MatrixXd; + + VectorXreal x(5); // the input vector x with 5 variables + x << 1, 2, 3, 4, 5; // x = [1, 2, 3, 4, 5] + double y[25]; // the output Jacobian as a flat array + + VectorXreal F; // the output vector F = f(x) evaluated together with Jacobian matrix below + Map J(y, 5, 5); // the output Jacobian dF/dx mapped onto the flat array + + jacobian(f, wrt(x), at(x), F, J); // evaluate the output vector F and the Jacobian matrix dF/dx + + std::cout << "F = \n" << F << std::endl; // print the evaluated output vector F + std::cout << "J = \n" << J << std::endl; // print the evaluated Jacobian matrix dF/dx + std::cout << "y = "; // print the flat array + for(int i = 0 ; i < 25 ; ++i) + std::cout << y[i] << " "; + std::cout << std::endl; +} \ No newline at end of file diff --git a/test/forward/11higher_order_der_of_scalar_func.cpp b/test/forward/11higher_order_der_of_scalar_func.cpp new file mode 100644 index 0000000..ba3196f --- /dev/null +++ b/test/forward/11higher_order_der_of_scalar_func.cpp @@ -0,0 +1,40 @@ +// NOTE：标量函数的高阶交叉导 +// C++ includes +#include + +// autodiff include +#include +using namespace autodiff; + +// The multi-variable function for which higher-order derivatives are needed (up to 4th order) +dual4th f(dual4th x, dual4th y, dual4th z) +{ + return 1 + x + y + z + x*y + y*z + x*z + x*y*z + exp(x/y + y/z); +} + +int main() +{ + dual4th x = 1.0; + dual4th y = 2.0; + dual4th z = 3.0; + + auto [u0, ux, uxy, uxyx, uxyxz] = derivatives(f, wrt(x, y, x, z), at(x, y, z)); + + std::cout << "u0 = " << u0 << std::endl; // print the evaluated value of u = f(x, y, z) + std::cout << "ux = " << ux << std::endl; // print the evaluated derivative du/dx + std::cout << "uxy = " << uxy << std::endl; // print the evaluated derivative d²u/dxdy + std::cout << "uxyx = " << uxyx << std::endl; // print the evaluated derivative d³u/dx²dy + std::cout << "uxyxz = " << uxyxz << std::endl; // print the evaluated derivative d⁴u/dx²dydz +} + +/*------------------------------------------------------------------------------------------------- +=== Note === +--------------------------------------------------------------------------------------------------- +In most cases, dual can be replaced by real, as commented in other examples. +However, computing higher-order cross derivatives has definitely to be done +using higher-order dual types (e.g., dual3rd, dual4th)! This is because real +types (e.g., real2nd, real3rd, real4th) are optimally designed for computing +higher-order directional derivatives. +对于计算高阶交叉导数,必须使用更高阶的 dual 类型(如 dual3rd、dual4th 等)。 +否则可以使用real +-------------------------------------------------------------------------------------------------*/ \ No newline at end of file diff --git a/test/forward/12higher_order_directional_der_of_scalar_func.cpp b/test/forward/12higher_order_directional_der_of_scalar_func.cpp new file mode 100644 index 0000000..b834c8e --- /dev/null +++ b/test/forward/12higher_order_directional_der_of_scalar_func.cpp @@ -0,0 +1,36 @@ +// NOTE: 标量函数的高阶方向导数 +// C++ includes +#include + +// autodiff include +#include +using namespace autodiff; + +// The multi-variable function for which higher-order derivatives are needed (up to 4th order) +real4th f(real4th x, real4th y, real4th z) +{ + return sin(x) * cos(y) * exp(z); +} + +int main() +{ + real4th x = 1.0; + real4th y = 2.0; + real4th z = 3.0; + + auto dfdv = derivatives(f, along(1.0, 1.0, 2.0), at(x, y, z)); // the directional derivatives of f along direction v = (1, 1, 2) at (x, y, z) = (1, 2, 3) + + std::cout << "dfdv[0] = " << dfdv[0] << std::endl; // print the evaluated 0th order directional derivative of f along v (equivalent to f(x, y, z)) + std::cout << "dfdv[1] = " << dfdv[1] << std::endl; // print the evaluated 1st order directional derivative of f along v + std::cout << "dfdv[2] = " << dfdv[2] << std::endl; // print the evaluated 2nd order directional derivative of f along v + std::cout << "dfdv[3] = " << dfdv[3] << std::endl; // print the evaluated 3rd order directional derivative of f along v + std::cout << "dfdv[4] = " << dfdv[4] << std::endl; // print the evaluated 4th order directional derivative of f along v +} + +/*------------------------------------------------------------------------------------------------- +=== Note === +--------------------------------------------------------------------------------------------------- +This example would also work if dual was used instead of real. However, real +types are your best option for directional derivatives, as they were optimally +designed for this kind of derivatives. +-------------------------------------------------------------------------------------------------*/ \ No newline at end of file diff --git a/test/forward/13higher_order_directional_der_of_vector_func.cpp b/test/forward/13higher_order_directional_der_of_vector_func.cpp new file mode 100644 index 0000000..15aa2bc --- /dev/null +++ b/test/forward/13higher_order_directional_der_of_vector_func.cpp @@ -0,0 +1,73 @@ +// NOTE: 向量的高阶方向导数 +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The vector function for which higher-order directional derivatives are needed (up to 4th order). +ArrayXreal4th f(const ArrayXreal4th& x, real4th p) +{ + return p * x.log(); +} + +int main() +{ + using Eigen::ArrayXd; + + ArrayXreal4th x(5); // the input vector x + x << 1.0, 2.0, 3.0, 4.0, 5.0; + + real4th p = 7.0; // the input parameter p = 1 + + ArrayXd vx(5); // the direction vx in the direction vector v = (vx, vp) + vx << 1.0, 1.0, 1.0, 1.0, 1.0; + + double vp = 1.0; // the direction vp in the direction vector v = (vx, vp) + + auto dfdv = derivatives(f, along(vx, vp), at(x, p)); // the directional derivatives of f along direction v = (vx, vp) at (x, p) + + std::cout << std::scientific << std::showpos; + std::cout << "Directional derivatives of f along v = (vx, vp) up to 4th order:" << std::endl; + std::cout << "dfdv[0] = " << dfdv[0].transpose() << std::endl; // print the evaluated 0th order directional derivative of f along v (equivalent to f(x, p)) + std::cout << "dfdv[1] = " << dfdv[1].transpose() << std::endl; // print the evaluated 1st order directional derivative of f along v + std::cout << "dfdv[2] = " << dfdv[2].transpose() << std::endl; // print the evaluated 2nd order directional derivative of f along v + std::cout << "dfdv[3] = " << dfdv[3].transpose() << std::endl; // print the evaluated 3rd order directional derivative of f along v + std::cout << "dfdv[4] = " << dfdv[4].transpose() << std::endl; // print the evaluated 4th order directional derivative of f along v + std::cout << std::endl; + + double t = 0.1; // the step length along direction v = (vx, vp) used to compute 4th order Taylor estimate of f + + ArrayXreal4th u = f(x + t * vx, p + t * vp); // start from (x, p), walk a step length t = 0.1 along direction v = (vx, vp) and evaluate f at this current point + + ArrayXd utaylor = dfdv[0] + t*dfdv[1] + (t*t/2.0)*dfdv[2] + (t*t*t/6.0)*dfdv[3] + (t*t*t*t/24.0)*dfdv[4]; // evaluate the 4th order Taylor estimate of f along direction v = (vx, vp) at a step length of t = 0.1 + + std::cout << "Comparison between exact evaluation and 4th order Taylor estimate:" << std::endl; + std::cout << "u(exact) = " << u.transpose() << std::endl; + std::cout << "u(taylor) = " << utaylor.transpose() << std::endl; +} + +/*------------------------------------------------------------------------------------------------- +=== Output === +--------------------------------------------------------------------------------------------------- +Directional derivatives of f along v = (vx, vp) up to 4th order: +dfdv[0] = +0.000000e+00 +4.852030e+00 +7.690286e+00 +9.704061e+00 +1.126607e+01 +dfdv[1] = +7.000000e+00 +4.193147e+00 +3.431946e+00 +3.136294e+00 +3.009438e+00 +dfdv[2] = -5.000000e+00 -7.500000e-01 -1.111111e-01 +6.250000e-02 +1.200000e-01 +dfdv[3] = +1.100000e+01 +1.000000e+00 +1.851852e-01 +3.125000e-02 -8.000000e-03 +dfdv[4] = -3.400000e+01 -1.625000e+00 -2.222222e-01 -3.906250e-02 -3.200000e-03 + +Comparison between exact evaluation and 4th order Taylor estimate: +u(exact) = +6.767023e-01 +5.267755e+00 +8.032955e+00 +1.001801e+01 +1.156761e+01 +u(taylor) = +6.766917e-01 +5.267755e+00 +8.032955e+00 +1.001801e+01 +1.156761e+01 +-------------------------------------------------------------------------------------------------*/ + +/*------------------------------------------------------------------------------------------------- +=== Note === +--------------------------------------------------------------------------------------------------- +This example would also work if dual was used instead of real. However, real +types are your best option for directional derivatives, as they were optimally +designed for this kind of derivatives. +-------------------------------------------------------------------------------------------------*/ \ No newline at end of file diff --git a/test/forward/14taylor_series_of_scalar_func_along_a_direction.cpp b/test/forward/14taylor_series_of_scalar_func_along_a_direction.cpp new file mode 100644 index 0000000..654fb20 --- /dev/null +++ b/test/forward/14taylor_series_of_scalar_func_along_a_direction.cpp @@ -0,0 +1,42 @@ +// NOTE: 标量函数沿某个方向的泰勒级数 +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The scalar function for which a 4th order directional Taylor series will be computed. +real4th f(const real4th& x, const real4th& y, const real4th& z) +{ + return sin(x * y) * cos(x * z) * exp(z); +} + +int main() +{ + real4th x = 1.0; // the input vector x + real4th y = 2.0; // the input vector y + real4th z = 3.0; // the input vector z + + auto g = taylorseries(f, along(1, 1, 2), at(x, y, z)); // the function g(t) as a 4th order Taylor approximation of f(x + t, y + t, z + 2t) + + double t = 0.1; // the step length used to evaluate g(t), the Taylor approximation of f(x + t, y + t, z + 2t) + + real4th u = f(x + t, y + t, z + 2*t); // the exact value of f(x + t, y + t, z + 2t) + + double utaylor = g(t); // the 4th order Taylor estimate of f(x + t, y + t, z + 2t) + + std::cout << std::fixed; + std::cout << "Comparison between exact evaluation and 4th order Taylor estimate of f(x + t, y + t, z + 2t):" << std::endl; + std::cout << "u(exact) = " << u << std::endl; + std::cout << "u(taylor) = " << utaylor << std::endl; +} + +/*------------------------------------------------------------------------------------------------- +=== Output === +--------------------------------------------------------------------------------------------------- +Comparison between exact evaluation and 4th order Taylor estimate of f(x + t, y + t, z + 2t): +u(exact) = -16.847071 +u(taylor) = -16.793986 +-------------------------------------------------------------------------------------------------*/ \ No newline at end of file diff --git a/test/forward/15taylor_series_of_vector_func_along_a_direction.cpp b/test/forward/15taylor_series_of_vector_func_along_a_direction.cpp new file mode 100644 index 0000000..9ad23c7 --- /dev/null +++ b/test/forward/15taylor_series_of_vector_func_along_a_direction.cpp @@ -0,0 +1,46 @@ +// NOTE: 向量函数沿某个方向的泰勒级数 +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The vector function for which a 4th order directional Taylor series will be computed. +ArrayXreal4th f(const ArrayXreal4th& x) +{ + return x.sin() / x; +} + +int main() +{ + using Eigen::ArrayXd; + + ArrayXreal4th x(5); // the input vector x + x << 1.0, 2.0, 3.0, 4.0, 5.0; + + ArrayXd v(5); // the direction vector v + v << 1.0, 1.0, 1.0, 1.0, 1.0; + + auto g = taylorseries(f, along(v), at(x)); // the function g(t) as a 4th order Taylor approximation of f(x + t·v) + + double t = 0.1; // the step length used to evaluate g(t), the Taylor approximation of f(x + t·v) + + ArrayXreal4th u = f(x + t * v); // the exact value of f(x + t·v) + + ArrayXd utaylor = g(t); // the 4th order Taylor estimate of f(x + t·v) + + std::cout << std::fixed; + std::cout << "Comparison between exact evaluation and 4th order Taylor estimate of f(x + t·v):" << std::endl; + std::cout << "u(exact) = " << u.transpose() << std::endl; + std::cout << "u(taylor) = " << utaylor.transpose() << std::endl; +} + +/*------------------------------------------------------------------------------------------------- +=== Output === +--------------------------------------------------------------------------------------------------- +Comparison between exact evaluation and 4th order Taylor estimate of f(x + t·v): +u(exact) = 0.810189 0.411052 0.013413 -0.199580 -0.181532 +u(taylor) = 0.810189 0.411052 0.013413 -0.199580 -0.181532 +-------------------------------------------------------------------------------------------------*/ diff --git a/test/forward/1single_variable_function.cpp b/test/forward/1single_variable_function.cpp new file mode 100644 index 0000000..7df5897 --- /dev/null +++ b/test/forward/1single_variable_function.cpp @@ -0,0 +1,25 @@ +// NOTE: 单变量函数的导数 + +// C++ includes +#include + +// autodiff include +#include +using namespace autodiff; + +// The single-variable function for which derivatives are needed +dual f(dual x) +{ + return 1 + x + x*x + 1/x + log(x); +} + +int main() +{ + dual x = 2.0; // the input variable x + dual u = f(x); // the output variable u + + double dudx = derivative(f, wrt(x), at(x)); // evaluate the derivative du/dx + + std::cout << "u = " << u << std::endl; // print the evaluated output u + std::cout << "du/dx = " << dudx << std::endl; // print the evaluated derivative du/dx +} \ No newline at end of file diff --git a/test/forward/2single_variable_function_with_custom_scalar.cpp b/test/forward/2single_variable_function_with_custom_scalar.cpp new file mode 100644 index 0000000..5bdc40a --- /dev/null +++ b/test/forward/2single_variable_function_with_custom_scalar.cpp @@ -0,0 +1,32 @@ +// NOTE: 使用自定义标量（复数）的单变量函数的导数 + +// C++ includes +#include +#include +using namespace std; + +// autodiff include +#include +using namespace autodiff; + +// Specialize isArithmetic for complex to make it compatible with dual +namespace autodiff::detail { + template + struct ArithmeticTraits>: ArithmeticTraits {}; +} // autodiff::detail + +using cxdual = Dual, complex>; + +cxdual f(cxdual x){ + return 1+x+x*x+1/x+log(x); +} + +int main(){ + cxdual x = 2.0; // the input variable x + cxdual u = f(x); // the output variable u + cxdual dudx = derivative(f, wrt(x), at(x)); // evaluate the derivative du/dx + + cout << "u = " << u << endl; + cout << "du/dx = " << dudx << endl; + return 0; +} \ No newline at end of file diff --git a/test/forward/3multi_varialbe_function.cpp b/test/forward/3multi_varialbe_function.cpp new file mode 100644 index 0000000..0f359da --- /dev/null +++ b/test/forward/3multi_varialbe_function.cpp @@ -0,0 +1,32 @@ +// NOTE: 多变量函数的导数 + +// C++ includes +#include + +// autodiff include +#include +using namespace autodiff; + +// The multi-variable function for which derivatives are needed +dual f(dual x, dual y, dual z) +{ + return 1 + x + y + z + x*y + y*z + x*z + x*y*z + exp(x/y + y/z); +} + +int main() +{ + dual x = 1.0; + dual y = 2.0; + dual z = 3.0; + + dual u = f(x, y, z); + + double dudx = derivative(f, wrt(x), at(x, y, z)); + double dudy = derivative(f, wrt(y), at(x, y, z)); + double dudz = derivative(f, wrt(z), at(x, y, z)); + + std::cout << "u = " << u << std::endl; // print the evaluated output u = f(x, y, z) + std::cout << "du/dx = " << dudx << std::endl; // print the evaluated derivative du/dx + std::cout << "du/dy = " << dudy << std::endl; // print the evaluated derivative du/dy + std::cout << "du/dz = " << dudz << std::endl; // print the evaluated derivative du/dz +} \ No newline at end of file diff --git a/test/forward/4multi_variable_function_with_params.cpp b/test/forward/4multi_variable_function_with_params.cpp new file mode 100644 index 0000000..dae543a --- /dev/null +++ b/test/forward/4multi_variable_function_with_params.cpp @@ -0,0 +1,59 @@ +// NOTE: 带参数的多变量函数的导数 + +// C++ includes +#include + +// autodiff include +#include +using namespace autodiff; + +// A type defining parameters for a function of interest +struct Params +{ + dual a; + dual b; + dual c; +}; + +// The function that depends on parameters for which derivatives are needed +dual f(dual x, const Params& params) +{ + return params.a * sin(x) + params.b * cos(x) + params.c * sin(x)*cos(x); +} + +int main() +{ + Params params; // initialize the parameter variables + params.a = 1.0; // the parameter a of type dual, not double! + params.b = 2.0; // the parameter b of type dual, not double! + params.c = 3.0; // the parameter c of type dual, not double! + + dual x = 0.5; // the input variable x + + dual u = f(x, params); // the output variable u + + double dudx = derivative(f, wrt(x), at(x, params)); // evaluate the derivative du/dx + double duda = derivative(f, wrt(params.a), at(x, params)); // evaluate the derivative du/da + double dudb = derivative(f, wrt(params.b), at(x, params)); // evaluate the derivative du/db + double dudc = derivative(f, wrt(params.c), at(x, params)); // evaluate the derivative du/dc + + std::cout << "u = " << u << std::endl; // print the evaluated output u + std::cout << "du/dx = " << dudx << std::endl; // print the evaluated derivative du/dx + std::cout << "du/da = " << duda << std::endl; // print the evaluated derivative du/da + std::cout << "du/db = " << dudb << std::endl; // print the evaluated derivative du/db + std::cout << "du/dc = " << dudc << std::endl; // print the evaluated derivative du/dc +} + +/* NOTE: + +This example would also work if real was used instead of dual. Should you +need higher-order cross derivatives, however, e.g.,: + + double d2udxda = derivative(f, wrt(x, params.a), at(x, params)); + +then higher-order dual types are the right choicesince real types are +optimally designed for higher-order directional derivatives. + +如果需要计算更高阶的交叉倒数（例如d2udxda）时，应该使用更高阶的dual类型，而不是real类型 + +*/ diff --git a/test/forward/5multi_varialbe_function_rely_on_analytical_der.cpp b/test/forward/5multi_varialbe_function_rely_on_analytical_der.cpp new file mode 100644 index 0000000..2b16d81 --- /dev/null +++ b/test/forward/5multi_varialbe_function_rely_on_analytical_der.cpp @@ -0,0 +1,95 @@ +// NOTE: 依赖于解析导数的多变量函数的导数 + +// C++ includes +#include + +// autodiff include +#include +using namespace autodiff; + +// Define functions A, Ax, Ay using double; analytical derivatives are available. +double A(double x, double y) { return x*y; } +double Ax(double x, double y) { return x; } +double Ay(double x, double y) { return y; } + +// Define functions B, Bx, By using double; analytical derivatives are available. +double B(double x, double y) { return x + y; } +double Bx(double x, double y) { return 1.0; } +double By(double x, double y) { return 1.0; } + +// Wrap A into Adual function so that it can be used within autodiff-enabled codes. +// A是一个double相关的函数，无法与autodiff混合使用；这里将A包装成Adual +dual Adual(dual const& x, dual const& y) +{ + dual res = A(x.val, y.val); + + if(x.grad != 0.0) + res.grad += x.grad * Ax(x.val, y.val); + + if(y.grad != 0.0) + res.grad += y.grad * Ay(x.val, y.val); + + return res; +} + +// Wrap B into Bdual function so that it can be used within autodiff-enabled codes. +dual Bdual(dual const& x, dual const& y) +{ + dual res = B(x.val, y.val); + + if(x.grad != 0.0) + res.grad += x.grad * Bx(x.val, y.val); + + if(y.grad != 0.0) + res.grad += y.grad * By(x.val, y.val); + + return res; +} + +// Define autodiff-enabled C function that relies on Adual and Bdual +dual C(dual const& x, dual const& y) +{ + const auto A = Adual(x, y); + const auto B = Bdual(x, y); + return A*A + B; +} + +int main() +{ + // --- begin of test --- + double tmp=0.0; + dual t1 = tmp; + dual t2 = 0.0; + std::cout << t1 << " " << t2 << std::endl; + // --- end of test --- + + dual x = 1.0; + dual y = 2.0; + + auto C0 = C(x, y); + + // Compute derivatives of C with respect to x and y using autodiff! + auto Cx = derivative(C, wrt(x), at(x, y)); + auto Cy = derivative(C, wrt(y), at(x, y)); + + // Compute expected analytical derivatives of C with respect to x and y + auto x0 = x.val; + auto y0 = y.val; + auto expectedCx = 2.0*A(x0, y0)*Ax(x0, y0) + Bx(x0, y0); + auto expectedCy = 2.0*A(x0, y0)*Ay(x0, y0) + By(x0, y0); + + std::cout << "C0 = " << C0 << "\n"; + + std::cout << "Cx(computed) = " << Cx << "\n"; + std::cout << "Cx(expected) = " << expectedCx << "\n"; + + std::cout << "Cy(computed) = " << Cy << "\n"; + std::cout << "Cy(expected) = " << expectedCy << "\n"; +} + +// Output: +// C0 = 7 +// Cx(computed) = 5 +// Cx(expected) = 5 +// Cy(computed) = 9 +// Cy(expected) = 9 \ No newline at end of file diff --git a/test/forward/6gradient_vector_of_scalar_func.cpp b/test/forward/6gradient_vector_of_scalar_func.cpp new file mode 100644 index 0000000..c02e15c --- /dev/null +++ b/test/forward/6gradient_vector_of_scalar_func.cpp @@ -0,0 +1,29 @@ +// NOTE: 标量函数的梯度向量 +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The scalar function for which the gradient is needed +real f(const ArrayXreal& x) +{ + return (x * x.exp()).sum(); // sum([xi * exp(xi) for i = 1:5]) +} + +int main() +{ + using Eigen::VectorXd; + + ArrayXreal x(5); // the input array x with 5 variables + x << 1, 2, 3, 4, 5; // x = [1, 2, 3, 4, 5] + + real u; // the output scalar u = f(x) evaluated together with gradient below + + VectorXd g = gradient(f, wrt(x), at(x), u); // evaluate the function value u and its gradient vector g = du/dx + + std::cout << "u = " << u << std::endl; // print the evaluated output u + std::cout << "g = \n" << g << std::endl; // print the evaluated gradient vector g = du/dx +} \ No newline at end of file diff --git a/test/forward/7gradient_vector_of_scalar_func_with_params.cpp b/test/forward/7gradient_vector_of_scalar_func_with_params.cpp new file mode 100644 index 0000000..a33ee86 --- /dev/null +++ b/test/forward/7gradient_vector_of_scalar_func_with_params.cpp @@ -0,0 +1,51 @@ +// NOTE: 带参数的标量函数的梯度向量 +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The scalar function for which the gradient is needed +real f(const ArrayXreal& x, const ArrayXreal& p, const real& q) +{ + return (x * x).sum() * p.sum() * exp(q); // sum([xi * xi for i = 1:5]) * sum(p) * exp(q) +} + +int main() +{ + using Eigen::VectorXd; + + ArrayXreal x(5); // the input vector x with 5 variables + x << 1, 2, 3, 4, 5; // x = [1, 2, 3, 4, 5] + + ArrayXreal p(3); // the input parameter vector p with 3 variables + p << 1, 2, 3; // p = [1, 2, 3] + + real q = -2; // the input parameter q as a single variable + + real u; // the output scalar u = f(x, p, q) evaluated together with gradient below + + VectorXd gx = gradient(f, wrt(x), at(x, p, q), u); // evaluate the function value u and its gradient vector gx = du/dx + VectorXd gp = gradient(f, wrt(p), at(x, p, q), u); // evaluate the function value u and its gradient vector gp = du/dp + VectorXd gq = gradient(f, wrt(q), at(x, p, q), u); // evaluate the function value u and its gradient vector gq = du/dq + VectorXd gqpx = gradient(f, wrt(q, p, x), at(x, p, q), u); // evaluate the function value u and its gradient vector gqpx = [du/dq, du/dp, du/dx] + + std::cout << "u = " << u << std::endl; // print the evaluated output u + std::cout << "gx = \n" << gx << std::endl; // print the evaluated gradient vector gx = du/dx + std::cout << "gp = \n" << gp << std::endl; // print the evaluated gradient vector gp = du/dp + std::cout << "gq = \n" << gq << std::endl; // print the evaluated gradient vector gq = du/dq + std::cout << "gqpx = \n" << gqpx << std::endl; // print the evaluated gradient vector gqpx = [du/dq, du/dp, du/dx] +} + +/*------------------------------------------------------------------------------------------------- +=== Note === +--------------------------------------------------------------------------------------------------- +This example would also work if dual was used instead. However, if gradient, +Jacobian, and directional derivatives are all you need, then real types are your +best option. You want to use dual types when evaluating higher-order cross +derivatives, which is not supported for real types. +如果只需要计算梯度、雅可比矩阵和方向导数，使用real类型 +是最佳的。如果需要计算更高阶的交叉倒数，则应该使用dual类型。 +-------------------------------------------------------------------------------------------------*/ diff --git a/test/forward/8jacobian_matrix_of_vector_func.cpp b/test/forward/8jacobian_matrix_of_vector_func.cpp new file mode 100644 index 0000000..d3f5451 --- /dev/null +++ b/test/forward/8jacobian_matrix_of_vector_func.cpp @@ -0,0 +1,32 @@ +// NOTE: 向量函数的雅可比矩阵 + +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The vector function for which the Jacobian is needed +VectorXreal f(const VectorXreal& x) +{ + return x * x.sum(); +} + +int main() +{ + using Eigen::MatrixXd; + + VectorXreal x(5); // the input vector x with 5 variables + x << 1, 2, 3, 4, 5; // x = [1, 2, 3, 4, 5] + + VectorXreal F; // the output vector F = f(x) evaluated together with Jacobian matrix below + MatrixXd J1 = jacobian(f, wrt(x), at(x), F); // evaluate the output vector F and the Jacobian matrix dF/dx + + std::cout << "F = \n" << F << std::endl; // print the evaluated output vector F + std::cout << "J1 = \n" << J1 << std::endl; // print the evaluated Jacobian matrix dF/dx + + MatrixXd J2 = jacobian(f, wrt(x), at(x)); + std::cout << "J2 = \n" << J2 << std::endl; // print the evaluated Jacobian matrix dF/dx +} diff --git a/test/forward/9jacobian_matrix_of_vector_func_with_params.cpp b/test/forward/9jacobian_matrix_of_vector_func_with_params.cpp new file mode 100644 index 0000000..6a53f39 --- /dev/null +++ b/test/forward/9jacobian_matrix_of_vector_func_with_params.cpp @@ -0,0 +1,40 @@ +// NOTE: 带参数的向量函数的雅可比矩阵 +// C++ includes +#include + +// autodiff include +#include +#include +using namespace autodiff; + +// The vector function with parameters for which the Jacobian is needed +VectorXreal f(const VectorXreal& x, const VectorXreal& p, const real& q) +{ + return x * p.sum() * exp(q); +} + +int main() +{ + using Eigen::MatrixXd; + + VectorXreal x(5); // the input vector x with 5 variables + x << 1, 2, 3, 4, 5; // x = [1, 2, 3, 4, 5] + + VectorXreal p(3); // the input parameter vector p with 3 variables + p << 1, 2, 3; // p = [1, 2, 3] + + real q = -2; // the input parameter q as a single variable + + VectorXreal F; // the output vector F = f(x, p, q) evaluated together with Jacobian below + + MatrixXd Jx = jacobian(f, wrt(x), at(x, p, q), F); // evaluate the function and the Jacobian matrix Jx = dF/dx + MatrixXd Jp = jacobian(f, wrt(p), at(x, p, q), F); // evaluate the function and the Jacobian matrix Jp = dF/dp + MatrixXd Jq = jacobian(f, wrt(q), at(x, p, q), F); // evaluate the function and the Jacobian matrix Jq = dF/dq + MatrixXd Jqpx = jacobian(f, wrt(q, p, x), at(x, p, q), F); // evaluate the function and the Jacobian matrix Jqpx = [dF/dq, dF/dp, dF/dx] + + std::cout << "F = \n" << F << std::endl; // print the evaluated output vector F + std::cout << "Jx = \n" << Jx << std::endl; // print the evaluated Jacobian matrix dF/dx + std::cout << "Jp = \n" << Jp << std::endl; // print the evaluated Jacobian matrix dF/dp + std::cout << "Jq = \n" << Jq << std::endl; // print the evaluated Jacobian matrix dF/dq + std::cout << "Jqpx = \n" << Jqpx << std::endl; // print the evaluated Jacobian matrix [dF/dq, dF/dp, dF/dx] +} diff --git a/test/forward/CMakeLists.txt b/test/forward/CMakeLists.txt new file mode 100644 index 0000000..d4fe459 --- /dev/null +++ b/test/forward/CMakeLists.txt @@ -0,0 +1,9 @@ + +file(GLOB CPPFILES "${CMAKE_CURRENT_SOURCE_DIR}/*.cpp") +foreach(item ${CPPFILES}) + get_filename_component(itemname ${item} NAME_WE) + add_executable(${itemname} ${item}) + message(STATUS "add_executable(${itemname} ${item})") + # target_link_libraries(${itemname} Eigen3::Eigen) +endforeach() +