Peter Benner

Peter Benner, Pawan Goyal, Serkan Gugercin

We investigate the optimal model reduction problem for large-scale quadratic-bilinear (QB) control systems. Our contributions are threefold. First, we discuss the variational analysis and the Volterra series formulation for QB systems. We then define the $\mathcal H_2$-norm for a QB system based on the kernels of the underlying Volterra series and also propose a truncated $\mathcal H_2$-norm. Next, we derive first-order necessary conditions for an optimal approximation, where optimality is measured in term of the truncated $\mathcal H_2$-norm of the error system. We then propose an iterative model reduction algorithm, which upon convergence yields a reduced-order system that approximately satisfies the newly derived optimality conditions. We also discuss an efficient computation of the reduced Hessian, using the special Kronecker structure of the Hessian of the system. We illustrate the efficiency of the proposed method by means of several numerical examples resulting from semi-discretized nonlinear partial differential equations and show its competitiveness with the existing model reduction schemes for QB systems such as moment-matching methods and balanced truncation.

Maximilian Behr, Peter Benner, Jan Heiland

The differential Sylvester equation and its symmetric version, the differential Lyapunov equation, appear in different fields of applied mathematics like control theory, system theory, and model order reduction. The few available straight-forward numerical approaches if applied to large-scale systems come with prohibitively large storage requirements. This shortage motivates us to summarize and explore existing solution formulas for these equations. We develop a unifying approach based on the spectral theorem for normal operators like the Sylvester operator $\mathcal S(X)=AX+XB$ and derive a formula for its norm using an induced operator norm based on the spectrum of $A$ and $B$. In view of numerical approximations, we propose an algorithm that identifies a suitable Krylov subspace using Taylor series and use a projection to approximate the solution. Numerical results for large-scale differential Lyapunov equations are presented in the last sections.

Peter Benner, Christian Himpe, Tim Mitchell

The identification of reduced-order models from high-dimensional data is a challenging task, and even more so if the identified system should not only be suitable for a certain data set, but generally approximate the input-output behavior of the data source. In this work, we consider the input-output dynamic mode decomposition method for system identification. We compare excitation approaches for the data-driven identification process and describe an optimization-based stabilization strategy for the identified systems.

Kapil Ahuja, Peter Benner, Eric de Sturler et al.

Krylov subspace recycling is a process for accelerating the convergence of sequences of linear systems. Based on this technique, the recycling BiCG algorithm has been developed recently. Here, we now generalize and extend this recycling theory to BiCGSTAB. Recycling BiCG focuses on efficiently solving sequences of dual linear systems, while the focus here is on efficiently solving sequences of single linear systems (assuming non-symmetric matrices for both recycling BiCG and recycling BiCGSTAB). As compared with other methods for solving sequences of single linear systems with non-symmetric matrices (e.g., recycling variants of GMRES), BiCG based recycling algorithms, like recycling BiCGSTAB, have the advantage that they involve a short-term recurrence, and hence, do not suffer from storage issues and are also cheaper with respect to the orthogonalizations. We modify the BiCGSTAB algorithm to use a recycle space, which is built from left and right approximate invariant subspaces. Using our algorithm for a parametric model order reduction example gives good results. We show about 40% savings in the number of matrix-vector products and about 35% savings in runtime.

Pawan Goyal, Igor Pontes Duff, Peter Benner

Scientific machine learning for inferring dynamical systems combines data-driven modeling, physics-based modeling, and empirical knowledge. It plays an essential role in engineering design and digital twinning. In this work, we primarily focus on an operator inference methodology that builds dynamical models, preferably in low-dimension, with a prior hypothesis on the model structure, often determined by known physics or given by experts. Then, for inference, we aim to learn the operators of a model by setting up an appropriate optimization problem. One of the critical properties of dynamical systems is stability. However, this property is not guaranteed by the inferred models. In this work, we propose inference formulations to learn quadratic models, which are stable by design. Precisely, we discuss the parameterization of quadratic systems that are locally and globally stable. Moreover, for quadratic systems with no stable point yet bounded (e.g., chaotic Lorenz model), we discuss how to parameterize such bounded behaviors in the learning process. Using those parameterizations, we set up inference problems, which are then solved using a gradient-based optimization method. Furthermore, to avoid numerical derivatives and still learn continuous systems, we make use of an integral form of differential equations. We present several numerical examples, illustrating the preservation of stability and discussing its comparison with the existing state-of-the-art approach to infer operators. By means of numerical examples, we also demonstrate how the proposed methods are employed to discover governing equations and energy-preserving models.

Harshit Kapadia, Lihong Feng, Peter Benner

When repeated evaluations for varying parameter configurations of a high-fidelity physical model are required, surrogate modeling techniques based on model order reduction are desired. In absence of the governing equations describing the dynamics, we need to construct the parametric reduced-order surrogate model in a non-intrusive fashion. In this setting, the usual residual-based error estimate for optimal parameter sampling associated with the reduced basis method is not directly available. Our work provides a non-intrusive optimality criterion to efficiently populate the parameter snapshots, thereby, enabling us to effectively construct a parametric surrogate model. We consider separate parameter-specific proper orthogonal decomposition (POD) subspaces and propose an active-learning-driven surrogate model using kernel-based shallow neural networks, abbreviated as ActLearn-POD-KSNN surrogate model. To demonstrate the validity of our proposed ideas, we present numerical experiments using two physical models, namely Burgers' equation and shallow water equations. Both the models have mixed -- convective and diffusive -- effects within their respective parameter domains, with each of them dominating in certain regions. The proposed ActLearn-POD-KSNN surrogate model efficiently predicts the solution at new parameter locations, even for a setting with multiple interacting shock profiles.

Pawan Goyal, Peter Benner

The engineering design process often relies on mathematical modeling that can describe the underlying dynamic behavior. In this work, we present a data-driven methodology for modeling the dynamics of nonlinear systems. To simplify this task, we aim to identify a coordinate transformation that allows us to represent the dynamics of nonlinear systems using a common, simple model structure. The advantage of a common simple model is that customized design tools developed for it can be applied to study a large variety of nonlinear systems. The simplest common model -- one can think of -- is linear, but linear systems often fall short in accurately capturing the complex dynamics of nonlinear systems. In this work, we propose using quadratic systems as the common structure, inspired by the lifting principle. According to this principle, smooth nonlinear systems can be expressed as quadratic systems in suitable coordinates without approximation errors. However, finding these coordinates solely from data is challenging. Here, we leverage deep learning to identify such lifted coordinates using only data, enabling a quadratic dynamical system to describe the system's dynamics. Additionally, we discuss the asymptotic stability of these quadratic dynamical systems. We illustrate the approach using data collected from various numerical examples, demonstrating its superior performance with the existing well-known techniques.

Xin Du, Peter Benner

This paper discusses model order reduction of LTI systems over limited frequency intervals within the framework of balanced truncation. Two new \emph{frequency-dependent balanced truncation} methods were developed, one is \emph{SF-type frequency-dependent balanced truncation} to copy with the cases that only a single dominating point of the operating frequency interval is pre-known, the other is \emph{interval-type frequency-dependent balanced truncation} to deal with the cases that both of the upper and lower bound of frequency interval are known \emph{a priori}. SF-type error bound and interval-type error bound are derived for the first time to estimate the desired approximation error over pre-specified frequency interval. We show that the new methods generally lead to good in-band approximation performance, at the same time, provide accurate error bounds under certain conditions. Examples are included for illustration.

Süleyman Yildiz, Pawan Goyal, Thomas Bendokat et al.

We present a framework for learning Hamiltonian systems using data. This work is based on a lifting hypothesis, which posits that nonlinear Hamiltonian systems can be written as nonlinear systems with cubic Hamiltonians. By leveraging this, we obtain quadratic dynamics that are Hamiltonian in a transformed coordinate system. To that end, for given generalized position and momentum data, we propose a methodology to learn quadratic dynamical systems, enforcing the Hamiltonian structure in combination with a weakly-enforced symplectic auto-encoder. The obtained Hamiltonian structure exhibits long-term stability of the system, while the cubic Hamiltonian function provides relatively low model complexity. For low-dimensional data, we determine a higher-dimensional transformed coordinate system, whereas for high-dimensional data, we find a lower-dimensional coordinate system with the desired properties. We demonstrate the proposed methodology by means of both low-dimensional and high-dimensional nonlinear Hamiltonian systems.

Peter Benner, Christian Himpe

A standard approach for model reduction of linear input-output systems is balanced truncation, which is based on the controllability and observability properties of the underlying system. The related dominant subspace projection model reduction method similarly utilizes these system properties, yet instead of balancing, the associated subspaces are directly conjoined. In this work we extend the dominant subspace approach by computation via the cross Gramian for linear systems, and describe an a-priori error indicator for this method. Furthermore, efficient computation is discussed alongside numerical examples illustrating these findings.

Peter Benner, Norman Lang

We consider the application of implicit and linearly implicit (Rosenbrock-type) peer methods to matrix-valued ordinary differential equations. In particular the differential Riccati equation (DRE) is investigated. For the Rosenbrock-type schemes, a reformulation capable of avoiding a number of Jacobian applications is developed that, in the autonomous case, reduces the computational complexity of the algorithms. Dealing with large-scale problems, an efficient implementation based on low-rank symmetric indefinite factorizations is presented. The performance of both peer approaches up to order 4 is compared to existing implicit time integration schemes for matrix-valued differential equations.

Petar Mlinarić, Takayuki Ishizaki, Aranya Chakrabortty et al.

We study nonlinear power systems consisting of generators, generator buses, and non-generator buses. First, looking at a generator and its bus' variables jointly, we introduce a synchronization concept for a pair of such joint generators and buses. We show that this concept is related to graph symmetry. Next, we extend, in two ways, the synchronization from a pair to a partition of all generators in the networks and show that they are related to either graph symmetry or equitable partitions. Finally, we show how an exact reduced model can be obtained by aggregating the generators and associated buses in the network when the original system is synchronized with respect to a partition, provided that the initial condition respects the partition. Additionally, the aggregation-based reduced model is again a power system.

Cleophas Kweyu, Lihong Feng, Matthias Stein et al.

The Poisson-Boltzmann equation (PBE) is a nonlinear elliptic PDE that arises in biomolecular modeling and is a fundamental tool for structural biology. It is used to calculate electrostatic potentials around an ensemble of fixed charges immersed in an ionic solution. Efficient numerical computation of the PBE yields a high number of degrees of freedom in the resultant algebraic system of equations. Coupled with the fact that in most cases the PBE requires to be solved multiple times for a large number of system configurations, this poses great computational challenges to conventional numerical techniques. To accelerate such computations, we here present the reduced basis method (RBM) which greatly reduces this computational complexity by constructing a reduced order model of typically low dimension. In this study, we employ a simple version of the PBE for proof of concept and discretize the linearized PBE (LPBE) with a centered finite difference scheme. The resultant linear system is solved by the aggregation-based algebraic multigrid method at different samples of ionic strength on a three-dimensional Cartesian grid. The discretized LPBE, which we call the high-fidelity full order model (FOM), yields solution as accurate as other LPBE solvers. We then apply the RBM to FOM. The discrete empirical interpolation method (DEIM) is applied to the Dirichlet boundary conditions which are nonaffine with the parameter, to reduce the complexity of the reduced order model (ROM). From the numerical results, we notice that the RBM reduces the model order from $\mathcal{N} = 2\times 10^{6}$ to $N = 6$ at an accuracy of $10^{-9}$ and reduces computational time by a factor of approximately $7,600$. DEIM, on the other hand, is also used in the offline-online phase of solving the ROM for different values of parameters which provides a speed-up of $20$ for a single iteration of the greedy algorithm.

Peter Benner, Sergey Dolgov, Akwum Onwunta et al.

Many problems in computational science and engineering are simultaneously characterized by the following challenging issues: uncertainty, nonlinearity, nonstationarity and high dimensionality. Existing numerical techniques for such models would typically require considerable computational and storage resources. This is the case, for instance, for an optimization problem governed by time-dependent Navier-Stokes equations with uncertain inputs. In particular, the stochastic Galerkin finite element method often leads to a prohibitively high dimensional saddle-point system with tensor product structure. In this paper, we approximate the solution by the low-rank Tensor Train decomposition, and present a numerically efficient algorithm to solve the optimality equations directly in the low-rank representation. We show that the solution of the vorticity minimization problem with a distributed control admits a representation with ranks that depend modestly on model and discretization parameters even for high Reynolds numbers. For lower Reynolds numbers this is also the case for a boundary control. This opens the way for a reduced-order modeling of the stochastic optimal flow control with a moderate cost at all stages.

Igor Pontes Duff, Sara Grundel, Peter Benner

In this paper, we propose new algebraic Gramians for continuous-time linear switched systems, which satisfy generalized Lyapunov equations. The main contribution of this work is twofold. First, we show that the ranges of those Gramians encode the reachability and observability spaces of a linear switched system. As a consequence, a simple Gramian-based criterion for reachability and observability is established. Second, a balancing-based model order reduction technique is proposed and, under some sufficient conditions, stability preservation and an error bound are shown. Finally, the efficiency of the proposed method is illustrated by means of numerical examples.

Peter Benner, Venera Khoromskaia, Boris N. Khoromskij

We introduce and analyze the new range-separated (RS) canonical/Tucker tensor format which aims for numerical modeling of the 3D long-range interaction potentials in multi-particle systems. The main idea of the RS tensor format is the independent grid-based low-rank representation of the localized and global parts in the target tensor which allows the efficient numerical approximation of $N$-particle interaction potentials. The single-particle reference potential like $1/\|x\|$ is split into a sum of localized and long-range low-rank canonical tensors represented on a fine 3D $n\times n\times n$ Cartesian grid. The smoothed long-range contribution to the total potential sum is represented on the 3D grid in $O(n)$ storage via the low-rank canonical/Tucker tensor. We prove that the Tucker rank parameters depend only logarithmically on the number of particles $N$ and the grid-size $n$. Agglomeration of the short range part in the sum is reduced to an independent treatment of $N$ localized terms with almost disjoint effective supports, calculated in $O(N)$ operations. Thus, the cumulated sum of short range clusters is parametrized by a single low-rank canonical reference tensor with a local support, accomplished by a list of particle coordinates and their charges. The RS canonical/Tucker tensor representations reduce the cost of multi-linear algebraic operations on the 3D potential sums arising in modeling of multi-dimensional data by radial basis functions, say, in computation of the electrostatic potential of a protein, in 3D integration and convolution transforms, computation of gradients, forces and the interaction energy of a many-particle systems, and in low parametric fitting of multi-dimensional scattered data by reducing all of them to 1D calculations.

Pawan Goyal, Peter Benner

Measurement noise is an integral part while collecting data of a physical process. Thus, noise removal is necessary to draw conclusions from these data, and it often becomes essential to construct dynamical models using these data. We discuss a methodology to learn differential equation(s) using noisy and irregular sampled measurements. In our methodology, the main innovation can be seen in the integration of deep neural networks with the neural ordinary differential equations (ODEs) approach. Precisely, we aim at learning a neural network that provides (approximately) an implicit representation of the data and an additional neural network that models the vector fields of the dependent variables. We combine these two networks by constraining using neural ODEs. The proposed framework to learn a model describing the vector field is highly effective under noisy measurements. The approach can handle scenarios where dependent variables are not available at the same temporal grid. Moreover, a particular structure, e.g., second-order with respect to time, can easily be incorporated. We demonstrate the effectiveness of the proposed method for learning models using data obtained from various differential equations and present a comparison with the neural ODE method that does not make any special treatment to noise.

Pawan Goyal, Süleyman Yıldız, Peter Benner

Discovering a suitable coordinate transformation for nonlinear systems enables the construction of simpler models, facilitating prediction, control, and optimization for complex nonlinear systems. To that end, Koopman operator theory offers a framework for global linearization for nonlinear systems, thereby allowing the usage of linear tools for design studies. In this work, we focus on the identification of global linearized embeddings for canonical nonlinear Hamiltonian systems through a symplectic transformation. While this task is often challenging, we leverage the power of deep learning to discover the desired embeddings. Furthermore, to overcome the shortcomings of Koopman operators for systems with continuous spectra, we apply the lifting principle and learn global cubicized embeddings. Additionally, a key emphasis is paid to enforce the bounded stability for the dynamics of the discovered embeddings. We demonstrate the capabilities of deep learning in acquiring compact symplectic coordinate transformation and the corresponding simple dynamical models, fostering data-driven learning of nonlinear canonical Hamiltonian systems, even those with continuous spectra.

Pawan Goyal, Igor Pontes Duff, Peter Benner

Machine-learning technologies for learning dynamical systems from data play an important role in engineering design. This research focuses on learning continuous linear models from data. Stability, a key feature of dynamic systems, is especially important in design tasks such as prediction and control. Thus, there is a need to develop methodologies that provide stability guarantees. To that end, we leverage the parameterization of stable matrices proposed in [Gillis/Sharma, Automatica, 2017] to realize the desired models. Furthermore, to avoid the estimation of derivative information to learn continuous systems, we formulate the inference problem in an integral form. We also discuss a few extensions, including those related to control systems. Numerical experiments show that the combination of a stable matrix parameterization and an integral form of differential equations allows us to learn stable systems without requiring derivative information, which can be challenging to obtain in situations with noisy or limited data.

Peter Benner, Pawan Goyal

We discuss model reduction for a particular class of quadratic-bilinear (QB) descriptor systems. The main goal of this article is to extend the recently studied interpolation-based optimal model reduction framework for QBODEs [Benner et al. '16] to a class of descriptor systems in an efficient and reliable way. Recently, it has been shown in the case of linear or bilinear systems that a direct extension of interpolation-based model reduction techniques to descriptor systems, without any modifications, may lead to poor reduced-order systems. Therefore, for the analysis, we aim at transforming the considered QB descriptor system into an equivalent QBODE system by means of projectors for which standard model reduction techniques for QBODEs can be employed, including aforementioned interpolation scheme. Subsequently, we discuss related computational issues, thus resulting in a modified algorithm that allows us to construct \emph{near}--optimal reduced-order systems without explicitly computing the projectors used in the analysis. The efficiency of the proposed algorithm is illustrated by means of a numerical example, obtained via semi-discretization of the Navier-Stokes equations.

Xin Du, Peter Benner

Balanced truncation is one of the most common model order reduction schemes. In this paper, we study finite-frequency model order reduction (FF-MOR) problems of linear continuous-time systems within the framework of balanced truncation method. Firstly, we construct a family of parameterized frequency-dependent (PFD) mappings which generate discrete-time PFD mapped systems and continuous-time PFD mapped systems of the given continuous-time system. The relationships between the maximum singular value of the given system over pre-specified frequency ranges and the maximum singular value of the PFD mapped systems over entire frequency range are established. By exploiting the properties of the discrete-time PFD mapped systems, a new parameterized frequency-dependent balanced truncation (PFDBT) method providing finite-frequency type error bound with respect to the maximum singular value of the approximation error systems are developed. Examples are included for illustration.

Tobias Damm, Peter Benner, Jan Hauth

The stochastic $H^\infty$-norm is defined as the $L^2$-induced norm of the input-output operator of a stochastic linear system. Like the deterministic $H^\infty$-norm it is characterised by a version of the bounded real lemma, but without a frequency domain description or a Hamiltonian condition. Therefore, we base its computation on a parametrised algebraic Riccati-type matrix equation.

Kirandeep Kour, Sergey Dolgov, Peter Benner et al.

High-dimensional data in the form of tensors are challenging for kernel classification methods. To both reduce the computational complexity and extract informative features, kernels based on low-rank tensor decompositions have been proposed. However, what decisive features of the tensors are exploited by these kernels is often unclear. In this paper we propose a novel kernel that is based on the Tucker decomposition. For this kernel the Tucker factors are computed based on re-weighting of the Tucker matrices with tuneable powers of singular values from the HOSVD decomposition. This provides a mechanism to balance the contribution of the Tucker core and factors of the data. We benchmark support tensor machines with this new kernel on several datasets. First we generate synthetic data where two classes differ in either Tucker factors or core, and compare our novel and previously existing kernels. We show robustness of the new kernel with respect to both classification scenarios. We further test the new method on real-world datasets. The proposed kernel has demonstrated a higher test accuracy than the state-of-the-art tensor train multi-way multi-level kernel, and a significantly lower computational time.

Peter Benner, Xiaobo Liu

We develop a mixed-precision iterative refinement framework for solving low-rank Lyapunov matrix equations $AX + XA^T + W =0$, where $W=LL^T$ or $W=LSL^T$. Via rounding error analysis of the algorithms we derive sufficient conditions for the attainable normwise residuals in different precision settings and show how the algorithmic parameters should be chosen. These conditions are independent of the choice of inner solver, provided that the prescribed residual accuracy is attained in the inner solves. Using the sign-function Newton iteration as the solver, we demonstrate that reduced precisions, such as half precision with unit roundoff $u_s$, can be used efficiently for Lyapunov equations with condition numbers of order $1/u_s$ without compromising the attainable solution quality. This provides an algorithmic framework towards exploiting native low-precision hardware to accelerate Lyapunov solvers without sacrificing accuracy.

Peter Benner, Boris N. Khoromskij, Venera Khoromskaia et al.

We introduce and analyze a mesh-free two-level hybrid Chebyshev-Tucker tensor representation for approximating multivariate functions, which combines tensor-product Chebyshev interpolation with the low-rank Tucker decomposition of the tensor of Chebyshev coefficients. This construction allows to avoid the expensive rank-structured grid-based approximation of function-related tensors on large spatial grids, while benefiting from the Tucker decomposition of the moderate-sized core tensor of Chebyshev coefficients. Thus, we can compute the nearly optimal Tucker decomposition of the 3D function with controllable accuracy $\varepsilon >0$ without discretizing the function on a full fine grid in the domain, but only using its values at a small set of Chebyshev nodes computed either from the explicit analytic expression of the target function or from its data-sparse representation in a rank-structured tensor format with moderate rank parameter. Finally, we can represent the function in the algebraic Tucker format with optimal $\varepsilon$-rank on an arbitrarily large 3D tensor grid in the computational domain by discretizing the Chebyshev polynomials on that grid. The rank parameters of the nonlinear Tucker-ALS decomposition of the coefficient tensor can be much smaller than the polynomial degrees of the initial Chebyshev linear interpolation in the function independent polynomial basis set. It is shown that our techniques can be gainfully applied to the long-range part of the singular electrostatic potential of multi-particle systems represented on a fine grid in the range-separated (RS) tensor format. We provide error and complexity estimates and demonstrate the computational efficiency of the proposed techniques on challenging examples, including the collective electrostatic potential for large bio-molecular systems and lattice-type compounds.

Ali Forootani, Pawan Goyal, Peter Benner

The discovery of governing equations from data has been an active field of research for decades. One widely used methodology for this purpose is sparse regression for nonlinear dynamics, known as SINDy. Despite several attempts, noisy and scarce data still pose a severe challenge to the success of the SINDy approach. In this work, we discuss a robust method to discover nonlinear governing equations from noisy and scarce data. To do this, we make use of neural networks to learn an implicit representation based on measurement data so that not only it produces the output in the vicinity of the measurements but also the time-evolution of output can be described by a dynamical system. Additionally, we learn such a dynamic system in the spirit of the SINDy framework. Leveraging the implicit representation using neural networks, we obtain the derivative information -- required for SINDy -- using an automatic differentiation tool. To enhance the robustness of our methodology, we further incorporate an integral condition on the output of the implicit networks. Furthermore, we extend our methodology to handle data collected from multiple initial conditions. We demonstrate the efficiency of the proposed methodology to discover governing equations under noisy and scarce data regimes by means of several examples and compare its performance with existing methods.

Mohammad S. Khorrami, Pawan Goyal, Jaber R. Mianroodi et al.

The purpose of the current work is the development of a so-called physics-encoded Fourier neural operator (PeFNO) for surrogate modeling of the quasi-static equilibrium stress field in solids. Rather than accounting for constraints from physics in the loss function as done in the (now standard) physics-informed approach, the physics-encoded approach incorporates or "encodes" such constraints directly into the network or operator architecture. As a result, in contrast to the physics-informed approach in which only training is physically constrained, both training and output are physically constrained in the physics-encoded approach. For the current constraint of divergence-free stress, a novel encoding approach based on a stress potential is proposed. As a "proof-of-concept" example application of the proposed PeFNO, a heterogeneous polycrystalline material consisting of isotropic elastic grains subject to uniaxial extension is considered. Stress field data for training are obtained from the numerical solution of a corresponding boundary-value problem for quasi-static mechanical equilibrium. This data is also employed to train an analogous physics-guided FNO (PgFNO) and physics-informed FNO (PiFNO) for comparison. As confirmed by this comparison and as expected on the basis of their differences, the output of the trained PeFNO is significantly more accurate in satisfying mechanical equilibrium than the output of either the trained PgFNO or the trained PiFNO.

Süleyman Yıldız, Pawan Goyal, Peter Benner

This paper presents an energy-preserving machine learning method for inferring reduced-order models (ROMs) by exploiting the multi-symplectic form of partial differential equations (PDEs). The vast majority of energy-preserving reduced-order methods use symplectic Galerkin projection to construct reduced-order Hamiltonian models by projecting the full models onto a symplectic subspace. However, symplectic projection requires the existence of fully discrete operators, and in many cases, such as black-box PDE solvers, these operators are inaccessible. In this work, we propose an energy-preserving machine learning method that can infer the dynamics of the given PDE using data only, so that the proposed framework does not depend on the fully discrete operators. In this context, the proposed method is non-intrusive. The proposed method is grey box in the sense that it requires only some basic knowledge of the multi-symplectic model at the partial differential equation level. We prove that the proposed method satisfies spatially discrete local energy conservation and preserves the multi-symplectic conservation laws. We test our method on the linear wave equation, the Korteweg-de Vries equation, and the Zakharov-Kuznetsov equation. We test the generalization of our learned models by testing them far outside the training time interval.

Celine Reddig, Pawan Goyal, Igor Pontes Duff et al.

Model reduction is an active research field to construct low-dimensional surrogate models of high fidelity to accelerate engineering design cycles. In this work, we investigate model reduction for linear structured systems using dominant reachable and observable subspaces. When the training set $-$ containing all possible interpolation points $-$ is large, then these subspaces can be determined by solving many large-scale linear systems. However, for high-fidelity models, this easily becomes computationally intractable. To circumvent this issue, in this work, we propose an active sampling strategy to sample only a few points from the given training set, which can allow us to estimate those subspaces accurately. To this end, we formulate the identification of the subspaces as the solution of the generalized Sylvester equations, guiding us to select the most relevant samples from the training set to achieve our goals. Consequently, we construct solutions of the matrix equations in low-rank forms, which encode subspace information. We extensively discuss computational aspects and efficient usage of the low-rank factors in the process of obtaining reduced-order models. We illustrate the proposed active sampling scheme to obtain reduced-order models via dominant reachable and observable subspaces and present its comparison with the method where all the points from the training set are taken into account. It is shown that the active sample strategy can provide us $17$x speed-up without sacrificing any noticeable accuracy.

Martin Köhler, Peter Benner

Although Fortran is almost 70 years old, the language continues to evolve in order to keep pace with developments in computer science. In particular, a flexible type system was introduced that allows developers to specify the sizes of floating-point numbers and integers. In the latest revisions of the Fortran standard, portable type variants for IEEE 754 binary64 (double precision, real64) and binary32 (single precision, real32) were added. However, the rapid development of AI toolkits and accelerator hardware has created a strong focus on floating-point types of lower precision and lower memory usage than binary32. While the IEEE 754-2019 standard defines the binary16 type for representing half-precision numbers, the Fortran standard does not provide the real16 variant in the type system. In contrast, most C compilers support such a data type. In numerical linear algebra, there is strong interest in exploiting the high performance of accelerator devices for core algorithms like matrix decompositions or iterative solvers. Especially when the performance ratio between double, single, and half precision is on the order of 1:2:20, as on current NVidia H100 accelerators, it becomes highly beneficial to use lower-precision types. Yet, before performance can be targeted, correctness and accuracy must be verified when operating below single precision. In this article, we present our Low Precision Fortran (LPF) library that enables the use of low-precision types -- binary16, bfloat16, fp8_e4m3, and fp8_e5m2 -- just like any other floating-point type in Fortran. Furthermore, we introduce extensions that support BLAS operations in low precision and show how easily existing routines can be rewritten to use these data types.

Igor Pontes Duff, Pawan Goyal, Peter Benner

This work primarily focuses on an operator inference methodology aimed at constructing low-dimensional dynamical models based on a priori hypotheses about their structure, often informed by established physics or expert insights. Stability is a fundamental attribute of dynamical systems, yet it is not always assured in models derived through inference. Our main objective is to develop a method that facilitates the inference of quadratic control dynamical systems with inherent stability guarantees. To this aim, we investigate the stability characteristics of control systems with energy-preserving nonlinearities, thereby identifying conditions under which such systems are bounded-input bounded-state stable. These insights are subsequently applied to the learning process, yielding inferred models that are inherently stable by design. The efficacy of our proposed framework is demonstrated through a couple of numerical examples.

Konrad Janik, Peter Benner

Measurement data is often sampled irregularly i.e. not on equidistant time grids. This is also true for Hamiltonian systems. However, existing machine learning methods, which learn symplectic integrators, such as SympNets [20] and HénonNets [4] still require training data generated by fixed step sizes. To learn time-adaptive symplectic integrators, an extension to SympNets, which we call TSympNets, was introduced in [20]. We adapt the architecture of TSympNets and extend them to non-autonomous Hamiltonian systems. So far the approximation qualities of TSympNets were unknown. We close this gap by providing a universal approximation theorem for separable Hamiltonian systems and show that it is not possible to extend it to non-separable Hamiltonian systems. To investigate these theoretical approximation capabilities, we perform different numerical experiments. Furthermore we fix a mistake in a proof of a substantial theorem [25, Theorem 2] for the approximation of symplectic maps in general, but specifically for symplectic machine learning methods.

Shuwen Sun, Lihong Feng, Peter Benner

Numerically solving a large parametric nonlinear dynamical system is challenging due to its high complexity and the high computational costs. In recent years, machine-learning-aided surrogates are being actively researched. However, many methods fail in accurately generalizing in the entire time interval $[0, T]$, when the training data is available only in a training time interval $[0, T_0]$, with $T_0<T$. To improve the extrapolation capabilities of the surrogate models in the entire time domain, we propose a new deep learning framework, where kernel dynamic mode decomposition (KDMD) is employed to evolve the dynamics of the latent space generated by the encoder part of a convolutional autoencoder (CAE). After adding the KDMD-decoder-extrapolated data into the original data set, we train the CAE along with a feed-forward deep neural network using the augmented data. The trained network can predict future states outside the training time interval at any out-of-training parameter samples. The proposed method is tested on two numerical examples: a FitzHugh-Nagumo model and a model of incompressible flow past a cylinder. Numerical results show accurate and fast prediction performance in both the time and the parameter domain.

Ion Victor Gosea, Luisa Peterson, Pawan Goyal et al.

In this work, we address the challenge of efficiently modeling dynamical systems in process engineering. We use reduced-order model learning, specifically operator inference. This is a non-intrusive, data-driven method for learning dynamical systems from time-domain data. The application in our study is carbon dioxide methanation, an important reaction within the Power-to-X framework, to demonstrate its potential. The numerical results show the ability of the reduced-order models constructed with operator inference to provide a reduced yet accurate surrogate solution. This represents an important milestone towards the implementation of fast and reliable digital twin architectures.

Pawan Goyal, Peter Benner

Learning dynamical models from data plays a vital role in engineering design, optimization, and predictions. Building models describing dynamics of complex processes (e.g., weather dynamics, or reactive flows) using empirical knowledge or first principles are onerous or infeasible. Moreover, these models are high-dimensional but spatially correlated. It is, however, observed that the dynamics of high-fidelity models often evolve in low-dimensional manifolds. Furthermore, it is also known that for sufficiently smooth vector fields defining the nonlinear dynamics, a quadratic model can describe it accurately in an appropriate coordinate system, conferring to the McCormick relaxation idea in nonconvex optimization. Here, we aim at finding a low-dimensional embedding of high-fidelity dynamical data, ensuring a simple quadratic model to explain its dynamics. To that aim, this work leverages deep learning to identify low-dimensional quadratic embeddings for high-fidelity dynamical systems. Precisely, we identify the embedding of data using an autoencoder to have the desired property of the embedding. We also embed a Runge-Kutta method to avoid the time-derivative computations, which is often a challenge. We illustrate the ability of the approach by a couple of examples, arising in describing flow dynamics and the oscillatory tubular reactor model.

Pawan Goyal, Peter Benner

Measurement noise is an integral part while collecting data of a physical process. Thus, noise removal is a necessary step to draw conclusions from these data, and it often becomes quite essential to construct dynamical models using these data. We discuss a methodology to learn differential equation(s) using noisy and sparsely sampled measurements. In our methodology, the main innovation can be seen in of integration of deep neural networks with a classical numerical integration method. Precisely, we aim at learning a neural network that implicitly represents the data and an additional neural network that models the vector fields of the dependent variables. We combine these two networks by enforcing the constraint that the data at the next time-steps can be given by following a numerical integration scheme such as the fourth-order Runge-Kutta scheme. The proposed framework to learn a model predicting the vector field is highly effective under noisy measurements. The approach can handle scenarios where dependent variables are not available at the same temporal grid. We demonstrate the effectiveness of the proposed method to learning models using data obtained from various differential equations. The proposed approach provides a promising methodology to learn dynamic models, where the first-principle understanding remains opaque.

Pawan Goyal, Peter Benner

Discovering dynamical models to describe underlying dynamical behavior is essential to draw decisive conclusions and engineering studies, e.g., optimizing a process. Experimental data availability notwithstanding has increased significantly, but interpretable and explainable models in science and engineering yet remain incomprehensible. In this work, we blend machine learning and dictionary-based learning with numerical analysis tools to discover governing differential equations from noisy and sparsely-sampled measurement data. We utilize the fact that given a dictionary containing huge candidate nonlinear functions, dynamical models can often be described by a few appropriately chosen candidates. As a result, we obtain interpretable and parsimonious models which are prone to generalize better beyond the sampling regime. Additionally, we integrate a numerical integration framework with dictionary learning that yields differential equations without requiring or approximating derivative information at any stage. Hence, it is utterly effective in corrupted and sparsely-sampled data. We discuss its extension to governing equations, containing rational nonlinearities that typically appear in biological networks. Moreover, we generalized the method to governing equations that are subject to parameter variations and externally controlled inputs. We demonstrate the efficiency of the method to discover a number of diverse differential equations using noisy measurements, including a model describing neural dynamics, chaotic Lorenz model, Michaelis-Menten Kinetics, and a parameterized Hopf normal form.

Pawan Goyal, Peter Benner

Mathematical modeling is an essential step, for example, to analyze the transient behavior of a dynamical process and to perform engineering studies such as optimization and control. With the help of first-principles and expert knowledge, a dynamic model can be built, but for complex dynamic processes, appearing, e.g., in biology, chemical plants, neuroscience, financial markets, this often remains an onerous task. Hence, data-driven modeling of the dynamics process becomes an attractive choice and is supported by the rapid advancement in sensor and measurement technology. A data-driven approach, namely operator inference framework, models a dynamic process, where a particular structure of the nonlinear term is assumed. In this work, we suggest combining the operator inference with certain deep neural network approaches to infer the unknown nonlinear dynamics of the system. The approach uses recent advancements in deep learning and possible prior knowledge of the process if possible. We also briefly discuss several extensions and advantages of the proposed methodology. We demonstrate that the proposed methodology accomplishes the desired tasks for dynamics processes encountered in neural dynamics and the glycolytic oscillator.

Peter Benner, Pawan Goyal, Jan Heiland et al.

Reduced-order modeling has a long tradition in computational fluid dynamics. The ever-increasing significance of data for the synthesis of low-order models is well reflected in the recent successes of data-driven approaches such as Dynamic Mode Decomposition and Operator Inference. With this work, we suggest a new approach to learning structured low-order models for incompressible flow from data that can be used for engineering studies such as control, optimization, and simulation. To that end, we utilize the intrinsic structure of the Navier-Stokes equations for incompressible flows and show that learning dynamics of the velocity and pressure can be decoupled, thus leading to an efficient operator inference approach for learning the underlying dynamics of incompressible flows. Furthermore, we show the operator inference performance in learning low-order models using two benchmark problems and compare with an intrusive method, namely proper orthogonal decomposition, and other data-driven approaches.

Anke Stoll, Peter Benner

Currently, the growth of material data from experiments and simulations is expanding beyond processable amounts. This makes the development of new data-driven methods for the discovery of patterns among multiple lengthscales and time-scales and structure-property relationships essential. These data-driven approaches show enormous promise within materials science. The following review covers machine learning applications for metallic material characterization. Many parameters associated with the processing and the structure of materials affect the properties and the performance of manufactured components. Thus, this study is an attempt to investigate the usefulness of machine learning methods for material property prediction. Material characteristics such as strength, toughness, hardness, brittleness or ductility are relevant to categorize a material or component according to their quality. In industry, material tests like tensile tests, compression tests or creep tests are often time consuming and expensive to perform. Therefore, the application of machine learning approaches is considered helpful for an easier generation of material property information. This study also gives an application of machine learning methods on small punch test data for the determination of the property ultimate tensile strength for various materials. A strong correlation between small punch test data and tensile test data was found which ultimately allows to replace more costly tests by simple and fast tests in combination with machine learning.

Pawan Goyal, Hussam Al Daas, Peter Benner

In this paper, we study the problem of image recovery from given partial (corrupted) observations. Recovering an image using a low-rank model has been an active research area in data analysis and machine learning. But often, images are not only of low-rank but they also exhibit sparsity in a transformed space. In this work, we propose a new problem formulation in such a way that we seek to recover an image that is of low-rank and has sparsity in a transformed domain. We further discuss various non-convex non-smooth surrogates of the rank function, leading to a relaxed problem. Then, we present an efficient iterative scheme to solve the relaxed problem that essentially employs the (weighted) singular value thresholding at each iteration. Furthermore, we discuss the convergence properties of the proposed iterative method. We perform extensive experiments, showing that the proposed algorithm outperforms state-of-the-art methodologies in recovering images.

Peter Benner, Pawan Goyal, Boris Kramer et al.

This work presents a non-intrusive model reduction method to learn low-dimensional models of dynamical systems with non-polynomial nonlinear terms that are spatially local and that are given in analytic form. In contrast to state-of-the-art model reduction methods that are intrusive and thus require full knowledge of the governing equations and the operators of a full model of the discretized dynamical system, the proposed approach requires only the non-polynomial terms in analytic form and learns the rest of the dynamics from snapshots computed with a potentially black-box full-model solver. The proposed method learns operators for the linear and polynomially nonlinear dynamics via a least-squares problem, where the given non-polynomial terms are incorporated in the right-hand side. The least-squares problem is linear and thus can be solved efficiently in practice. The proposed method is demonstrated on three problems governed by partial differential equations, namely the diffusion-reaction Chafee-Infante model, a tubular reactor model for reactive flows, and a batch-chromatography model that describes a chemical separation process. The numerical results provide evidence that the proposed approach learns reduced models that achieve comparable accuracy as models constructed with state-of-the-art intrusive model reduction methods that require full knowledge of the governing equations.

Kirandeep Kour, Sergey Dolgov, Martin Stoll et al.

An increasing amount of collected data are high-dimensional multi-way arrays (tensors), and it is crucial for efficient learning algorithms to exploit this tensorial structure as much as possible. The ever-present curse of dimensionality for high dimensional data and the loss of structure when vectorizing the data motivates the use of tailored low-rank tensor classification methods. In the presence of small amounts of training data, kernel methods offer an attractive choice as they provide the possibility for a nonlinear decision boundary. We develop the Tensor Train Multi-way Multi-level Kernel (TT-MMK), which combines the simplicity of the Canonical Polyadic decomposition, the classification power of the Dual Structure-preserving Support Vector Machine, and the reliability of the Tensor Train (TT) approximation. We show by experiments that the TT-MMK method is usually more reliable computationally, less sensitive to tuning parameters, and gives higher prediction accuracy in the SVM classification when benchmarked against other state-of-the-art techniques.

Nicat Aliyev, Peter Benner, Emre Mengi et al.

We deal with the minimization of the ${\mathcal H}_\infty$-norm of the transfer function of a parameter-dependent descriptor system over the set of admissible parameter values. Subspace frameworks are proposed for such minimization problems where the involved systems are of large order. The proposed algorithms are greedy interpolatory approaches inspired by our recent work [Aliyev et al., SIAM J. Matrix Anal. Appl., 38(4):1496--1516, 2017] for the computation of the ${\mathcal H}_\infty$-norm. In this work, we minimize the ${\mathcal H}_\infty$-norm of a reduced-order parameter-dependent system obtained by two-sided restrictions onto certain subspaces. Then we expand the subspaces so that Hermite interpolation properties hold between the full and reduced-order system at the optimal parameter value for the reduced order system. We formally establish the superlinear convergence of the subspace frameworks under some smoothness assumptions. The fast convergence of the proposed frameworks in practice is illustrated by several large-scale systems.

Peter Benner, Pawan Goyal

In this work, we investigate a model order reduction scheme for polynomial parametric systems. We begin with defining the generalized multivariate transfer functions for the system. Based on this, we aim at constructing a reduced-order system, interpolating the defined generalized transfer functions at a given set of interpolation points. Furthermore, we provide a method, inspired by the Loewner approach for linear and (quadratic-)bilinear systems, to determine a good-quality reduced-order system in an automatic way. We also discuss the computational issues related to the proposed method and a potential application of CUR matrix approximation in order to further speed-up simulations of reduced-order systems. We test the efficiency of the proposed methods via several numerical examples.

Nicat Aliyev, Peter Benner, Emre Mengi et al.

We are concerned with the computation of the ${\mathcal L}_\infty$-norm for an ${\mathcal L}_\infty$-function of the form $H(s) = C(s) D(s)^{-1} B(s)$, where the middle factor is the inverse of a meromorphic matrix-valued function, and $C(s),\, B(s)$ are meromorphic functions mapping to short-and-fat and tall-and-skinny matrices, respectively. For instance, transfer functions of descriptor systems and delay systems fall into this family. We focus on the case where the middle factor is large-scale. We propose a subspace projection method to obtain approximations of the function $H$ where the middle factor is of much smaller dimension. The ${\mathcal L}_\infty$-norms are computed for the resulting reduced functions, then the subspaces are refined by means of the optimal points on the imaginary axis where the ${\mathcal L}_\infty$-norm of the reduced function is attained. The subspace method is designed so that certain Hermite interpolation properties hold between the largest singular values of the original and reduced functions. This leads to a locally superlinearly convergent algorithm with respect to the subspace dimension, which we prove and illustrate on various numerical examples.

Peter Benner, Tobias Damm, Yolanda Rocio Rodriguez Cruz

We consider two approaches to balanced truncation of stochastic linear systems, which follow from different generalizations of the reachability Gramian of deterministic systems. Both preserve mean-square asymptotic stability, but only the second leads to a stochastic H-infinity-type bound for the approximation error of the truncated system.