Kristjan Jonasson
Abstract
A set of Fortran subroutines for reverse mode algorithmic (or automatic) differentiation of the basic linear algebra subprograms (BLAS) is presented. This is preceded by a description of the mathematical tools used to obtain the formulae of these derivatives, with emphasis on special matrices supported by the BLAS: triangular, symmetric, and band. All single and double precision BLAS derivatives have been implemented, together with the Cholesky factorization from Linear Algebra Package (LAPACK). The subroutines are written in Fortran 2003 with a Fortran 77 interface to allow use from C and C++, as well as dynamic languages such as R, Python, Matlab, and Octave. The subroutines are all implemented by calling BLAS, thereby attaining fast runtime. Timing results show derivative runtimes that are about twice those of the corresponding BLAS, in line with theory. The emphasis is on reverse mode because it is more important for the main application that we have in mind, numerical optimization. Two examples are presented, one dealing with the least squares modeling of groundwater, and the other dealing with the maximum likelihood estimation of the parameters of a vector autoregressive time series. The article contains comprehensive tables of formulae for the BLAS derivatives as well as for several non-BLAS matrix operations commonly used in optimization.
Supplemental Material
Available for Download
Software for Fortran Subroutines for Reverse Mode Algorithmic Differentiation of BLAS Matrix Operations
- 2017. Intel Math Kernel Library. Retrieved October 13, 2017 from https://software.intel.com/en-us/mkl.Google Scholar
- 2017. OpenBLAS: An optimized BLAS library. Retrieved December 4, 2017 from http://www.openblas.net/.Google Scholar
- Ai Barraud. 1977. A numerical algorithm to solve AT XA − X = Q. IEEE Transactions on Automatic Control 22, 5 (1977), 883--885.Google ScholarCross Ref
- Thomas F. Coleman and Arun Verma. 1998. ADMAT: An automatic differentiation toolbox for MATLAB. In Proceedings of the SIAM Workshop on Object Oriented Methods for Inter-Operable Scientific and Engineering Computing, SIAM, Philadelphia, PA, Vol. 2.Google Scholar
- Thomas F. Coleman and Wei Xu. 2016. Automatic Differentiation in MATLAB Using ADMAT with Applications. SIAM.Google Scholar
- Jack J. Dongarra, Jeremy Du Croz, Sven Hammarling, and Iain S. Duff. 1990. A set of level 3 basic linear algebra subprograms. ACM Transactions on Mathematical Software (TOMS) 16, 1 (1990), 1--17.Google ScholarDigital Library
- Jack J. Dongarra and David W. Walker. 1995. Software libraries for linear algebra computations on high performance computers. SIAM Review 37, 2 (1995), 151--180.Google ScholarDigital Library
- Paul S. Dwyer and M. S. MacPhail. 1948. Symbolic matrix derivatives. The Annals of Mathematical Statistics 19, 4 (1948), 517--534.Google ScholarCross Ref
- Shaun A. Forth. 2006. An efficient overloaded implementation of forward mode automatic differentiation in MATLAB. ACM Transactions on Mathematical Software (TOMS) 32, 2 (2006), 195--222.Google ScholarDigital Library
- Ralf Giering and Thomas Kaminski. 1998. Recipes for adjoint code construction. ACM Transactions on Mathematical Software (TOMS) 24, 4 (1998), 437--474.Google ScholarDigital Library
- M. B. Giles. 2008. An Extended Collection of Matrix Derivative Results for Forward and Reverse Mode Automatic Differentiation. Technical Report. Oxford University Computing Laboratory.Google Scholar
- Mike B. Giles. 2008. Collected matrix derivative results for forward and reverse mode algorithmic differentiation. In Advances in Automatic Differentiation, Christian H. Bischof, H. Martin Bücker, Paul D. Hovland, Uwe Naumann, and J. Utke (Eds.). Springer, Berlin, 35--44.Google Scholar
- Kazushige Goto and Robert van de Geijn. 2008. High-performance implementation of the level-3 BLAS. ACM Transactions on Mathematical Software (TOMS) 35, 1 (2008), 4.Google ScholarDigital Library
- Andreas Griewank. 1989. On automatic differentiation. In Mathematical Programming, Masao Iri and K. Tanabe (Eds.). Kluwer Academic Publishers, Dordrecht, 83--108.Google Scholar
- Andreas Griewank. 2003. A mathematical view of automatic differentiation. In Acta Numerica. Vol. 12. Cambridge University Press, Cambridge, 321--398.Google Scholar
- Andreas Griewank. 2012. Who invented the reverse mode of differentiation? Documenta Mathematica, Extra Volume ISMP (2012), 389--400.Google Scholar
- Andreas Griewank. 2015. Automatic differentiation, Mark R. Dennis, Paul Glendinning, Paul A. Martin, Fadil Santosa, and Jared Tanner (Eds.). Princeton University Press, 749--752.Google Scholar
- Andreas Griewank, David Juedes, and Jean Utke. 1996. Algorithm 755: ADOL-C: A package for the automatic differentiation of algorithms written in C/C++. ACM Transactions on Mathematical Software 22, 2 (1996), 131--167.Google ScholarDigital Library
- Andreas Griewank and Andrea Walther. 2008. Evaluating Derivatives: Principles and Techniques of Algorithmic Differentiation (2nd ed.). Number 105 in Other Titles in Applied Mathematics. SIAM, Philadelphia, PA.Google ScholarCross Ref
- Masao Iri. 1991. History of automatic differentiation and rounding error estimation. In Automatic Differentiation of Algorithms: Theory, Implementation, and Application, Andreas Griewank and George F. Corliss (Eds.). SIAM, Philadelphia, 1--16.Google Scholar
- Kristjan Jonasson. 2008. Algorithm 878: Exact VARMA likelihood and its gradient for complete and incomplete data with Matlab. ACM Transactions on Mathematical Software (TOMS) 35, 1 (2008), 6:1--6:11.Google Scholar
- Kristjan Jonasson. 2009. Algorithm 892: DISPMODULE, a Fortran 95 module for pretty-printing matrices. ACM Transactions on Mathematical Software (TOMS) 36, 1 (2009), 6:1--6:7.Google ScholarDigital Library
- Gershon Kedem. 1980. Automatic differentiation of computer programs. ACM Trans. Math. Software 6, 2 (June 1980), 150--165.Google ScholarDigital Library
- Chuck L. Lawson, Richard J. Hanson, David R. Kincaid, and Fred T. Krogh. 1979. Basic linear algebra subprograms for Fortran usage. ACM Transactions on Mathematical Software (TOMS) 5, 3 (1979), 308--323.Google ScholarDigital Library
- Uwe Naumann. 2008. Optimal Jacobian accumulation is NP-complete. Mathematical Programming 112, 2 (2008), 427--441.Google ScholarDigital Library
- Uwe Naumann. 2012. The Art of Differentiating Computer Programs: An Introduction to Algorithmic Differentiation. Siam, Philadelphia, PA.Google ScholarDigital Library
- S. P. Smith. 1995. Differentiation of the Cholesky algorithm. Journal of Computational and Graphical Statistics 4, 2 (1995), 134--147.Google Scholar
- Arun Verma. 1998. Structured Automatic Differentiation. Ph.D. Dissertation. Department of Computer Science, Cornell University, Ithaca, NY.Google Scholar
- Arun Verma. 1999. ADMAT: Automatic differentiation in MATLAB using object oriented methods. In Object Oriented Methods for Interoperable Scientific and Engineering Computing: Proceedings of the 1998 SIAM Workshop, M. E. Henderson, C. R. Anderson, and S. L. Lyons (Eds.). SIAM, Philadelphia, 174--183.Google Scholar
- Sebastian F. Walter and Lutz Lehmann. 2013. Algorithmic differentiation in Python with AlgoPy. Journal of Computational Science 4, 5 (2013), 334--344.Google ScholarCross Ref
- Qian Wang, Xianyi Zhang, Yunquan Zhang, and Qing Yi. 2013. AUGEM: Automatically generate high performance dense linear algebra kernels on x86 CPUs. In Proceedings of the International Conference on High Performance Computing, Networking, Storage and Analysis. ACM, 25:1--25:12.Google Scholar
- Matthew J. Weinstein and Anil V. Rao. 2017. Algorithm 984: ADiGator, a toolbox for the algorithmic differentiation of mathematical functions in matlab using source transformation via operator overloading. ACM Transactions on Mathematical Software (TOMS) 44, 2 (2017), 21:1--21:25.Google ScholarDigital Library
- R. Clint Whaley, Antoine Petitet, and Jack J. Dongarra. 2001. Automated empirical optimizations of software and the ATLAS project. Parallel Computing 27, 1 (2001), 3--35.Google ScholarCross Ref
Index Terms
- Algorithm 1005: Fortran Subroutines for Reverse Mode Algorithmic Differentiation of BLAS Matrix Operations
Recommendations
Practical Use of Some Krylov Subspace Methods for Solving Indefinite and Nonsymmetric Linear Systems
The main purpose of this paper is to develop stable versions of some Krylov subspace methods for solving linear systems of equations $Ax = b$. As in the case of Paige and Saunders's SYMMLQ [SIAM J. Numer. Anal., 12 (1975), pp. 617–624], our algorithms ...
Row Modifications of a Sparse Cholesky Factorization
Given a sparse, symmetric positive definite matrix C and an associated sparse Cholesky factorization LDL$\tr$, we develop sparse techniques for updating the factorization after a symmetric modification of a row and column of C. We show how the ...
Solving unsymmetric sparse systems of linear equations with PARDISO
Special issue: Selected numerical algorithmsSupernode partitioning for unsymmetric matrices together with complete block diagonal supernode pivoting and asynchronous computation can achieve high gigaflop rates for parallel sparse LU factorization on shared memory parallel computers. The progress ...
Comments