[M-GI29] Data assimilation: A fundamental approach in geosciences

[MGI29-P03] Symplectic-adjoint-based uncertainty quantification method for large-scale data assimilation problems

*Shin-ichi Ito¹, Takeru Matsuda³, Yuto Miyatake² (1.The University of Tokyo, 2.Osaka University, 3.RIKEN)

Keywords:Uncertainty quantification, Second-order adjoint method, Symplectic structure

Computing the derivative of an objective function whose explanatory variables are constrained by differential equations is an essential task in various statistical-related fields, such as Hamiltonian Monte Carlo method, neural networks, and four-dimensional variational data assimilation (4DVarDA). In particular, the Hessian matrix, which is the second-order derivative of such an objective function, is important since its inverse matrix elements are directly linked to the uncertainty quantification of the optimal solution. However, as its computational costs are expensive, there is a need for a framework of fast and accurate computations of the Hessian matrix even for large-scale problems.
To date, we have promoted the development of algorithms based on the second-order adjoint (SOA) method used in 4DVarDA and its applied research to conduct high-speed/high-accuracy computations of the Hessian matrix and its inverse matrix elements. The SOA model used in the SOA method is given by a set of ordinary differential equations, of which numerical integrations allow us fast evaluations of the Hessian matrix elements. However, the accuracy of the Hessian matrix elements largely depends on the selection of the numerical integrators to be applied to the SOA model, and there existed no mathematical guideline for the selection of the integrators to give optimal accuracy.
For this background, we propose a construction way to obtain optimal numerical integrators that allow us to suppress the numerical errors up to the machine errors, based on the theory of symplectic geometry. Our method constructs numerical integrators that preserve the invariant inherent in the set of ordinary differential equations appearing in the SOA method even after discretization to enable exact Hessian matrix computations. The integrators obtained by our method also ensure to achieve optimal memory efficiency, not only the optimal accuracy of the Hessian. Through numerical experiments using initial-value estimation problems and parameter estimation problems of reaction-diffusion systems and wave equation systems, we verified that the numerical integrator proposed by our method drastically suppressed numerical errors included in the Hessian matrix compared to that of the conventionally-used numerical integrators.