In this article we present a new approach for automatic adjoint differentiation
(AAD) with a special focus on computations where derivatives ∂F(X)
∂X are required for multiple instances of vectors X. In practice, the presented approach
is able to calculate all the differentials faster than the primal (original) C++
program for F.publishe

Goloubentsev, Dmitri

Lakshtanov, Evgeny

Wilmott

English

Repositório Institucional da Universidade de Aveiro

AAD: Breaking the Primal BarrierDmitri Goloubentsev∗, Evgeny Lakshtanov†Here, we present a new approach for Auto-matic Adjoint Differentiation with a special fo-cus on computations where derivatives ∂F (X)∂Xarerequired for multiple instances of vectors X. Inpractice, the presented approach is able to cal-culate all the differentials faster than the primal(original) C++ program for F . Major applicationareas are• Gradient methods for optimisation problems,including global model calibration, speechrecognition, deblurring of images and ma-chine learning in general• Derivatives of mathematical expectation• Path-wise sensitivities of SDECode Transformation vs Operator Overload-ing Currently, two main approaches are used forthe Automatic AD tools:Code Transformation (CT): Analises the com-puter program which implements function F toproduce a code of the AD method.Operator Overloading (OO): All mathematicaloperations are overloaded in such a way that theinformation about a computational graph of F issaved in the data structure called Tape1. Tape isused afterwards to process the backward pass ofthe AD method.There are rather succesful CT AAD tools, how-ever, they limit the available language featuresand make the build system more complex. Theneed for such a tool is further reflected in the fact∗MathLogic LTD, London, UK., dmitri@matlogica.com†CIDMA, Department of Mathematics, University ofAveiro, Aveiro 3810-193, Portugal and MathLogic LTD,London, UK., lakshtanov@matlogica.com1its choice depends on a concrete AAD toolthat there is currently no CT AAD available forC++. The OO approach usually demonstratesweaker speed performance due to a runtime over-head in each iteration. Let us enter into moredetails at this point.Consider a schematic processing of multiplesamples Xi using a standard OO AAD library:Loop i =0. .N−1BeginTapeRecording ( ) ;Y[ i ] = F(X[ i ] ) ;StopTapeRecording ( ) ;dX[ i ] = Reverse (dY[ i ] , tape ) ;Next i ;Each overloaded operator collects informationon valuations, after which the second backwardspass updates the adjoint variables. This approachcomes with a number of disadvantages:• OO or/and tape interpretation runtime over-head for each iteration• Multithreading is only possible if the origi-nal function F(X[i]) is thread-safe• Relying on compiler for CPU vectorizationof a scalar primal function means that vec-torization is used very sparsely• F () may perform unnecessary operations thatdon’t depend on X, and yet are executedat each iteration (e.g. mathematical opera-tions, virtual function calls, dictionary lookupsetc)• Additional memory is required to store theTape data-structure. Its size is proportionalto the number of operations of the primalfunction F , which may be prohibitive in somecases.1OO×CT Our innovative idea is to cross2 bothapproaches, namely to use the Overloaded Oper-ators to auto-generate AD-version of the primalfunction at run-time. The created AD-functionscan be used for different Xi instead of performingclassic OO AAD approach on each F (Xi).During the first run of the original function,every overloaded function - or operator - will gen-erate instructions for forward and reverse AD-pass.For instance, consider the function f(a, b, c) =a ∗ b + c. The first column of the following ta-ble lists the consecutive calls of atomic valuationsduring the execution of the function f(a, b, c).Valuation → Forward() → Reverse()Initialization v0=a; v1=b;Initialization v2=c;operator * v3=v0*v1; d1+ = d3 * v0;d0 += d3 * v1;operator + v4=v3+v2; d2 += d4;d3 += d4;Initialization f=v4;Initialization d0=0; d1=0;Initialization d2=0; d3=0;Concatenated entries of the second column formthe body of the Forward pass, while the body ofthe Reverse pass is formed by the entries in thethird column in reverse order:void AD Function( double a , b , c ,& f ,&d0,&d1,&d2 , d4 ) {double v3 , d3 ( 0 ) ;v0=a ; v1=b ; v2=c ;v3=v0∗v1 ;v4=v3+v2 ;f=v4 ;// Reversed2 += d4 ;d3 += d4 ;d0 += d3 ∗ v1 ;d1 += d3 ∗ v0 ;}Now all differentials of the function f and itsvalue at the point (a, b, c) can be computed bycallingdouble d0 ( 0 ) , d1 ( 0 ) , d2 ( 0 ) , d4 ( 1 ) , f ;2Not to be confused with AAD of the mixed type, wherethe tool does an analysis first and then takes the decisionon what approach should be used in each particular situa-tionAD function (a , b , c , f , d0 , d1 , d2 , d4 )Within the proposed framework, the OO isused for only one sample X, after which the result-ing program can process multiple input data. Letus list the immediate benefits of this approach:• Unlike the classic OO AAD, the AD-functiondoes not change from one iteration to thenext one. Hence there is no any OO ortape interpretation run-time overhead perXi sample.• The AD-function is completely segregatedfrom the user program. All user data areencapsulated within the AD-function, andits memory state is limited to vectors v andd. Thus, multiple samples of X can be pro-cessed safely in the multithreaded mode.• The AD-function can be generated to con-sistently utilize native CPU vectorization toprocess 4(8)-double chunks of user data (AVX2\ AVX512 speed up x4-x8).• That can lead to a final acceleration of order8 × #Cores compared to the modern AADtools.• Highly efficient, i.e. operations that don’tdepend on X are not included in the con-structed AD-function.• Although additional memory is also requiredto store the AD-function, its code remainsstatic and can be shared between CPU cores.Combining all the stated benefits one can achievea fantastic performance of 0.4 at one core, com-paring to the original program implemented us-ing a standard double arithmetics. For anybodywho wishes to play around we prepared a ”pro-totype” C++ implementation available at www.matlogica.com free of charge. This prototypeimplementation is simplified and cannot work ef-fectively with large computations. In addition toa two-pass compilation, the relative compilation2time and the volume of memory involved are un-acceptable for practical use3. Part of the reasonfor this is that the size of AAD functions is propor-tional to the tape i.e. to the linearized (unfolded)version of the primal algorithm. Below we discusshow we can address these drawbacks.A prototype C++ implementation Thecode of the prototype library is short and self ex-planatory. Nevertheless, we found it reasonableto provide some comments in case it needs fur-ther clarificaiton. According to the canonical AD,all variables can be distinguished by their roles asinputs, intermediates, or outputs. The header fileNaiveAADLib.h defines the active class dagdoublewhich transforms each double variable into a nodeof a future calculation Directed Acyclic Graph.c l a s s dagdouble {pub l i c :dagdouble ( ) {}dagdouble ( const double& val ,bool i s I npu t = f a l s e ): va l ( va l ) {indx = getNextVarCounter ( ) ;i f ( ! i s I npu t ) {aadAssignConst ( indx , va l ) ;} e l s e {inputIndex . i n s e r t ( indx ) ;}}dagdouble ( const dagdouble& other ) :va l ( other . va l ) {indx = getNextVarCounter ( ) ;aadAssign ( indx , other . indx ) ;}dagdouble& operator =( const dagdouble& other ) {va l = other . va l ;indx = getNextVarCounter ( ) ;aadAssign ( indx , other . indx ) ;r e turn ∗ t h i s ;}void markAsOutput ( ) {outputIndex . i n s e r t ( indx ) ;}double va l ;i n t indx ;} ;As an example, we supplied an overloaded ver-sion of the multiplication operator. Here, we pro-vide only the scalar version of the code. For thevector version, the reader is encouraged to consultthe Prototype library source code.dagdouble operator ∗( const dagdouble& a ,const dagdouble& b) {dagdouble r e s ;3See the relative compilation time benchmarks for pro-totype and professional AAD-Compilersr e s . va l = a . va l ∗ b . va l ;r e s . indx = getNextVarCounter ( ) ;aadMult ( r e s . indx , a . indx , b . indx ) ;r e turn r e s ;}void aadMult ( i n t r e s indx , i n t a indx ,i n t b indx ) {s t r i ng s t r eam f s t r , r s t r ;i f ( codeVers ion == ScalarCode ) {f s t r << ”v” << r e s i ndx << ”=” << ”v”<< a indx << ”∗v” << b indx << ” ; ” ;r s t r << ”d” << a indx << ”+=” << ”d”<< r e s i ndx << ”∗v” << b indx << ” ;”<< ”d” << b indx << ”+=” << ”d”<< r e s i ndx << ”∗v” << a indx << ” ; ” ;}aad func fwd . push back ( f s t r . s t r ( ) ) ;aad func rev . push back ( r s t r . s t r ( ) ) ;}It creates a new variable (getNextVarCounter())and writes the corresponding code of the Forwardand Reverse passes (applying aadMult()).The distribution contains a basic example main.cppwhere the library‘s functioning is carried out bydefining a preprocessor variable#define USE GENRATED AAD FUNCTION. De-pending on its value, the user program either gen-erates a file which contains a newly created ADfunction, or uses the earlier generated one.The included synthetic benchmark test (exam-ple.cpp) produce the following results:clang Absolute time Relative factorAVX2 Primal 83ms 1AVX2 Adjoint 64ms 0.73AVX512 Primal 84ms 1AVX512 Adjoint 39ms 0.46All of the aforementioned disadvantages of theprototype approach can be addressed by gener-ating a binary code directly from OO. This ideawas implemented in AAD-Compiler by MathLogicLTD.A JIT AAD-Compiler is a professional ver-sion of the prototype library. It represents a com-pletely enabled OO AAD library and offers thefollowing features:• Optimized machine binary code generationoptimized for both run-time performance andquick code generation (patent pending)• Streaming Compilation (patent pending)3• Incremental Checkpointing• Proprietary AD-Code-Folding Compression• Support for multiple platforms and C++compilers• Support for AVX2 and AVX512 vectoriza-tion (patent pending)• Enables multithreaded valuations even whenunderlying user program isn’t multithreadsafe (patent pending)The AAD-Compiler has been extensively testedand documented. Licensing terms distribution canbe found at www.matlogica.com.AAD-Compiler speed-benchmark resultsThe benchmark results for the AAD-Compiler arebased on tests of different nature, including ran-dom synthetic tests, previously published tests suchas LW, Toon or GMM, as well as the standard Fi-nancial models like LMM or Stochastic VolatilityCalibration model [1],[2],[4].The first test we consider is based on the opensource benchmark from the article F.Srajer et al.(2018). Authors test various AAD-tools for theGaussian Mixture Model with 2.5M iterations andget the following results (absolute time of adjointsin seconds):# variables 5.36e+4 4.29e+5Manual 3.89e+2 6.16e+3Finite diff.Adept C++ 4.09e+3 3.99e+4ADOLC C++ 1.04e+4Ceres C++Tapenade C 1.32e+3 1.59e+4DiffSharp F#MuPAD Matlab • •Julia-F JuliaJulia-F (vect) Julia • •Autograd PythonTheano Python • •Theano (vect) Python • •Authors note that ”The bullet symbolizes thata tool crashed and no entry means that a tooldid not finish in the time limit. Only tools thatcould compute at least one problem instance areshown.”We execute the Adept‘s and AAD-Compiler‘sversions of the GMM with 429 thousands argu-ments and 2.5M iterations and get the followingabsolute time of adjoints (in seconds). The AAD-Compiler has optionally activated Code compres-sor, so both cases are presented here.Tool Compressor AVX2 AVX512Adept - 2.89e+4 N/AAADC On 1.3e+3 1.17e+3AADC Off 2.86e+3 2.86e+3Note that the huge number of required differ-entials (0.5M) makes the GMM a memory-boundproblem and memory bandwidth is almost maxedout at AVX2. For smaller sized problems perfor-mance scales well with the AVX vector size.AAD-Compiler Relative time of adjointsfor real world models The following table demon-strates relative time of the full gradient to the ex-ecution time of a primal algorithm.Performance coeff: AVX2 AVX512Toon 0.25 0.2Heston 0.37 0.22LMM 0.35 0.21Toon benchmark is assumed to run over multipleinput data.The behaviour of the performance relative timeon the number of CPU cores can be found in thefollowing table:LMM 1 Core 4 Cores 8 CoresAVX2 0.35 0.095 0.05AVX512 0.21 0.06 0.035AAD-Compiler Relative compilation timefor various real world and toy models Fromthe following table one can conclude that compi-lation time is equivalent to around 400 executionsof the primal algorithm.Model Relative Compilation timeHeston 970LMM 650GMM 700Toon 240AAD-Compiler memory-benchmark re-4sultsFor memory benchmarks we also used the con-ventional Lax&Wendroff and Toon tests from [1].Tool Comp. LW Toon GMMAdept 3mb 24.5mb 23mbAADC ON 42kb 11.5mb 7mbAADC OFF 16mb 59.5mb 15mbAADC’s values are normalized to a single sam-ple. To obtain the real memory consumption, oneneeds to multiply the values by an AVX vectorlength.Forward function. The Forward functionis useful on its own and constitutes a vectorizedreplication of the primal algorithm. Similarly tothe complete AD function, it is MT safe even if theprimal algorithm is not. This replication provesextremely useful if one wishes to use Finite Dif-ferences, or just accelerate any complex computa-tions. Given the current trend in using multipleCPU cores to accelerate computations, this fea-ture is useful in its own right.References[1] Hogan, R. J. 2014. Fast reverse-mode auto-matic differentiation using expression tem-plates in C++. ACM Transactions on Math-ematical Software (TOMS), 40(4), 26.[2] Srajer, F., Kukelova, Z. and Fitzgibbon, A.2018. A benchmark of selected algorithmicdifferentiation tools on some problems incomputer vision and machine learning. Op-timization Methods and Software, 33(4-6),889-906.[3] AAD-Compiler prototype library. www.github.com[4] Goloubentcev, D. and Lakshtanov, E. 2019.Remarks on stochastic automatic adjointdifferentiation and financial models calibra-tion. arXiv preprint arXiv:1901.04200.[5] AAD-Compiler by MathLogics LTD. www.matlogica.com5

AAD: breaking the primal barrier

https://ria.ua.pt/bitstream/10773/26588/1/wilmott.pdf

AAD: breaking the primal barrier

Abstract

Similar works

Full text

Available Versions

Repositório Institucional da Universidade de Aveiro