Code Monkey home page Code Monkey logo

clad's Introduction

Clad

Conda-Forge Anaconda-Server Badge Anaconda-Server Badge Anaconda-Server Badge

Linux & Osx Status Coverity Scan Build Status codecov Binder

Clad is a source-transformation automatic differentiation (AD) library for C++,
implemented as a plugin for the Clang compiler.

About Clad

Clad enables automatic differentiation (AD) for C++. It is based on LLVM compiler infrastructure and is a plugin for Clang compiler. Clad is based on source code transformation. Given C++ source code of a mathematical function, it can automatically generate C++ code for computing derivatives of the function. It supports both forward-mode and reverse-mode AD.Clad has extensive coverage of modern C++ features and a robust fallback and recovery system in place.

How to use Clad

Clad provides five API functions:

API functions are used to label an existing function for differentiation. Both functions return a functor object containing the generated derivative which can be called via .execute method, which forwards provided arguments to the generated derivative function.

For a guide on compiling your clad-based programs, look here.

Forward mode - clad::differentiate

For a function f of several inputs and single (scalar) output, forward mode AD can be used to compute (or, in case of Clad, create a function) a directional derivative of f with respect to a single specified input variable. Derivative function created by the forward-mode AD is guaranteed to have at most a constant factor (around 2-3) more arithmetical operations compared to the original function.

clad::differentiate(f, ARGS) takes 2 arguments:

  1. f is a pointer to a function or a method to be differentiated
  2. ARGS is either:
  • a single numerical literal indicating an index of independent variable (e.g. 0 for x, 1 for y)
  • a string literal with the name of independent variable (as stated in the definition of f, e.g. "x" or "y"), and if the variable is an array the index needs to be specified, e.g. "arr[1]"

Generated derivative function has the same signature as the original function f, however its return value is the value of the derivative. Example:

#include "clad/Differentiator/Differentiator.h"
#include <iostream>

double f(double x, double y) { return x * y; }

int main() {
  // Call clad to generate the derivative of f wrt x.
  auto f_dx = clad::differentiate(f, "x");
  // Execute the generated derivative function.
  std::cout << f_dx.execute(/*x=*/3, /*y=*/4) << std::endl;
  // Dump the generated derivative code to standard output.
  f_dx.dump();
}

Reverse mode - clad::gradient

Reverse-mode AD allows computing the gradient of f using at most a constant factor (around 4) more arithmetical operations compared to the original function. While its constant factor and memory overhead is higher than that of the forward-mode, it is independent of the number of inputs. E.g. for a function having N inputs and consisting of T arithmetical operations, computing its gradient takes a single execution of the reverse-mode AD and around 4*T operations, while it would take N executions of the forward-mode, this requiring up to N*3*T operations.

clad::gradient(f, /*optional*/ ARGS) takes 1 or 2 arguments:

  1. f is a pointer to a function or a method to be differentiated
  2. ARGS is either:
  • not provided, then f is differentiated w.r.t. its every argument
  • a string literal with comma-separated names of independent variables (e.g. "x" or "y" or "x, y" or "y, x")

Since a vector of derivatives must be returned from a function generated by the reverse mode, its signature is slightly different. The generated function has void return type and same input arguments. The function has additional n arguments (where n refers to the number of arguments whose gradient was requested) of type T*, where T is the type of the corresponding original variable. Each of these variables stores the derivative of the elements as they appear in the orignal function signature. The caller is responsible for allocating and zeroing-out the gradient storage. Example:

auto f_grad = clad::gradient(f);
double dx = 0, dy = 0;
// After this call, dx and dy will store the derivatives of x and y respectively.
f_grad.execute(x, y, &dx, &dy);
std::cout << "dx: " << dx << ' ' << "dy: " << dy << std::endl;

// Same effect as before.
auto f_dx_dy = clad::gradient(f, "x, y"); 
auto f_dy_dx = clad::gradient(f, "y, x");

// The same effect can be achieved by using an array instead of individual variables.
double result2[2] = {};
f_dy_dx.execute(x, y, /*dx=*/&result2[0], /*dy=*/&result2[1]);
// note that the derivatives are mapped to the "result" indices in the same order as they were specified in the argument:
std::cout << "dy: " << result2[0] << ' ' << "dx: " << result2[1] << std::endl;

Hessian mode - clad::hessian

Clad can produce the hessian matrix of a function using its forward and reverse mode capabilities. Its interface is similar to reverse mode but differs when arrays are involved. It returns the matrix as a flattened vector in row major format.

clad::hessian(f, /*optional*/ ARGS) takes 1 or 2 arguments:

  1. f is a pointer to a function or a method to be differentiated
  2. ARGS is either:
    • not provided, then f is differentiated w.r.t. its every argument except in the case of arrays where it needs to be provided
    • a string literal with comma-separated names of independent variables (e.g. "x" or "y" or "x, y" or "y, x" or in case of arrays "x[0:2]")

The generated function has void return type and same input arguments. The function has an additional argument of type T*, where T is the return type of f. This variable stores the hessian matrix. The caller is responsible for allocating and zeroing-out the hessian storage. Example:

#include "clad/Differentiator/Differentiator.h"
#include <iostream>

double f(double x, double y) { return x * y; }
double g(double x, double y[2]) { return x * y[0] * y[1]; }

int main() {
    // Since we are differentiating variables that are not arrays the interface
    // is same as in reverse mode
    auto f_hess = clad::hessian(f);
    // The size of the resultant matrix should be the square of the 
    // number of independent variables
    double mat_f[4] = {0};
    
    // Execute the hessian function
    f_hess.execute(/*x=*/3, /*y=*/4, mat_f);
    std::cout << "[" << mat_f[0] << ", " << mat_f[1] << "\n  " 
                     << mat_f[2] << ", " << mat_f[3] << "]";
    
    // When arrays are involved the array indexes that are to be differentiated needs to be specified
    // even if the whole array is being differentiated
    auto g_hess = clad::hessian(g, "x, y[0:1]");
    // The rest of the steps are the same.
}

Jacobian mode - clad::jacobian

Clad can produce the jacobian of a function using its reverse mode. It returns the jacobian matrix as a flattened vector in row major format.

clad::jacobian(f, /*optional*/ ARGS) takes 1 or 2 arguments:

  1. f is a pointer to a function or a method to be differentiated
  2. ARGS is either:
    • not provided, then f is differentiated w.r.t. its every argument
    • a string literal with comma-separated names of independent variables (e.g. "x" or "y" or "x, y" or "y, x")

The generated function has void return type and same input arguments. The function has an additional argument of type T *, where T is the pointee type of the output (the last variable) of f. This variable stores the jacobian matrix. The caller is responsible for allocating and zeroing-out the jacobian storage. Example:

#include "clad/Differentiator/Differentiator.h"
#include <iostream>

void h(double a, double b, double output[]) {
    output[0] = a * a * a;
    output[1] = a * a * a + b * b * b;
    output[2] = 2 * (a + b);
}

int main() {
    // This sets all the input variables (i.e a and b) as independent variables 
    auto h_jac = clad::jacobian(h);
    
    // The jacobian matrix size should be the number of 
    // independent variables * the number of outputs of the original function
    // In this case it is 2 * 3 = 6
    double jac[6] = {0};
    double output[3] = {0};
    h_jac.execute(/*a=*/3, /*b=*/4, output, jac);
    
    std::cout << jac[0] << " " << jac[1] << std::endl
              << jac[2] << " " << jac[3] << std::endl
              << jac[4] << " " << jac[5] << std::endl;
}

Floating-point error estimation - clad::estimate_error

Clad is capable of annotating a given function with floating point error estimation code using the reverse mode of AD. An interface similar to clad::gradient(f) is provided as follows:

clad::estimate_error(f) takes 1 argument:

  1. f is a pointer to the function or method to be annotated with floating point error estimation code.

The function signature of the generated code is the same as from clad::gradient(f) with the exception that it has an extra argument at the end of type double&, which returns the total floating point error in the function by reference. For a user function double f(double, double) example usage is described below:

// Generate the floating point error estimation code for 'f'.
auto df = clad::estimate_error(f);
// Print the generated code to standard output.
df.dump();
// Declare the necessary variables.
double x, y, d_x, d_y, final_error = 0;
// Finally call execute on the generated code.
df.execute(x, y, &d_x, &d_y, final_error);
// After this, 'final_error' contains the floating point error in function 'f'.

The above example generates the the error code using an in-built taylor approximation model. However, clad is capable of using any user defined custom model, for information on how to use you own custom model, please visit this demo.

More detail on the APIs can be found under clad's user documentation.

Compiling and executing your code with clad

Using Jupyter Notebooks

xeus-cling provides a Jupyter kernel for C++ with the help of the C++ interpreter Cling and the native implementation of the Jupyter protocol xeus. Within the xeus-cling framework, Clad can enable automatic differentiation (AD) such that users can automatically generate C++ code for their computation of derivatives of their functions.

To set up your environment, use:

mamba create -n xeus-clad -c conda-forge clad xeus-cling jupyterlab

conda activate xeus-clad

Next, running jupyter notebook will show 3 new kernels for C++ 11/14/17 with Clad attached.

Try out a Clad tutorial interactively in your browser through binder:

Binder

Using as a plugin for Clang

Since Clad is a Clang plugin, it must be properly attached when Clang compiler is invoked. First, the plugin must be built to get libclad.so (or .dylib).

To compile SourceFile.cpp with Clad enabled, use the following commands:

  • Clang++: clang++ -std=c++11 -I /full/path/to/include/ -fplugin=/full/path/to/lib/clad.so Sourcefile.cpp
  • Clang: clang -x c++ -std=c++11 -I /full/path/to/include/ -fplugin=/full/path/to/lib/clad.so SourceFile.cpp -lstdc++ -lm

Clad also provides certain flags to save and print the generated derivative code:

  • To save the Clad generated derivative code to Derivatives.cpp: -Xclang -plugin-arg-clad -Xclang -fgenerate-source-file
  • To print the Clad generated derivative: -Xclang -plugin-arg-clad -Xclang -fdump-derived-fn

How to install

At the moment, LLVM/Clang 8.0.x - 18.1.x are supported.

Conda Installation

Clad is available using conda:

conda install -c conda-forge clad

If you have already added conda-forge as a channel, the -c conda-forge is unnecessary. Adding the channel is recommended because it ensures that all of your packages use compatible versions:

conda config --add channels conda-forge
conda update --all

Building from source (example was tested on Ubuntu 20.04 LTS)

#sudo apt install clang-11 libclang-11-dev llvm-11-tools llvm-11-dev
sudo bash -c "$(wget -O - https://apt.llvm.org/llvm.sh)"
sudo -H pip install lit
git clone https://github.com/vgvassilev/clad.git clad
mkdir build_dir inst; cd build_dir
cmake ../clad -DClang_DIR=/usr/lib/llvm-11 -DLLVM_DIR=/usr/lib/llvm-11 -DCMAKE_INSTALL_PREFIX=../inst -DLLVM_EXTERNAL_LIT="$(which lit)"
make && make install

NOTE: On some Linux distributions (e.g. Arch Linux), the LLVM and Clang libraries are installed at /usr/lib/cmake/llvm and /usr/lib/cmake/clang. If compilation fails with the above provided command, ensure that you are using the correct path to the libraries.

Building from source (example was tested on macOS Big Sur 11.6)

brew install llvm@12
brew install python
python -m pip install lit
git clone https://github.com/vgvassilev/clad.git clad
mkdir build_dir inst; cd build_dir
cmake ../clad -DLLVM_DIR=/opt/homebrew/opt/llvm@12/lib/cmake/llvm -DClang_DIR=/opt/homebrew/opt/llvm@12/lib/cmake/clang -DCMAKE_INSTALL_PREFIX=../inst  -DLLVM_EXTERNAL_LIT="`which lit`"
make && make install
make check-clad

Developer Environment - Build LLVM, Clang and Clad from source:

pip3 install lit

Clone the LLVM project and checkout the required LLVM version (Currently supported versions 8.x - 18.x)

git clone https://github.com/llvm/llvm-project.git
cd llvm-project
git checkout llvmorg-16.0.0

Build Clang:

mkdir build && cd build
cmake -DLLVM_ENABLE_PROJECTS="clang" -DCMAKE_BUILD_TYPE="DEBUG" -DLLVM_TARGETS_TO_BUILD=host -DLLVM_INSTALL_UTILS=ON ../llvm
cmake --build . --target clang --parallel $(nproc --all)
make -j8 check-clang # this installs llvm-config required by lit
cd ../..

Clone and build Clad:

git clone https://github.com/vgvassilev/clad.git
cd clad
mkdir build && cd build
cmake -DLLVM_DIR=PATH/TO/llvm-project/build -DCMAKE_BUILD_TYPE=DEBUG -DLLVM_EXTERNAL_LIT="$(which lit)" ../
make -j8 clad

Run the Clad tests:

make -j8 check-clad

Further reading

What can be differentiated

Clad is based on compile-time analysis and transformation of C++ abstract syntax tree (Clang AST). This means that Clad must be able to see the body of a function to differentiate it (e.g. if a function is defined in an external library there is no way for Clad to get its AST).

Note: Clad currently differentiates types such as int/char/boolean as any real type (such as float, double, etc.) would be differentiated. Users should keep in mind that Clad does not warn against lossy casts, which on differentiation may result in incorrect derivatives.

Note: If for any reason clad is unable to algorithmically differentiate a function, it automatically switches to numerically differentiating the same. To disable this behavior, please compile your programs with the -DCLAD_NO_NUM_DIFF flag. The numerical differentiation functionality can also be used standalone over a wide range of function signatures with minimal user intervention. This presentation provides more information on what can be numerically differentiated. For a comprehensive demo on using custom user defined types with numerical differentiation, you can check out this demo.

Specifying custom derivatives

Sometimes Clad may be unable to differentiate your function (e.g. if its definition is in a library and source code is not available). Alternatively, an efficient/more numerically stable expression for derivatives may be know. In such cases, it is useful to be able to specify a custom derivatives for your function.

Clad supports that functionality by allowing to specify your own derivatives in namespace custom_derivatives. For a function named FNAME you can specify:

  • a custom derivative w.r.t I-th argument by defining a function FNAME_dargI inside namespace custom_derivatives
  • a custom gradient w.r.t every argument by defining a function FNAME_grad inside namespace custom_derivatives

When Clad will encounter a function FNAME, it will first do a lookup inside the custom_derivatives namespace to try to find a suitable custom function, and only if none is found will proceed to automatically derive it.

Example:

  • Suppose that you have a function my_pow(x, y) which computes x to the power of y. However, Clad is not able to differentiate my_pow's body (e.g. it calls an external library or uses some non-differentiable approximation):
double my_pow(double x, double y) { // something non-differentiable here... }

However, you know analytical formulas of its derivatives, and you can easily specify custom derivatives:

namespace custom_derivatives {
  double my_pow_darg0(double x, double y) { return y * my_pow(x, y - 1); }
  double my_pow_darg1(dobule x, double y) { return my_pow(x, y) * std::log(x); }
}

You can also specify a custom gradient:

namespace custom_derivatives {
  void my_pow_grad(double x, double y, array_ref<double> _d_x, array_ref<double> _d_y) {
     double t = my_pow(x, y - 1);
     *_d_x = y * t;
     *_d_y = x * t * std::log(x);
   }
}

Whenever Clad will encounter my_pow inside differentiated function, it will find and use provided custom functions instead of attempting to differentiate it.

Note: Clad provides custom derivatives for some mathematical functions from <cmath> inside clad/Differentiator/BuiltinDerivatives.h.

Details on custom derivatives, other supported C++ syntax (already supported or in-progress) and further resources can be found over at clad's user documentation.

Citing Clad

% Peer-Reviewed Publication
%
% 16th International workshop on Advanced Computing and Analysis Techniques
% in physics research (ACAT), 1-5 September, 2014, Prague, The Czech Republic
%
@inproceedings{Vassilev_Clad,
  author = {Vassilev,V. and Vassilev,M. and Penev,A. and Moneta,L. and Ilieva,V.},
  title = {{Clad -- Automatic Differentiation Using Clang and LLVM}},
  journal = {Journal of Physics: Conference Series},
  year = 2015,
  month = {may},
  volume = {608},
  number = {1},
  pages = {012055},
  doi = {10.1088/1742-6596/608/1/012055},
  url = {https://iopscience.iop.org/article/10.1088/1742-6596/608/1/012055/pdf},
  publisher = {{IOP} Publishing}
}

Founders

Founder of the project is Vassil Vassilev as part of his research interests and vision. He holds the exclusive copyright and other related rights, described in Copyright.txt.

License

clad is an open source project, licensed by GNU LESSER GENERAL PUBLIC LICENSE (see License.txt). If there is module with different that LGPL license it will be explicitly stated in the License.txt in the module's source code folder.

Please see License.txt for further information.

How to Contribute

We are always looking for improvements to the tool, as such open-source developers are greatly appreciated! If you are interested in getting started with contributing to clad, make sure you checkout our contribution guide.

clad's People

Contributors

alexander-penev avatar arora-vidushi avatar blide avatar davidlange6 avatar dependabot[bot] avatar efremale avatar gojakuch avatar grimmmyshini avatar hahnjo avatar ioanaif avatar jacklqiu avatar kchristin22 avatar krishna-13-cyber avatar maximusron avatar mcbarton avatar mfoco avatar mihailmihov avatar mvassilev avatar nirhar avatar oshadura avatar parth-07 avatar pcanal avatar petrozarytskyi avatar phrygiangates avatar reikdas avatar ris-bali avatar rohanjulka19 avatar sudo-panda avatar vaithak avatar vgvassilev avatar

Stargazers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

Watchers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

clad's Issues

Clad picks up system LLVM

Even if -DCLAD_PATH_TO_LLVM_BUILD is passed, the system llvm is discovered first.

 /cvmfs/sft.cern.ch/lcg/releases/CMake/3.11.1-773ff/x86_64-centos7-gcc62-opt/bin/cmake  -DCMAKE_BUILD_TYPE=Release -DCMAKE_C_COMPILER=/cvmfs/projects.cern.ch/intelsw/psxe/linux/x86_64/2018/compilers_and_libraries_2018.2.199/linux/bin/intel64/icc -DCMAKE_C_FLAGS= -restrict -wd1572 -m64 -wd279 -wd873 -wd2536 -wd597 -wd1098 -wd1292 -wd1478 -wd3373 -pthread -fPIC -Werror=date-time -w -ffunction-sections -fdata-sections -DCMAKE_CXX_COMPILER=/cvmfs/projects.cern.ch/intelsw/psxe/linux/x86_64/2018/compilers_and_libraries_2018.2.199/linux/bin/intel64/icc -DCMAKE_CXX_FLAGS= -wd1476 -wd1572 -m64 -wd279 -wd873 -wd2536 -wd597 -wd1098 -wd1292 -wd1478 -wd3373 -pthread -std=c++11  -fPIC -fvisibility-inlines-hidden -Werror=date-time -std=c++11 -w -ffunction-sections -fdata-sections -fno-common -Woverloaded-virtual -fno-strict-aliasing -Wno-nested-anon-types -Wno-covered-switch-default -Wno-unused-local-typedef -fno-rtti -DCMAKE_INSTALL_PREFIX=/data/sftnight/workspace/root-benchmark/BUILDTYPE/Release/COMPILER/icc18/LABEL/performance-sandy-cc7/build/etc/cling/plugins/ -DCLAD_PATH_TO_LLVM_BUILD=/data/sftnight/workspace/root-benchmark/BUILDTYPE/Release/COMPILER/icc18/LABEL/performance-sandy-cc7/build/interpreter/llvm/src -G"Unix Makefiles" /data/sftnight/workspace/root-benchmark/BUILDTYPE/Release/COMPILER/icc18/LABEL/performance-sandy-cc7/build/interpreter/llvm/src/tools/cling/tools/plugins/clad/clad-prefix/src/clad
-- Could NOT find Subversion (missing: Subversion_SVN_EXECUTABLE) 
CMake Error at /opt/llvm-5.0.1/lib64/cmake/llvm/LLVMExports.cmake:960 (message):
  The imported target "LLVMDemangle" references the file

     "/opt/llvm-5.0.1/lib64/libLLVMDemangle.a"

  but this file does not exist.  Possible reasons include:

  * The file was deleted, renamed, or moved to another location.

  * An install or uninstall procedure did not complete successfully.

  * The installation package was faulty and contained

     "/opt/llvm-5.0.1/lib64/cmake/llvm/LLVMExports.cmake"

  but not all the files it references.

Call Stack (most recent call first):
  /opt/llvm-5.0.1/lib64/cmake/llvm/LLVMConfig.cmake:188 (include)
  CMakeLists.txt:8 (find_package)


-- Configuring incomplete, errors occurred!

Clad picks up wrong clang version

When embedded in ROOT we want to build clad with external LLVM but internal clang. Currently there is not way of detecting this and we pick up the wrong header files which manifest into a very hard-to-find bug.

Allow to turn on/off clad derivative scanning

Such mechanism would help clients to selectively enable the derivative computation. Eg:

#pragma clad on
// synthesized code which we know will have a call to clad::*
#pragma clad off

Add support for vector of points

Support double f(double*, unsigned numElements); The user should be able to specify the independent variable from a given array of variables.

Derivative name convention

Hi,
The suffix of the derivatives should be darg0, darg1, ... dargN, because we can have scenarios where the users rely on using builtin derivatives. For example:
Builtins.h

float builtin_f(float x, double y); // in builtin derivatives becomes
float builtin_fdx(float x, double y) {}

And the user file:

float builtin_f(float a, double b);
float f(float t, float u) {
  // Here clad won't be able to find builtin_fdx, because there is not enough info in the call
  return t * u * builtin_f(t, u); 
}

Vassil

Repeated expressions in the derived function result in repeated statements in the derivative.

Example:

double f(double x, double y) {
  double t = (x - y) * (x - y);
  return t;
}

->

double f_darg0(double x, double y) {
  double _t0 = (x - y);
  double _t1 = (x - y);
  double _d_t = (1. - 0.) * _t1 + _t0 * (1. - 0.);
  double t = _t0 * _t1;
  return _d_t;
}

Variables _t0 and _t1 are redundant. Since the original function uses (x - y) on several places, we get several temporaries for it.

We do not check the original function for repeated expressions, it would add quite a huge computational overhead. Moreover, since it was acceptable for the original function to use the same expression repeatedly, it is presumably acceptable to do so in the derivative.

Emit hard error when differentiation visits unsupported statement

Currently, if we Visit unsupported statement (e.g. try-catch, range for loop), it is simply cloned without changes and a warning is emitted.

Should we instead emit a hard error and not produce the derivative?

This is somewhat harder to achieve now.
How should the error be promoted in the visitor?
If we return null StmtDiff whenever we meet an unsupported statement, we would have to add a lot of checks for null.
We cannot throw an exception sice clang disables them.
We can set some bool flag and keep visiting function and check for the error flag in the end, but it is inefficient and error-prone.

Assert message and error must be more descriptive

Assertion failed: (lhs_derived->getType() == rhs_derived->getType() && "Must be the same types."), function VisitBinaryOperator, file /Users/alexanderpenev/clad/src/tools/clad/lib/Differentiator/DerivativeBuilder.cpp, line 453.
This message does not help the user. Better text should contains information on function, which can not differentiate, if possible row and column of the user program. Best messages or warnings in style clang.

Implement simple constant folder

The chained rule of differentiation produces a lot of trivial constants to fold. For example:
return 1 + (0) + (1) + (0) + (1).

Methods and Virtual methods do not differentiated

If we have

class A {
public:
  int f(int x) {
    return x;
  }
  virtual float vm(float x, float y) {
    return x + y;
  }
}

class B : public A {
public:
  float vm(float x, float y) override {
    return x*x + y*y;
  }
}

int main() {
  A a;
  B b;
  clad::differentiate(&A::f, 0);
  clad::differentiate(&A::vm, 0);
  printf("Result is = %f\n", a.vm_dx(2,3)); // Result is = 2.F
  clad::differentiate(&B::vm, 0);
  printf("Result is = %f\n", a.vm_dx(2,3)); // Result is = 4.F
  return 0;
}

clad throw error and assert message

error: no member named 'vm_dx' in 'A'
  printf("Result is = %f\n", a.vm_dx(2,3));
                             ~ ^

...

Assertion failed: (Access != AS_none && "Access specifier is AS_none inside a record decl"), function AccessDeclContextSanity, file /Users/alexanderpenev/clad/src/tools/clang/lib/AST/DeclBase.cpp, line 701.

Avoid storing simple expressions like -x in temporaries

Example:

double f(double x) {
  double t = 2 * (-x);
  return t;
}

->

double f_darg0(double x) {
    double _t0 = -x;
    double _d_t = 0 * _t0 + 2 * -1.;
    double t = 2 * _t0;
    return _d_t;
}

In this case, double _t0 = -x seems useless so it makes sense to add such simple unary operators to a list of cases for which temporaries are not generated (now it is literals and references). On the other hand, results of unary operators like ++x must be stored to avoid repeated increment.

Built-in derivatives do not work

There are problem with internal function like sqrt. For example

float func(float x, float y) {
  return sqrt(x * x + y * y) - y;
}

after clad::differentiate(func, 1) is

float func_derived_x(float x, float y) {
    return sqrt(x * x + y * y) - (0.F);
}

but it must something like that:

float func_derived_x(float x, float y) {
    return ((1.F * x + x * 1.F) + ((0.F * y + y * 0.F))) * (1.F / (2.F * sqrt(x * x + y * y)) - (0.F));
}

Clad generate "Must be the same types" assert when one of params is struct

struct Vec {
  float x, y, z;
}

float sphere_implicit_func(float x, float y, float z, Vec &p, float r) {
  return (x-p.x)*(x-p.x) + (y-p.y)*(y-p.y) + (z-p.z)*(z-p.z) - r*r;
}

Fires "Assertion failed: (lhs_derived->getType() == rhs_derived->getType() && "Must be the same types."), function VisitBinaryOperator" during compilation,

but

float sphere_implicit_func(float x, float y, float z, float px, float py, float pz, float r) {
  return (x-px)*(x-px) + (y-py)*(y-py) + (z-pz)*(z-pz) - r*r;
}

is OK.

Would it be possible to implement forward mode by simply replacing all numerical types with dual-number types?

The idea is to replace each type RealT with the type dual<RealT, RealT>, where dual has overloaded operators such that the first element is the same as the original value and the second element is the value of its derivative. This is easy to achieve by employing basic differentiation rules for operators.

On the first sight, it appears to be a much easier way of supporting more general C++ constructs as they need not be treated separately anymore.

Similar idea is in "Automatic Differentiation in 10 minutes with Julia": https://www.youtube.com/watch?v=vAp6nUMrKYg&t=683s

Structure fields do not differentiated correctly

Differentiate

struct Vec { float x,y,z; };
float f_const_args_func_4(float x, float y, const Vec v) {
  return x * x + y * y - v.x;
}

return

float f_const_args_func_4_dx(const float x, const float y, const Vec v) {
  return (1.F * x + x * 1.F) + ((0.F * y + y * 0.F)) - (v.x);
}

but we expect

float f_const_args_func_4_dx(const float x, const float y, const Vec v) {
  return (1.F * x + x * 1.F) + ((0.F * y + y * 0.F)) - (0.F);
}

Recursive functions are not differentiated correctly

A call to itself inside a function body is not treated correctly.

For example, if we differentiate

double r(double x) {
  if (x > 0)
    return r(x-1);
  else
    return x;

the call r(x-1) is treated as a call to some other function and we get the warning: function 'r' was not differentiated because it is not declared in namespace 'custom_derivatives' attempted differention of function 'r', which does not have a definition.

The following derivative is generated:

double r_darg0(double x) {
  if (x > 0)
    return r(x-1);
  else
    return 1.0;

but the correct derivative would be

double r_darg0(double x) {
  if (x > 0)
    return r_darg0(x-1) * 1.0;
  else
    return 1.0;

Add an option to specify independent variables for gradients

Now we can only compute gradients w.r.t. to every argument of a function.

For efficiency reasons, it seems to be a good idea to be able to specify indices required of independent variables, e.g.:

clad::gradient(f, 0, 2, 4)

or

clad::gradient(f, {0, 2, 4})

or even

clad::gradient(f, std::index_sequence<0, 2, 4>{})

Add clearer indication for results of differentiation of functions for which no derivative in custom_derivatives exist

Currently, if differentiation process meets a function which has no defined derivative in custom_derivatives namespace, its derivative is set to be 0.

Example:

double f(double x) {
  printf("%d", x);
  return x;
}

->

double f_darg0(double x) {
  0;
  printf("%d", x);
  return 1;
}

Since there exist no custom_derivatives::printf_darg0, we set the derivative of printf to be zero.
In principle, we could detect it and return null instead, but this would break cases when printf is used in expressions, e.g. int i = x + printf(...); -> int _d_i = _d_x + 0;. This is an obscure case, but it makes sense to use 0 in such expressions (say, if it is not printf, but some non-differentiable math function, we cannot determine that).

Should we indicate more clearly where does this 0 come from? E.g.

double f_darg0(double x) {
  int printf_darg0_non_differentiable = 0;
  printf_darg0_non_differentiable;
  printf("%d", x);
  return 1;
}

"Must be the same types" strange case

This function is differentiate without problems:

float sphere_implicit_func(float x, float y, float z, float px, float py, float pz, float r) {
  return (x-px)*(x-px) + (y-py)*(y-py) + (z-pz)*(z-pz) - r*r;
}

but this

float sphere_distance_func(float x, float y, float z, float px, float py, float pz, float r) {
  return (x-px)*(x-px) + (y-py)*(y-py) + (z-pz)*(z-pz) - r;
}

generate compile time assertion:

Assertion failed: (lhs_derived->getType() == rhs_derived->getType() && "Must be the same types."), function VisitBinaryOperator, file /Users/alexanderpenev/clad/src/tools/clad/lib/Differentiator/DerivativeBuilder.cpp, line 453.

This workaround

float sphere_distance_func(float x, float y, float z, float px, float py, float pz, float r) {
  return (x-px)*(x-px) + (y-py)*(y-py) + (z-pz)*(z-pz) - r*1.0f;
}

also is OK.

x/=r; does not differentiate properly

double hyperbolic_octahedron_func(double x, double y, double z, const Vec &p, double r) {
  x/=r; y/=r; z/=r;                        
  return pow(x, 2/3) + pow(y, 2/3) + pow(z, 2/3) - 1;
}

generate:

double hyperbolic_octahedron_func_dx(double x, double y, double z, const Vec &p, double r) {
    1. /= 0.;
    0. /= 0.;
    0. /= 0.;
    return pow_dx(x, 2 / 3);
}

but

double hyperbolic_octahedron_func(double x, double y, double z, const Vec &p, double r) {
  x=x/r; y=y/r; z=z/r;                        
  return pow(x, 2/3) + pow(y, 2/3) + pow(z, 2/3) - 1;
}

is OK, generating:

double hyperbolic_octahedron_func_dx(double x, double y, double z, const Vec &p, double r) {
    x = r / (r * r);
    y = 0.;
    z = 0.;
    return pow_dx(x, 2 / 3);
}

Wrong virtual method override called.

class A {
public:
  virtual float vm(float x, float y) {
    return x + y;
  }
};
class B : public A {
public:
  float vm(float x, float y) override {
    return x*x + y*y;
  }
};
int main() {
  auto vm_darg0_A = clad::differentiate(&A::vm, 0);
  printf("Result is = %f\n", vm_darg0_A.execute(a, 2, 3)); // CHECK-EXEC: Result is = 1.0000
  auto vm_darg0_B = clad::differentiate(&B::vm, 0);
  printf("Result is = %f\n", vm_darg0_B.execute(b, 2, 3)); // CHECK-EXEC: Result is = 4.0000
}

For some reason vm_darg0_A.execute(a, 2, 3) calls the B::vm. If the execute()-s get swapped it works as expected.

Position of output parameter in gradient calls

Now we pass the output _result array to gradient functions as a last parameter.
This causes some problems, e.g.,
default parameter values cannot be used as the _result parameter is last and has no default value.
This leads to generation of functions like

void f_grad(double x, double y = default_y, double* _result);

where there is no way to actually call it with default_y value.

Some problems may also arise when f or f_grad is declared as variadic template function.

Are there better solutions?
Does it make sense to pass _result as a first parameter?

Clad work right only when include "Differentiator.h" is first include in program

This work:

// Necessary for clad to work include
#include "clad/Differentiator/Differentiator.h"
#include <math.h>
#include <stdlib.h>
#include <stdio.h>

This generate compile time assert ("Cannot find builtin derivatives!"):

#include <math.h>
#include <stdlib.h>
#include <stdio.h>
// Necessary for clad to work include
#include "clad/Differentiator/Differentiator.h"

Testsuite fails when clad is built against Debug clang

Failing Tests (3):
    clad :: FirstDerivative/TemplateFunction.C
    clad :: Gradient/Gradients.C
    clad :: Gradient/TestAgainstDiff.C

  Expected Passes    : 18
  Expected Failures  : 3
  Unexpected Failures: 3
lldb -- /Users/vvassilev/workspace/builds/llvm-root-debug/obc/bin/clang-5.0 -cc1 -triple x86_64-apple-macosx10.12.0 -Wdeprecated-objc-isa-usage -Werror=deprecated-objc-isa-usage -emit-obj -mrelax-all -disable-free -main-file-name Gradients.C -mrelocation-model pic -pic-level 2 -mthread-model posix -mdisable-fp-elim -masm-verbose -munwind-tables -faligned-alloc-unavailable -target-cpu penryn -target-linker-version 305 -dwarf-column-info -debugger-tuning=lldb -resource-dir /Users/vvassilev/workspace/builds/llvm-root-debug/obc/bin/../../../../../llvm-root-debug/lib/clang/5.0.0/ -I /Users/vvassilev/workspace/sources/root-llvm/tools/clad/test/Gradient/../../include -stdlib=libc++ -std=c++11 -fdeprecated-macro -fdebug-compilation-dir /Users/vvassilev/workspace/builds/llvm-root-debug/obc -ferror-limit 19 -fmessage-length 117 -stack-protector 1 -fblocks -fobjc-runtime=macosx-10.12.0 -fencode-extended-block-signature -fcxx-exceptions -fexceptions -fmax-type-align=16 -fdiagnostics-show-option -fcolor-diagnostics -add-plugin clad -plugin-arg-clad -fdump-derived-fn -load /Users/vvassilev/workspace/builds/llvm-root-debug/obc/./lib/clad.dylib -o /var/folders/cp/g4fhftbd5c3fsp08mc8nr1pr0000gn/T/Gradients-848bce.o -x c++ /Users/vvassilev/workspace/sources/root-llvm/tools/clad/test/Gradient/Gradients.C

* thread #1, queue = 'com.apple.main-thread', stop reason = signal SIGABRT
    ❲ 0❳ __pthread_kill ❮libsystem_kernel.dylib❯
    ❲ 1❳ pthread_kill ❮libsystem_pthread.dylib❯
    ❲ 2❳ abort ❮libsystem_c.dylib❯
    ❲ 3❳ __assert_rtn ❮libsystem_c.dylib❯
    ❲ 4❳ void clang::CodeGen::CodeGenFunction::EmitCallArgs<clang::FunctionProtoType>(this=0x00007fff5fbf6738, Args=0x00007fff5fbf5588, CallArgTypeInfo=0x000000010d04dc10, ArgRange=iterator_range<clang::Stmt::ConstExprIterator> @ 0x00007fff5fbf4dd0, AC=AbstractCallee @ 0x00007fff5fbf4dc8, ParamsToSkip=0, Order=Default) at CodeGenFunction.h:3770 ❮clang-5.0❯
    ❲ 5❳ clang::CodeGen::CodeGenFunction::EmitCall(this=0x00007fff5fbf6738, CalleeType=QualType @ 0x00007fff5fbf5490, OrigCallee=0x00007fff5fbf5950, E=0x000000010d04df90, ReturnValue=ReturnValueSlot @ 0x00007fff5fbf5480, Chain=0x0000000000000000) at CGExpr.cpp:4468 ❮clang-5.0❯
  * ❲ 6❳ clang::CodeGen::CodeGenFunction::EmitCallExpr(this=0x00007fff5fbf6738, E=0x000000010d04df90, ReturnValue=ReturnValueSlot @ 0x00007fff5fbf59e0) at CGExpr.cpp:4089 ❮clang-5.0❯
    ❲ 7❳ (anonymous namespace)::AggExprEmitter::VisitCallExpr(this=0x00007fff5fbf5ba8, E=0x000000010d04df90) at CGExprAgg.cpp:782 ❮clang-5.0❯
    ❲ 8❳ clang::StmtVisitorBase<clang::make_ptr, (anonymous namespace)::AggExprEmitter, void>::Visit(this=0x00007fff5fbf5ba8, S=0x000000010d04df90) at StmtNodes.inc:329 ❮clang-5.0❯
    ❲ 9❳ (anonymous namespace)::AggExprEmitter::Visit(this=0x00007fff5fbf5ba8, E=0x000000010d04df90) at CGExprAgg.cpp:104 ❮clang-5.0❯
    ❲10❳ clang::CodeGen::CodeGenFunction::EmitAggExpr(this=0x00007fff5fbf6738, E=0x000000010d04df90, Slot=AggValueSlot @ 0x00007fff5fbf5c10) at CGExprAgg.cpp:1548 ❮clang-5.0❯
    ❲11❳ clang::CodeGen::CodeGenFunction::EmitAnyExpr(this=0x00007fff5fbf6738, E=0x000000010d04df90, aggSlot=AggValueSlot @ 0x00007fff5fbf5cf0, ignoreResult=true) at CGExpr.cpp:180 ❮clang-5.0❯
    ❲12❳ clang::CodeGen::CodeGenFunction::EmitIgnoredExpr(this=0x00007fff5fbf6738, E=0x000000010d04df90) at CGExpr.cpp:159 ❮clang-5.0❯
    ❲13❳ clang::CodeGen::CodeGenFunction::EmitStmt(this=0x00007fff5fbf6738, S=0x000000010d04df90) at CGStmt.cpp:107 ❮clang-5.0❯
    ❲14❳ clang::CodeGen::CodeGenFunction::EmitCompoundStmtWithoutScope(this=0x00007fff5fbf6738, S=0x000000010d0512d8, GetLast=false, AggSlot=AggValueSlot @ 0x00007fff5fbf5fa0) at CGStmt.cpp:381 ❮clang-5.0❯
    ❲15❳ clang::CodeGen::CodeGenFunction::EmitCompoundStmt(this=0x00007fff5fbf6738, S=0x000000010d0512d8, GetLast=false, AggSlot=AggValueSlot @ 0x00007fff5fbf60f0) at CGStmt.cpp:371 ❮clang-5.0❯
    ❲16❳ clang::CodeGen::CodeGenFunction::EmitSimpleStmt(this=0x00007fff5fbf6738, S=0x000000010d0512d8) at CGStmt.cpp:344 ❮clang-5.0❯
    ❲17❳ clang::CodeGen::CodeGenFunction::EmitStmt(this=0x00007fff5fbf6738, S=0x000000010d0512d8) at CGStmt.cpp:53 ❮clang-5.0❯
    ❲18❳ clang::CodeGen::CodeGenFunction::EmitCompoundStmtWithoutScope(this=0x00007fff5fbf6738, S=0x000000010d05eb78, GetLast=false, AggSlot=AggValueSlot @ 0x00007fff5fbf6370) at CGStmt.cpp:381 ❮clang-5.0❯
    ❲19❳ clang::CodeGen::CodeGenFunction::EmitFunctionBody(this=0x00007fff5fbf6738, Args=0x00007fff5fbf6600, Body=0x000000010d05eb78) at CodeGenFunction.cpp:1035 ❮clang-5.0❯
    ❲20❳ clang::CodeGen::CodeGenFunction::GenerateCode(this=0x00007fff5fbf6738, GD=GlobalDecl @ 0x00007fff5fbf6588, Fn=0x000000010c81c718, FnInfo=0x000000010c81d5e0) at CodeGenFunction.cpp:1205 ❮clang-5.0❯
    ❲21❳ clang::CodeGen::CodeGenModule::EmitGlobalFunctionDefinition(this=0x000000010d819c00, GD=GlobalDecl @ 0x00007fff5fbf6730, GV=0x000000010c81c718) at CodeGenModule.cpp:3207 ❮clang-5.0❯
    ❲22❳ clang::CodeGen::CodeGenModule::EmitGlobalDefinition(this=0x000000010d819c00, GD=GlobalDecl @ 0x00007fff5fbf7b88, GV=0x0000000000000000) at CodeGenModule.cpp:2037 ❮clang-5.0❯
    ❲23❳ clang::CodeGen::CodeGenModule::EmitGlobal(this=0x000000010d819c00, GD=GlobalDecl @ 0x00007fff5fbf7d70) at CodeGenModule.cpp:1812 ❮clang-5.0❯
    ❲24❳ clang::CodeGen::CodeGenModule::EmitTopLevelDecl(this=0x000000010d819c00, D=0x000000010d04cc50) at CodeGenModule.cpp:3937 ❮clang-5.0❯
    ❲25❳ clang::CodeGeneratorImpl::HandleTopLevelDecl(this=0x000000010ce095b0, DG=DeclGroupRef @ 0x00007fff5fbf97d8) at ModuleBuilder.cpp:322 ❮clang-5.0❯
    ❲26❳ clang::BackendConsumer::HandleTopLevelDecl(this=0x000000010ce093d0, D=DeclGroupRef @ 0x00007fff5fbf9870) at CodeGenAction.cpp:136 ❮clang-5.0❯

(lldb) f 4
❲ 4❳ void clang::CodeGen::CodeGenFunction::EmitCallArgs<clang::FunctionProtoType>(this=0x00007fff5fbf6738, Args=0x00007fff5fbf5588, CallArgTypeInfo=0x000000010d04dc10, ArgRange=iterator_range<clang::Stmt::ConstExprIterator> @ 0x00007fff5fbf4dd0, AC=AbstractCallee @ 0x00007fff5fbf4dc8, ParamsToSkip=0, Order=Default) at CodeGenFunction.h:3770 ❮clang-5.0❯
   3767	                E = CallArgTypeInfo->param_type_end();
   3768	           I != E; ++I, ++Arg) {
   3769	        assert(Arg != ArgRange.end() && "Running over edge of argument list!");
-> 3770	        assert((isGenericMethod ||
   3771	                ((*I)->isVariablyModifiedType() ||
   3772	                 (*I).getNonReferenceType()->isObjCRetainableType() ||
   3773	                 getContext()
(lldb) p isGenericMethod
(bool) $1 = false
(lldb) up
❲ 5❳ clang::CodeGen::CodeGenFunction::EmitCall(this=0x00007fff5fbf6738, CalleeType=QualType @ 0x00007fff5fbf5490, OrigCallee=0x00007fff5fbf5950, E=0x000000010d04df90, ReturnValue=ReturnValueSlot @ 0x00007fff5fbf5480, Chain=0x0000000000000000) at CGExpr.cpp:4468 ❮clang-5.0❯
   4465	    }
   4466	  }
   4467	
-> 4468	  EmitCallArgs(Args, dyn_cast<FunctionProtoType>(FnType), E->arguments(),
   4469	               E->getDirectCallee(), /*ParamsToSkip*/ 0, Order);
   4470	
   4471	  const CGFunctionInfo &FnInfo = CGM.getTypes().arrangeFreeFunctionCall(
(lldb) p E->dump()
CallExpr 0x10d04df90 'CladFunction<false, void, double, double, double *>':'class clad::CladFunction<false, void, double, double, double *>'
|-ImplicitCastExpr 0x10d04df78 'CladFunction<false, void, double, double, double *> (*)(double (*)(double, double), const char *)' <FunctionToPointerDecay>
| `-DeclRefExpr 0x10d04dee0 'CladFunction<false, void, double, double, double *> (double (*)(double, double), const char *)' lvalue Function 0x10d04dd40 'gradient' 'CladFunction<false, void, double, double, double *> (double (*)(double, double), const char *)' (FunctionTemplate 0x10d03cf80 'gradient')
|-ImplicitCastExpr 0x10d066070 'void (*)(double, double, double *)' <FunctionToPointerDecay>
| `-DeclRefExpr 0x10d066048 'void (double, double, double *)' lvalue Function 0x10d065cc0 'f_add1_grad' 'void (double, double, double *)'
`-ImplicitCastExpr 0x10d0661a0 'const char *' <ArrayToPointerDecay>
  `-StringLiteral 0x10d066108 'const char [107]' lvalue "void f_add1_grad(double x, double y, double *_result) {\n    _result[0UL] += 1.;\n    _result[1UL] += 1.;\n}\n"

It looks like some of our newly generated expressions do not set some relevant properties.

@efremale, could you take a look?

Simple constant folder do not work in some cases

For example:

double sphere_implicit_func(double x, double y, double z, const Vec &p, double r) {
  return (x-p.x)*(x-p.x) + (y-p.y)*(y-p.y) + (z-p.z)*(z-p.z) - r*r;
}

generate:

double sphere_implicit_func_dx(double x, double y, double z, const Vec &p, double r) {
    return ((x - p.x) + (x - p.x));
}
double sphere_implicit_func_dy(double x, double y, double z, const Vec &p, double r) {
    return 0. + ((y - p.y) + (y - p.y));
}
double sphere_implicit_func_dz(double x, double y, double z, const Vec &p, double r) {
    return 0. + ((z - p.z) + (z - p.z));
}

Reduce cloning complexity in forward mode

Each call to Clone() has complexity of O(n), where n is the size of the Stmt's subtree, since it recursively calls Clone() on every substatement until leaves are reached.

When differentiating a function via Forward/ReverseModeVisitor, we potentially call Clone() on every node of the AST tree (which recursively clones all subnodes, we are cloning same nodes multiple times) and have the total complexity of O(n^2).

Wouldn't it be more efficient to Clone() the whole tree just once before calling Visit() on it? The complexity should be just 2xO(n).

Recommend Projects

  • React photo React

    A declarative, efficient, and flexible JavaScript library for building user interfaces.

  • Vue.js photo Vue.js

    🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.

  • Typescript photo Typescript

    TypeScript is a superset of JavaScript that compiles to clean JavaScript output.

  • TensorFlow photo TensorFlow

    An Open Source Machine Learning Framework for Everyone

  • Django photo Django

    The Web framework for perfectionists with deadlines.

  • D3 photo D3

    Bring data to life with SVG, Canvas and HTML. 📊📈🎉

Recommend Topics

  • javascript

    JavaScript (JS) is a lightweight interpreted programming language with first-class functions.

  • web

    Some thing interesting about web. New door for the world.

  • server

    A server is a program made to process requests and deliver data to clients.

  • Machine learning

    Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.

  • Game

    Some thing interesting about game, make everyone happy.

Recommend Org

  • Facebook photo Facebook

    We are working to build community through open source technology. NB: members must have two-factor auth.

  • Microsoft photo Microsoft

    Open source projects and samples from Microsoft.

  • Google photo Google

    Google ❤️ Open Source for everyone.

  • D3 photo D3

    Data-Driven Documents codes.