Argonne - Automatic Differentiation
The goal of the computational differentiation project is to enable the generation of efficient
derivative code for computer models described in Fortran or C, with minimal human effort. We
employ approaches based on automatic differentiation which provide derivatives (in the engineering
context often also referred to as "sensitivities") that are accurate up to working precision, and
usually less compute-intensive than derivative approximations derived from perturbed model runs.
This project is pursuing three main lines of research:
I."Black-box" automatic differentiation
tools. One of our objectives is to create
automatic differentiation applications that are
simple to use. To this end these tools require, in
essence, no knowledge about the nature of the
source code being transformed. We provide
* ADIFORfor Fortran 77 (developed in
collaboration with the Center for Research
on Parallel Computation at Rice University).
http://www-unix.mcs.anl.gov/autodiff/ADIFOR/
http://www-unix.mcs.anl.gov/autodiff/ADIFOR/adifor_demo/
* ADIC for ANSI C.
* SparsLinC library to support sparse first-order automatic differentiation.
"Black-box" automatic differentiation tools.
Given a computer program for f:Rn -> Rm, automatic differentiation produces a program for exactly evaluating the derivatives
by applying the chain rule to elementary instructions.
Traditionally, there are two ways for propagating derivatives:
The Forward Mode
propagates derivative working from independant (input) variables to dependant (output) variables. The time and space
requirements are of the order O(n*cost(f)).
The Reverse Mode
propagates derivatives working from dependent, function variables to independant (input) variables. The time
requirements are of the order O(m*time(f)), but space requirements may be of the order O(#flops(f)), since one must
be able to reverse program flow and remember (or recompute) any intermediate value that nonlinearly impacts the final
result. Snapshotting schemes are a way to address this issue. (Argonne preprint 228) **P228.ps**
These two methods are particular instances of accumulating derivatives, and derivatives may be computed in many different
ways (Argonne preprint 267). **P267.ps** However, no matter which accumulation strategy was employed, AD produces
derivatives accurate to machine precision, without truncation error.
We are working on the following tools to provide automatic differentiation functionality for programs written in Fortran and C:
ADIFOR: Automatic Differentiation of Fortran Codes.
ADIC: Automatic Differentiation of C Codes.
The paper ADIC: An Extensible Automatic Differentiation Tool for ANSI-C ***NA***
***describes the system issues
that need to be addressed in automatic differentiation, illustrates ADIC's particular approach to AD,
provides a short example, and reports on applications of ADIC to a semiconductor device simulator, 3-D
CFD grid generator, vehicle simulator, and neural network code.
The application of ADIC to the CSCMDO 3-D volume grid generator
******************************
http://geolab.larc.nasa.gov/CSCMDO/
The Coordinate and Sensitivity Calculator for Multi-disciplinary Design Optimization CSCMDO is a general purpose
multi-block three-dimensional volume grid generator which is suitable for Multi-disciplinary Design Optimization. The code is
timely, robust, highly automated, and written in ANSI "C" for platform independence. Algebraic techniques are used to
generate and/or modify block face and volume grids to reflect geometric changes resulting from design optimization. Volume
grids are generated/modified in a batch environment and controlled via an ASCII user input deck. This allows the code to be
incorporated directly into the design loop. Volume grids have been successfully generated/modified for a wide variety of
configurations. CSCMDO is being supported and developed by the GEOmetry LABoratory at the NASA Langley Research
Center.
******************************
is decribed in detail in the paper ***P512.ps***
"Experiences with the Application of the ADIC Automatic Differentiation Tool to the CSCMDO 3-D
Volume Grid Generation Code". The application of ADIC to the FCAP circuit simulator from HP Labs is
described in the paper ***P698.ps*** "Statistical On-Chip Interconnect Modelling: An Application of Automatic
Differentiation".
ADIC currently provides for the generation of both the first and second order derivative code.
Second-order capabilities are in described in ***P636.ps*** "Algorithms and Design for a Second-Order Automatic
Differentiaton Module". ADIC also provides some coverage of C++ constructs.
The ADIC project was conceived and initiated by Christian Bischof, and principally developed and
currently led by Lucas Roh.
SparsLinC: A Library to support sparse vector accumulation, they key component in first-order forward-mode based
differentiation for large sparse Jacobians or gradients of partially separable functions (this class includes all functions
with a sparse Hessian).
SparsLinC: A Library for Computing Sparse Linear Combinations
***********************************
http://www-unix.mcs.anl.gov/autodiff/sparslinc.html
The main purpose of SparsLinC is to support the computation of linear combinations of sparse vectors, i.e. for
single and double precision real and complex types. This operation is the key kernel in first-order derivative propagation, with
vectors corresponding to derivative objects, and, ultimately, to rows of the Jacobian to be computed. We note that, in general,
no intermediate derivative object can have more nonzeros than any row of the final Jacobian.
However, even for gradients, which are dense, sparsity can potentially be exploited in the gradient computation. If
has a sparse Hessian, it is a partially separable function and can be represented as where each only
depends on a few of the , and hence each of is sparse. Hence, whether is expressed in this partially separable form or
not, there is "sparsity 'under the rug".
So by using dynamic data structures, we can
*reduce flops/storage for derivative computations, and
*automatically derive the sparsity pattern of the Jacobian.
SparsLinC supports sparse vector accumulation employing a polyalgorithm that transparently switches between three different
representations of sparse vectors, provides a buffered memory-management scheme, and special support for the type of
vector accumulations arising in the computation of gradients of partially separable functions. SparsLinC is written in portable
ANSI C and supports single, double, complex, and double complex arithmetic.
Employing ADIFOR/SParsLinC on some large-scale problems arising in large-scale optimization, we obtained up to three
orders of magnitude improvement in both runtime and storage compared to the performance of ADIFOR without
SParsLinC (which usually is already 2-4 times faster than numerical approximations of derivatives).
***********************************
Purse: A tool to provide flexible handling of Fortran intrinsics such as sqrt(x) for which the derivative at x=0 is
undefined and thus requires special attention.
II.Strategies to exploit program structurefor more efficient use of black-box tools.
A small amount of knowledge about the source code can be leveraged to provide a large
decrease in computational complexity of the derivative code.
While derivative codes generated by automatic differentiation tools such as ADIFOR or ADIC usually outperform
divided-difference approximations, the code they produce is often much less efficient than one produced by hand. This
shortcoming may be overcome by utilizing a programmer's insight into the high-level structure of the computation.
Let us consider the situation where we have two codes computing: F : x(n) -> y(m) and
G : y(m) -> z(k) and we want to compute dz/dx. If we simply employ ADIFOR or ADIC, then we first compute F' = dy/dx
and then the directional derivative dy/dx · F'. The complexity of the generated derivative code is O(n · [cost(F) + cost(G)]).
On the other hand, we can separately compute dy/dx and dz/dy and then multiply the two matrices to arrive at dz/dx . The
complexityof this approach is O(n·cost(F) + m·cost(G)). We see that the latter approach is preferable when m « n, a
situation we call interface contraction, since the amount of information flowing into F is greater than the amount of
information passed on to G.
III.Mathematics relating automatic differentation and "do-what-I-mean."
For example, do the derivatives computed in an autodifferentiated version of an iterative
solver converge, and if so, what do they converge to?
********************************************************
http://www-unix.mcs.anl.gov/autodiff/advsdwim.html
The derivatives computed by automatic differentation are those defined by the statements that were executed by a particular
program run, and, in that sense, are very well-defined. One should ask oneself the question then whether these derivatives are
also those that are "desired", and whether the convergence of the original code implies the convergence of the derivative code.
Consider for example the solution of a system of nonlinear equations, where R* is implicitly defined by R(z,x*)=0
The iterative process employed to solve this nonlinear equation is
for m=1,... do
if ((||R(zm,x*)|| <= tolR or (||zm-zm-1|| <= tolz)) then stop
Compute a preconditioner Pm .
zm+1 <- zm - PmR(zm,x*)
end for
Domain-specific knowledge is employed in defining Pm and computing its action onto R(zm,x*) The preconditiner Pm is
usually not formed explicitly.
When we want to compute z*', (' =ð/ðx) , then the application of automatic differentiation to the previous iteration results in
the following code (except for the code in astrisks):
for m=1,... do
if ((||R(zm,x*)|| <= tolR or (||zm-zm-1|| <= tolz)) then stop
** if ((||R'(zm,x*)|| <= tol ðR or (||z'm-z'm-1|| <= tolðz)) then stop**
Compute Pm and P'm
zm+1 <- zm - PmRm.
z'm+1 = z'm - Pm R'm - P'mRm.
end for
We make a few observations:
The implicit function theorem states that
R'(z*,x*) = ðR/ðz · z*' + ðR/ðx = 0
that is, only first-order derivatives are needed to determine z*'.
The univariate derivative R'(z*,x*) can be evaluated cheaply (at approximately the cost of R ).
The computation of P'm is probably computationally dominant. Is it useful in light of the implicit function theorem?
AD does not generate the augmented stopping criterion shown in the box.
As explained in detail in (see Derivative Convergence for Iterative Solvers, A. Griewank, C.Bischof, A. Carle, G. Corliss,
and K. Williamson, ANL-MCS-P333-1192.),
***P333.ps***
we can prove that for a large class of iterative scheme, the derivatives are
likely to converge. However, it is also shown that the derivatives z' may converge slower than the z themselves, necessitating
the augmented stopping criterion based on R' which is NOT generated automatically by AD. Hence, the process is not
altogether automatic any more.
We summarize the key results of this paper:
zm' -> z*' for large class of iterations, with or without Pm'.
Ignoring Pm' term allows z use of the same solver as in the iteration for , and is likely to avoid numerical difficulties
(e.g., Pm' -> ).
The convergence of z' may lag behind that of z and may be of a lesser order. Hence an augmented stopping criterion
monitoring ||R'|| is critical.
No matter how derivative approximation zm' was generated, if vanishes, then zm' is
the tangent of the perturbed solution set
*********************************************************
For an overview of the way this work impacts science and engineering:
Highlights of automatic differentiation applications using ADIFOR and ADIC.
ADIFOR and ADIC have been applied to application codes from various domains of science and engineering.
Atmospheric Chemistry http://www-unix.mcs.anl.gov/autodiff/highlights/matchm.html
On-Chip Interconnect Modeling ***P698.ps*** (not available here)
Mesoscale Weather Modeling ***P532.ps***
CFD Analysis of the High-Speed Civil Transport
Automated multidisciplinary design of aircraft requires the optimization of comlex performance objectives with respect to a
number of design parameters and constraints. The effect of these independent design variables on the system performance
criteria can be quantified in terms of sensitivity derivatives for the indivudal disipline simulation codes. Typical advanced CFD
codes do not provide such derivatives as part of a flow solution. These derivatives are expensive to obtain by
divided-differences from perturbed solutions, and may be unreliable, particularly for noisy functions. In collaboration with
researchers at NASA Langley, we have applied ADIFOR to advanced 3-D CFD codes.
To see a picture of the sensitivity of pressure with respect to root chord on a HSCT24E blended wing-body configuration,
click here (28k) (Plot courtesy of Vamshi Mohan Korivi, mohan@tab00.larc.nasa.gov). Details of this investigation are
described in detail in "Aerodynamic Optimization Studies using a 3-D Supersonic Euler Code with Efficient Calculation of
Sensitivity Derivatives", V. M. Korivi, P. Newman, A. Taylor, Proc. 5th AIAA/NASA/USAF/ISSMO Symp. on
Multidisciplinary Analysis and Optimization, AIAA-94-4270-CP).
Rotorcraft Flight
Research Objective: Automatic Differentiation (AD) is a technique which transforms an existing computer program into a
new program which performs sensitivity analysis in addition to the original analysis. If the original FORTRAN program
calculates a set of dependent (output) variables from a set of independent (input) variables, the new FORTRAN program
calculates the derivatives of the dependent variables with respect to the independent variables. AD is not an automatic
implementation of finite differencing which would produce approximate derivatives. Nor is it related to symbolic manipulation
which requires reprogramming in a special purpose language and which results in convoluted expressions for the derivatives.
Rather, AD is a systematic implementation of the chain-rule of differentiation which produces derivatives to machine accuracy
at a cost generally lower than finite differencing. The objective of this research is to apply AD to typical engineering analysis
codes and evaluate the resulting sensitivity analysis codes for use in design and optimization.
Approach: To implement AD, one modifies any ANSI standard FORTRAN code by inserting special lines of code which
identify the independent and dependent variables. The program is then further modified automatically by a preprocessor which
augments the code to calculate derivatives. The modified program is compiled conventionally, linked with the ADIFOR library
and executed. The ADIFOR preprocessor must execute on a SUN workstation, however the new FORTRAN code can be
compiled and executed on any computer.
Accomplishment Description: The initial phase of the research demonstrated the use of the preprocessor code ADIFOR,
developed by Rice University and Argonne National Laboratory. It was applied to the rotorcraft codes HOVT (a hover
analysis code) and WOPWOP ( a rotor acoustics code) and a variety of other engineering analysis codes. The derivatives
produced by the preprocessed versions of HOVT and WOPWOP were similar to the finite difference derivatives but are not
dependent upon the proper step size choice (see figure) .
Significance: Engineering design and optimization require sensitivity information. For some disciplines, analytical techniques
exist to generate those sensitivities; for most disciplines, however, those derivatives have not been developed. AD is a viable
alternative for generating those sensitivities. It generates exact derivatives at a cost which is competitive with the finite
difference approximations.
Future Plans: NASA is supporting the development of ADIFOR. The next version will fully implement complex arithmetic
and will recognize a larger subset of common FORTRAN syntax. This will facilitate ADIFOR processing of a comprehensive
rotor analysis code like CAMRAD/JA and the testing of ADIFOR derivatives in the IADS procedure.
3-D Groundwater Contaminant Transport
3-D Grid Generation for the High-Speed Civil Transport
The Coordinate and Sensitivity Calculator for Multidisciplinary Design Optimization (CSCMDO) code is a three-dimensional volume grid generator. It provides a
rapid and highly automated grid generation capability which produces changes in the surface and volumetric grids to reflect the perturbations of the baseline system.
Algebraic techniques are used to generate and/or modify block face and volume grids to reflect geometric changes resulting from design optimization. Details are
described in "A Grid Generation System for Multidisciplinary Design Optimization", William T. Jones and Jamshid Samareh-Abolhassani, NASA Langley Internal
Memorandum, 1994. The experience of using the prototype ADIC tool is described in "Experiences with the Application of the ADIC Automatic Differentiation
Tool to the CSCMDO 3-D Volume Grid Generation Code." ***P512.ps*** (in folder)
For sensitivity images of the most recent version of CSCMDO, generated via ADIC, press HERE.
Sensitivity of Volume X wrt Sweep Angle. http://geolab.larc.nasa.gov/CSCMDO/Pictures/GIF/ad1.gif
Sensitivity of Volume Y wrt Sweep Angle. http://geolab.larc.nasa.gov/CSCMDO/Pictures/GIF/ad2.gif
Sensitivity of Volume X wrt Bevel Angle. http://geolab.larc.nasa.gov/CSCMDO/Pictures/GIF/ad3.gif
Sensitivity of Volume Y wrt Bevel Angle. http://geolab.larc.nasa.gov/CSCMDO/Pictures/GIF/ad4.gif
Sensitivity of Volume X wrt Length. http://geolab.larc.nasa.gov/CSCMDO/Pictures/GIF/ad5.gif
Sensitivity of Volume Y wrt Length. http://geolab.larc.nasa.gov/CSCMDO/Pictures/GIF/ad6.gif
A Numerically Complicated Statistical Function -- the Log-Likelihood for log-F distribution (LLDRLF).
A Comparison of ADIFOR-generated versus Hand-Generated Derivatives for a Complicated Statistical Function
Authors: Barry W. Brown and Kathy Russell (University of Texas, M.D. Anderson Cancer Center)
Code Processed by: Alan Carle and Mike Fagan (Center for Reserach on Parallel Computaton (Rice University)
Goals
We applied the ADIFOR 2.0 automatic differentiation system to produce code to compute derivatives for a complicated statistical function: the log-likelihood for
log-F distribution (LLDRLF). The ADIFOR system generated code to compute first and second derivative functions from the statistical function code. The authors
of the statistical function code also ``hand-coded'' a version of the derivatives.
The likelihood function calculation is exceptionally difficult to calculate numerically due to an exceptionally wide input domain. Consequently, the strategy for
computing this function is complicated. One of five different methods (Extreme Value Method, Bratio Method, Asymptotic Method, Continued Fraction Method,
Series Method) is used to compute the function approximation within regions of the input domain. Thus, numerical accuracy is maintained over the entire space of
input values. The hand-coded derivatives follow the same basic strategy -- one of five different methods (Extreme Value Method, Asymptotic Method, Continued
Fraction Method, Series Method, Direct Differentiation Method) is used to compute the derivative approximation. Although four of the function approximation
methods (all but the Bratio Method) correspond directly to four of the derivative approximation methods (all but the Direct Differentiation Method), the authors of
LLDRLF partitioned the input domain for the function and derivative computations differently.
Results
The initial application of ADIFOR to LLDRLF resulted in first derivatives that were accurate to 6 figures in 99.9% of the tested cases, and second derivatives that
were accurate in 97.0% of the cases, in double precision. In quad precision, first derivatives were accurate to 6 figures in 99.9% of the tested cases, and second
derivatives were accurate in 97.7% of the cases. Table 1 and Table 2 show the initial results of the comparison for double and quad precision. In each case where
there was a discreprency between the hand-code derivatives and the ADIFOR-generated derivatives, the presence of the discrepency could be determined by
looking at the output of the ADIFOR-generated code -- ``SQRT'' and ``ABS'' exception handler messages were generated for the expression discrepancies, and
IEEE ``NAN'' and ``INF'' values were returned for the numerical discrepancies. In addition, by invoking ADIFOR with the ``AD_TRACE_RESULTS'' option
enabled, it was easy to detect where each of the numerical discrepancies occurred in the code. To deal with the discrepencies, we made two one-line changes to the
LLDRLF source code and achieved double precision first derivatives and quad precision first and second derivatives accurate to 6 figures in 100% of the tested
cases, and double precision second derivatives that were accurate in 99.4% of the cases. In all of the cases where the ADIFOR-generated second derivatives had
fewer than 6 figures of accuracy, the ADIFOR-generated code provided at least 1 digit of accuracy, and averaged 3.9 digits of accuracy. Table 3 and Table 4 show
the final results of the comparison for double and quad precision.
Lessons
Based on the results of this study, users of ADIFOR are advised to take a close look at any messages generated by ADIFOR's exception handler. The cause of
IEEE exceptions which appear when the ADIFOR-generated code is executed, but which did not appear in the original function evaluation, should also be examined
carefully.
This study also appears to indicate that ADIFOR should probably be more careful about how it generates derivative code for ``log.'' Unfortunately, to maximize the
accuracy of derivative computations, ADIFOR would need to generate several different versions of the derivative code and then select the best version to execute
based on the value of the argument to ``log'' and the values of the derivatives of the argument with respect to the independent variables.
For more information:
Table 1: http://www.caam.rice.edu/~carle/AdiforApplications/LLDRLF/report/P001.html
Table 2: http://www.caam.rice.edu/~carle/AdiforApplications/LLDRLF/report/P002.html
Table 3: http://www.caam.rice.edu/~carle/AdiforApplications/LLDRLF/report/P003.html
Table 4: http://www.caam.rice.edu/~carle/AdiforApplications/LLDRLF/report/P004.html
A somewhat more in-depth description of this effort is described in CRPC-TR95526 available from the CRPC Online Technical Report Database.
http://softlib.rice.edu/CRPC/softlib/TRs_online.html
The LLDRLF code is available by ftp at odin.mda.uth.tmc.edu in /pub/accflf.
Technical Reports of the computational differentiation project, available via anonymous ftp.
http://www-unix.mcs.anl.gov/autodiff/tech_reports.html
A Bibliography of Publications by usersof our automatic differentiation tools.
http://www-unix.mcs.anl.gov/autodiff/testimonials.html
On-Line Documentation
http://www-unix.mcs.anl.gov/autodiff/online.html
Other
http://www.ticam.utexas.edu/~abate/hessian/index.html
http://www-dsed.llnl.gov/files/programs/unix/latex2html/manual/
http://www.ticam.utexas.edu/~abate/papers/ADIPARS/ADIPARS.html
http://www-unix.mcs.anl.gov/alice/
http://www.ticam.utexas.edu/~abate/slides/ADOverview/index.htm
http://www.ticam.utexas.edu/~abate/slides/ADOverview/sld001.htm
http://www-fp.mcs.anl.gov/division/publications/
Automatic differentiation (AD):
AD is a technique for computing numerical values of derivatives which are accurate to machine precision, based on the chain-rule of calculus. There are
several areas of research that I am persuing:
Techniques for computing second-order derivatives via AD:
I have studied several techniques for computing various forms of Hessians using automatic differentiation, and have extended both the ADIFOR and
ADIC tools to support second-order derivatives.
Source transformation AD techniques for Fortran 90:
As part of my Ph.D. thesis, I am studying methods for augmenting Fortran 90 code to compute derivatives. This involves study of how to correctly
handle the new constructs in Fortran 90 and what performance improvements can be obtained using these features. It is expected that the use of these
higher-level constructs should make it easier to produce efficient derivative codes than the corresponding Fortran 77 code. As this work progresses, it
will be documented here. As part of this work, I will develop a prototype AD tool for Fortran 90. It will not support the full language (in particular, lots
of legacy Fortran 77 features are being left out). The supported features most closely match those of F.
The development of the Automatic differentiation Intermediate Form (AIF):
AIF is a language-neutral representation of programs which is intended to glue together various pieces of AD tools, which allows components to be
written once and used many times. Currently, support is in place for augmentation at the statement level. The first AIF module is the Hessian module
that I wrote. We are currenlty working on extending AIF to higher-level program structures such as basic blocks and loops. This will allow the
augmentation algorithms to exploit these program structures and produce more efficient derivative codes.
For more information on AD, see the Computational Differentiation web page at Argonne.
Purpose
Automatic differentiation (AD) tools make it easy to differentiate complicated simulation codes and have been applied to programs such as CFD solvers for 3-D
turbulent flow and 3-D semiconductor device models.
Typically, these programs do not compute the mathematical function of interest, but rather a computational approximation to that function. Numerical analysis and
computational experimentation are employed to show that the approximations thus computed are close to the mathematical model in a suitable norm. A fundamental
question is, then, how the derivatives of such approximate computer models (as they are computed by automatic differentiation and, to some extent, divided
differences) relate to the derivatives of the mathematical function of interest. Specific issues include the following:
1.What is the convergence behavior of derivatives of programs employing an iterative scheme? We have experienced instances where AD applied to
iterative solvers resulted in derivative oscillation. We also have experienced different behavior depending on the preconditioners that were differentiated.
2.How do internal feedback loops such as adaptive stepsize control in numerical integrators affect program derivatives? Experiments show that
other numerical artifacts in integrators also play a role that is not fully understood.
3.Are derivatives of discretized equations good approximations to mathematical derivatives when one is interested in modeling shocks, for
example? Other researchers are considering extended algebraic systems such as interval arithmetic or remainder differential algebras that have a greater
mathematical richness than the usual IEEE floating-point arithmetic. Are derivatives ``safer to compute'' in such contexts?
4.What are suitable methods of differentiating programs with discontinuities? In a computer program, potentially every branch introduces a ``kink''.
Researchers in numerical methods also must grapple with this issue in order to define sensible derivatives for systems exhibiting shocks or discontinuities.
5.The derivatives of a discretized function can be more or less ``smooth'' than the derivatives of the original function. One consequence is that
discretizing, then differentiating may not produce the same results as differentiating, then discretizing. Furthermore, due to the limitations of standard computer
arithmetic, both may be poor approximations to the ``true'' derivatives. How can we overcome these limitations? Is it always better to apply AD to a
discretized function than to discretize the derivatives?
By bringing together an international group of researchers who have investigated these issues from a perspective of automatic differentiation and numerical modeling,
we hope to develop a better understanding of the ``program derivatives'' in numerical modeling and their relation to function derivatives derived in continuous and
often mathematically richer spaces.