Introducing a Julia package to efficiently differentiate an objective function defined implicitly by the solution of a large PDE system.
Steady-state systems of nonlinear partial differential equations (PDEs) are common in engineering and the biogeosciences. These systems are typically controlled by parameters that can be estimated efficiently using second-order optimization algorithms. However, computing the gradient vector and Hessian matrix of a given objective function defined implicitly by the solution of large PDE systems is seldom economical.
A fast and easy-to-use algorithm is introduced for computing the gradient and Hessian of an objective function implicitly constrained by a steady-state PDE system. The algorithm, which is based on the use of hyperdual numbers, is called the F-1 algorithm, because it requires only one factorization of the constraint-equation Jacobian. Careful examination of the relationships that arise from differentiating the PDE system reveal analytical shortcuts that the F-1 algorithm leverages. Benchmarks of the F-1 algorithm against five numerical differentiation schemes are shown in the context of optimizing a global steady-state model of the marine phosphorus cycle that depends explicitly on m=6 parameters. In this context, the F-1 algorithm computes the Hessian 16 to 100 times faster than other algorithms, allowing for the entire optimization procedure to be performed 4 to 26 times faster. This is because other algorithms require O(m) to O(m²) factorizations, which suggests even greater speedups for larger problems.
A live demonstration of using the F-1 algorithm, which is implemented as a Julia package, is given.