Complex step: a different perspective on solving over-determined non-linear systems of equations

Abstract

One of the main challenges in various areas of engineering and science is solving over-determined non-linear systems of equations, largely due to the high computational and memory demands required to compute the Jacobian matrix. This work focuses on a framework that employs the non-symmetric LSQR (ns-LSQR) algorithm. The core of our proposal addresses the approximation of the Jacobian-vector product (JVP) and the Jacobian-transposed-vector product (JTVP). For the JVP, the Complex Step method is used, demonstrating its superior accuracy and numerical stability compared to finite differences. For the JTVP, it is shown that naive approximation strategies can lead to the collapse of the Krylov subspace. To overcome this, the Random Heuristic is proposed and validated, comparing it against the Always Best Heuristic and Vanilla Heuristic. This method uses a random subset of canonical or wavelet basis vectors to construct a robust JTVP approximation. Numerical experiments validate that the Random Heuristic prevents the collapse of the Krylov subspace, enabling the ns-LSQR method to converge even with a limited computational budget. Furthermore, it is demonstrated that using a wavelet basis can accelerate convergence in specific problems by providing a sparser representation of the JTVP.

Uno de los principales retos en varias áreas de la ingeniería y ciencias es resolver sistemas de ecuaciones no lineales sobredeterminados, en gran parte debido a las altas demandas de cómputo y memoria necesarias para calcular la matriz Jacobiana. Este trabajo se centra en un marco de trabajo que utiliza el algoritmo no-simétrico LSQR (ns-LSQR). El núcleo de nuestra propuesta aborda la aproximación de los productos Jacobiano-vector (JVP) y Jacobiano-traspuesto-vector (JTVP). Para el JVP, se emplea el método de Complex Step, demostrando su superioridad en precisión y estabilidad numérica frente a las diferencias finitas. Para el JTVP, se evidencia que las estrategias de aproximación triviales pueden llevar al colapso del subespacio de Krylov, por lo que se propone y valida la Heurística Aleatoria, comparándola con las Heurísticas Siempre Mejor y Vanilla. Este método utiliza un subconjunto aleatorio de vectores base canónicos o de bases wavelet para construir una aproximación robusta del JTVP. Los experimentos numéricos validan que la Heurística Aleatoria previene el colapso del subespacio de Krylov, permitiendo la convergencia del método ns-LSQR incluso con un presupuesto computacional limitado. Además, se demuestra que el uso de una base wavelet puede acelerar la convergencia en problemas específicos al ofrecer una representación más dispersa del JTVP.