Half-Inverse Gradients for Physical Deep Learning
Authors: Patrick Schnell
Philipp Holl
Nils Thuerey
Venue: ICLR 2022 Spotlight (2022)
Keywords: physical simulation
partial differential equations
physical loss functions
optimization

Abstract: Recent works in deep learning have shown that integrating differentiable physics simulators into the training process can greatly improve the quality of results. Although this combination represents a more complex optimization task than usual neural network training, the same gradient-based optimizers are used to minimize the loss function. However, the integrated physics solvers have a profound effect on the gradient flow as manipulating scales in magnitude and direction is an inherent property of many physical processes. Consequently, the gradient flow is often highly unbalanced and creates an environment in which existing gradient-based optimizers perform poorly. In this work, we analyze the characteristics of both physical and neural network optimizations separately to derive a new method based on a half-inversion of the Jacobian. Our approach combines principles of both classical network and physics optimizers to solve the combined optimization task. Compared to state-of-the-art neural network optimizers, our method converges more quickly and to better solutions, which we demonstrate on three complex learning problems involving nonlinear oscillators, the Schroedinger equation and the Poisson problem.

Published as a conference paper at ICLR 2022
HALF-INVERSE GRADIENTS FOR PHYSICAL DEEP LEARNING
Patrick Schnell, Philipp Holl and Nils Thuerey
Department of Informatics Technical University of Munich
Bolzzmannstr. 3, 85748 Garching, Germany
{patrick. Patrick.schnelll scwnelhh, pelll,phiipppe how., nils. mtw.oleeeestttteeeetetttt.t ttue
ABSTRACT
Recent works in deep learning have shown that integrating differentiable physics simulators into the training process greatly improve the quality of results. Although this combination represents a more complex optimization task than supervised neural network training, the same gradient-based optimizers are typically employed to minimize the loss function. However, the integrated physics solvers have a profound effect the gradient flow as manipulating scales in magnitude and direction is an inherent property of many physical processes. Consequently, the gradient flow is often highly unbalanced and creates an enviromment in which existing gradient-based optimizers perform poorly. In this work, we analyze the characteristics of both physical and neural network optimizations to derive a new method that does not suffer from this phenomenon. Our method is based on half-inversion of the Jacobian and combines principles of both classical network and physics optimizers solve the combined optimization task. Compared to state-of-the-art neural network optimizers, our method converges more quickly and yields better solutions, which we demonstrate on three complex learning problems involving nonlinear oscillators, the Schrödinger equation and the Poisson problem.
INTRODUCTION
The groundbreaking successes of deep learning (Krizhevsky et al. 2012 Sutskever et al., 2014 Silver et al., 2017) have led to ongoing efforts to study the capabilities of neural networks across all scientilic disciplines. In the area of physical simulation, neural networks have been used in various ways, such as creating accurate aeduced-order models (Morton et al., 2018), inferring improved discretization stencils (Bar-Sinai al. 2019), suppressing numerical errors (Um et al. 2020). The long-term goal of these methods is to exceed classical simulations in terms of accuracy and speed, which has been achieved, e.g., for rigid bodies (de Avila Belbute-Peres et al. (2018), physical inverse problems (Holl et al. 2020), and two-dimensional turbulence (Kochkov et al. 2021).
The successful application of deep learning to physical systems naturally hinges on the training setup. In recent years, the use of physical loss functions has proven beneficial for the training procedure, yielding substantial improvements over purely supervised training approaches (Tompson) et al. 2017, Wu & Tegmark [2019] Greydanus et al. 2019). These improvements were shown to stem from three aspects (Battaglia al. [216 Holl et 2020): (i) Incorporating prior knowledge from physical principles facilitates the learning process (ii) the ambiguities of multimodal cases are resolved naturally, and (iii) simulating the physics at training time can provide more realistic data distributions than pre-computed data sets. Approaches for training with physical losses can be divided into two categories. On the one hand, cquation-focused approaches that introduce physical residuals (Tompson al. 2017 Raissi et al. 2019), and on the other hand, solver-focused approaches that additionally integrate well wetl-established numerical procedures into training (Um et al, 2020 Kochkov et al. 2021).
From a mathematical point of view, training a neural network with a physical loss function bears the difficulties of both network training and physics optimization. In order to obtain satisfying

---

Published as a conference paper at ICLR 2022
results, it is vital to treat flat regions of the optimization landscapes effectively. In learning, the challenging loss landscapes are addressed using gradient- optimizers with data-based normalizing schemes, such as Adam (Kingma & Ba 2015), whereas in physies. the optimizers of choice are higher-order techniques, such as Newton's 's method (Gill & Murray 1978), which inherently make use of inversion processes. However, Holl et al. (2021) found that these approaches can not effectively handle the joint optimization of network and physics.
optimizers suffer from vanishing or exploding gradients, preventing effective convergence, while higher-order methods do not generally scale to the high-dimensional parameter spaces required by deep learning (Goodfellow et al., 2016).
Inspired by the insight that inversion is crucial for physics problems in learning from Holl et al. (202)) we focus on an inversion-based approach but propose a new method for joint physics and network optimization which we refer to as half-inverse gradients. At its core lies a partial matrix inversion, which we derive from the interaction between network and physies both formally and geometrically, An important property of our method is that its runtime scales linearly with the number of network parameters. To demonstrate the wide-ranging and practical applicability of our method, show that it yields significant improvements in terms of convergence speed and final loss values over existing methods. These improvements are measured both in terms of absolute accuracy well as wall-clock time. We evaluate a diverse set of physical systems, such the Schrödinger equation, a nonlinear chain system and the Poisson problem.
2 GRADIENTS BASED ON HALF-INVERSE JACOBIANS
Optimization on continuous spaces can be effectively performed with derivative-basod methods, the simplest of which is gradient descent. For target function L(0) to be minimized of several variables 0, using bold symbols for vector-valued quantitics in this section, and learning rate n,. gradient descent proceeds by repeatedly applying updates
ATC - (ar. 0T
00
For quadratic objectives, this algorithm convergences linearly with the rate of convergence depending on the condition number À of the Hessian matrix (Lax) 2014). In the ill-conditioned case À >> 1, flat regions in the optimization landscape can significantly slow down the optimization progress. This is ubiquitous problem in non-convex optimization tasks of the generic form:
L(0) = Σ (v. 0:) = E (f(x); 9))
Here (x,, 9.) denotes the ith data points from chosen set of measurements, f is function parametrized by 0 to be optimized to model the relationship between the data points y; (0) = [# 0), and 7 denotes loss function measuring the optimization progress. In the following, assume the most common case of lyy, ÿt) _ lly; ÿ: 1l being the squared 1.2-loss.
Physics Optimization. Simulating a physical system consists of two steps: (i) mathematically modeling the system by a differential equation, and (ii) discretizing its differential operators to obtain solver for computer. Optimization tasks occur for instance when manipulating physical system through an external force to reach a given configuration, for which we have to solve an inverse problem of form 2 In such control task, the sum reduces to single data point (æx, ÿ) with æ being the initial state, ý the target state and 0 the external force we want to find. The physical solver corresponds to the function f representing time evolution y(0) f fa; 0). This single data point sum still includes summation over vector components of y in the L2-loss. Sensitive behavior of the physical system arising from its high-frequency modes is present in the physical solver f. and produces small singular values in its Jacobian. This leads to an ill-conditioned Jacobian and flat regions in the optimization landscape when minimizing 2 This is addressed by using methods that incorporate more information than only the gradient. Prominent examples are Newton's method or the latter one is based on the Jacobian of f
the Gauss-Newton 's algorithm (Gill & Murray 1978): and the loss gradient:
AOCx = -0-)) (yy
-1
AL
(3)
00
Ty
2

---

Published as a conference paper at ICLR 2022
Here the inversion of the Jacobian is calculated with the pseudoinverse. The Gauss-Newtes update maps the steepest descent direction in y-space to the parameter space 0. Therefore, to first order, the resulting update approximates gradient descent steps in y-space, further details are given in appendix A.2 An advantage of such higher-order methods is that the update steps in y-space are invariant under arbitrary rescaling of the parameters 0, which cancels inherent scales in f and ensures quick progress in the optimization landscape.
Neural Network Training. For f representing a neural network in equation 2 the optimization matches the typical supervised learning task. In this context, the problem of flat regions in the optimization landscape is also referred as pathological curvature (Martens] 2010). Solving this problem with higher-order methods is considered to be too expensive given the large number of parameters O. For learning tasks, popular optimizers, such as Adam, instead use gradient information from earlier update steps, for instance in the form of momentum or adaptive learning rate terms, thereby improving convergence speed at little additional computational cost. Furthermore, the updates are computed mini-batches instead of the full data set, which saves computational resources and benefits generalization (Goodfellow et al. 2016).
Neural Network Training with Physics Objectives. For the remainder of the paper, we consider joint optimization problems, where f denotes composition of a neural network parameterized by O and physics solver. Using classical network optimizers for minimizing equation 2 is inefficient in this case since data normalization in the network output space is not possible and the classical initialization schemes cannot normalize the effects of the physics solver. As such, they unsuited to capture the strong coupling between optimization parameters typically encountered in physics applications. While tends to result in large overshoots in the updates when the involved physics solver is ill-conditioned. As we will demonstrate, this leads oversaturation of neurons, hampering the learning capability of the neural network.
seems promising for these cases, the involved Jacobian inversion
2.1 AN ILL-CONDITONNED TOY EXAMPLE
To illustrate the argumentation far, we consider a data set sampled from ((x) - (6in(x)),, cos (9r))
for x E [-1, 1]: We train a neural network to describe this data set by using the loss function:
(y ÿ: 7: (00 0²)² 2
1 [ ("))
(4)
Here, we denote vector components by superscripts. For a scale factor of y = 1, we receive the well-conditioned mean squared error loss. However, 1 becomes increasingly ill-conditioned as is decreased, imitating the effects of physics solver. For real-world physies solvers, the situation would be even more complex since these scales usually vary strongly in direction and magnitude across different data points and optimization steps. We use a small neural network with a single hidden layer with 7 neurons and a tanh activation. We then compare training with the well-conditional y 1 loss against ill-conditioned y 0.01 loss. In both cases, we train the network using both Adam and Geaus-Neetton as representatives of gradient- based and higher-order optimizers, respectively. The results are shown in figure 0
In the well-conditioned case, Adam and Gauss-Newton behave similarly, decreasing the loss by about three orders of magnitude. However, in the ill-conditioned case, both optimizers fail to minimize the objective beyond a certain point. To explain this observation, we first illustrate the behavior from the physics viewpoint by considering the trajectory of the network output f(x) for single value æ during training (figure T] right). For y == 1, Adam optimizes the network to accurately predict (xx) while for 0.011 the updates neglect the second component preventing Adam move efficiently along the small-scale coordinate (blue curve in figure [[, right). To illustrate the situation from the viewpoint of the network, we consider the variance the outputs of specific neurons over different x (figure [T. middle). When 7 = 1, all neurons process information by producing different outcomes for different x. However, for 7 = 0.01, Gauss-Newton's 's inversion of the small- scale component y² results in large updates, leading to an oversaturation of neurons (red curve in figure Пp. middle). These neurons stop processing information, reducing the effective capacity of the network and preventing the network from accurately fitting ÿ. Facing these problems. natural questions arises: Is it possible to construct an algorithm that can successfully process the inherently different scales of a physics solver while training a neural network at the same time?
3

---

Published as a conference paper at ICLR 2022
a
Loss curve
Neuron saturation
Physical trajectory
Addm
GN
10"
1.0
Well-condittonell
Neuron value
Loss
10
0.5
10
-05-
10
20
b)
10-
GN loep variares
10a
Ill-conditioned
Neuron value
Loss
- GX a Acandoin mock the oeger saa
10
-10-
Start
Target
26
Wall clock time [sec]
Wall clock time [sec]
y4
Figure I: Results of the learning problem of section 2.1 Optimization is performed with a) a well- conditioned loss and b) ill-conditioned loss. Plots show loss curves over training time (left), data set mean and standard deviation of the output of a neuron output over training time (middle), and the training trajectory of data point (right).
2.2 UPDATES BASED ON HALF-INVERSE JACOBIANS
We propose a novel method for optimizing neural networks with physics objectives. Since pure physics or neural network optimization can be thought of as special cases of the joint optimization, we analogously look for a potential method in the continuum of optimization methods between gradient descent and Gauss-Newton. We consider both of them to be the most elementary algorithms representing network and physics optimizers, respectively. The following equation describes updates that lie between the two.
Oy 00
OL
A k) k)= -
(5)
By
Here, the exponent of the Jacobian denotes the following procedure defined with the aid of the singular value decomposition J UAVT:
;= V
When x = 1, equation reduces to the well-known form of gradient descent. Likewise, the case x = -1 yields Gauss-Newton since the result of the Jacobian exponentiation then gives the pseudoinverse of the Jacobian. Unlike other possible interpolations between gradient descent and Gauss-Newtte exponentiation by in equation significantly affects the scales inherent in the Jacobian. This is highly important to appropriately influence physics and neural network scales.
To determine K, recall our goal to perform update steps which are optimal in both O- and y- space. However, since any update AO and its corresponding effect on the solver output Ay are connected by the inherent scales encoded in the Jacobian, no single к exists that normalizes both the same time. Instead, we distribute the burden equally between network and physics by choosing 1/2. From geometric viewpoint, the resulting update can regarded as a steepest descent step when the norm to measure distance is chosen accordingly. This alternative way to approach our method is explained in the appendix (A.2p and summarized in table
For batch size b and learning rate 7,, we define the following update step for our method by stacking dy; and loss gradients By. 3
network-solver Jacobians 50 æ;
8L of different data points (æ,, 9:):
4

---

Published as a conference paper at ICLR 2022
Table 1: Optimization algorithms viewed steepest descent algorithm the given 1.2-norms.
Optimization Method
performs Steepest Descent:
Norm (0-space)
Norm (y-space)
Gradient Descent
in Parameter Space
Gauss-Newton Ours
in Physics Space in Intermediate Space
JO
12/40
11y-
Pu
1/2
0L
dyi
22
duu 8e
0L
A =
dya #2.92
aL
80
22%
Oye
Besides batch size b and learning rate 1), we specify truncation parameter t additional hyperparameter enabling us to suppress numerical noise during the half-inversion process in equation As with the computation of the pseudoinverse via SVD, set the result of the - exponentiation of every singular value smaller than T to O.
The use of half-inversion instead of full inversion helps prevent exploding updates of network parameters while still guaranteeing substantial progress in directions of low curvature. With the procedure outlined above, we arrived at balanced method that combines the advantages of optimization methods from deep learning and physics. As our method uses half-inverse Jacobians multiplied with gradients we refer them in short as half-inverse gradients (HIGs).
Half-inverse Gradients in the Toy Example. With the definition of HIGs, we optimize the toy example introduced in section 2.1 The results in figure show that for y 1, HIGs minimize the objective as well as Adam and Gauss-Newton's is method. More interestingly, HIGs achieve a 0.01. On the hand, the physics trajectory
better result than the other two methods for (figure [b, right) highlights that HIGs can process information along the small-scale component y well and successfully progress along this direction. On the other hand, by checking neuron saturation (figure Db, middle), we see that HIGs - in contrast to Gauss Newton avoid oversaturating neurons.
2.3 PRACTICAL CONSIDERATIONS
Computational Cost. A HIG update step consists of constructing the stacked Jacobian and computing the half-inversion. The first step can be efficiently parallelized on modern GPUs, and therefore induces runtime cost comparable to regular backpropagation at the expense of higher memory requirements. In situations where the computational cost of the HIG step dominated by the half-inversion, memory requirements can be further reduced by parallelizing the Jacobian computation only partially. At the heart the half-inversion lies a divide and conquer algorithm for the singular value decomposition (Trefethen & Bau 1997). Hence, the cost of HIG step scales as ((10) ([])))-)(()(()) i.e. is linear in the number of network parameters 0 and quadratic in the batch size b and the dimension of the physical state ly). Concrete numbers for memory requirements and duration of a HIG step are listed in the appendix.
Hyperparameters. Our method depends on several hyperparameters. First, we need a suitable choice of the learning rate. The normalizing effects of HIGs allow for larger learning rates than commonly used gradient descent variants. We are able to use Second, the batch size b affects the number of data points included in the half-inversion process. should be noted that the way the feedback of individual data points is processed is fundamentally different from the standard gradient optimizers: Instead the averaging procedure of individual gradients of mini batch, batch. Consequently, the quality of updates increases with higher batch size. However, overly large batch sizes can cause the Jacobian become increasingly ill-conditioned and destabilize the learning progress. In appendix([] we discuss the remaining parameters 7 and x with several ablation experiments to illustrate their effects in detail.
for many of our experiments.
approach constructs an update that is optimal for the complete
5

---

Published as a conference paper at ICLR 2022
a)
b)
Asam
SGD
Adam bel2e
0M
Adam
2-32
Adaddelia
10:2
RMSprop
10"
Loss
Loss
10.4
10-4
-1
10:4
t
0
1500
300G
1500
30000
Wall clock time [sec]
Wall clack time [sec]
Nonlinear oscillator system: Time evolution controlled by HIG-trained neural network. Its inferred output is shown in blue in the background. b) Loss curves for different optimization methods. e) Loss curves for Adam, GN, and HIG with different batch sizes b.
Figure 2:
EXPERIMENTS
We evaluate method on three physical systems: controlling nonlinear oscillators, the Poisson problem, and the quantum dipole problem. Details of the numerical setups are given in the appendix show results with
along with results for broad range of hyperparameters. For fair comparison, the best set of hyperparameters for each of the methods below and plot the loss against wall clock time measured in seconds. All learning curves are recorded on previously unseen data set.
3.1 CONTROL OF NONLINEAR OSCILLATORS
First, consider control task for system of coupled oscillators with nonlinear interaction term. This system is practical importance in many areas physics, such as solid state physics (Ibach & Lüth 2003). Its equations of motions governed by the Hamiltonian
1 p 2 + (z: x441)* + u(t) xq :
H( Ps: 1) Σ
where z; and p; denote the Hamiltonian conjugate variables of oscillator i, the interaction strength, and the vector e specifies how to scalar-valued control function u(t) is applied. In our setup, we train neural network to learn the control signal u(t) that transforms given initial state into given target state with 96 time steps integrated by a 4th order Runge-Kutte scheme. We use dense neural network with three hidden layers totalling 2956 trainable parameters and ReLU activations.
The Mean-S Soinnn-Suuur----Eee Error loss is used to quantify differences between predicted and target state. A visualization of this control task shown in figure 2
Optimizer comparison. The goal of our first experiments is to give broad comparison of the proposed HIGs with commonly used optimizers. This includes stochastic gradient descent (SGD), Adagrad Duchi al., 2011). Adadelta (Zeiler 2012). RMSprop (Hinton et al. 2012). Adam (Kingma & Ba 2015). and Gauss-Newton (GN) applied to mini batches. The results are shown in figure 2p where all curves show the best for each optimizer with suitable hyperparameters independently selected, as explained the appendix. We find that the state-of-the-an optimizers stagnate early, with Adam achieving the best result with final loss value of 10-4 In comparison, our method and GN converge faster, exceeding Adam's 's accuracy after about three minutes. While GN exhibits stability problems, the best stable run from our hyperparameter search reaches loss value of 10-6, HIGs, the other hand, yield the best result with loss value of 10-7. These results clearly show the potential of our method process different scales of the physies solver more accurately and robustly. They also make clear that the poor result of the widely-used network optimizers cannot be attributed to simple numerical issues as HIG converges to better levels of accuracy with otherwise identical setup.
Role of the batch size. b central parameter of our method. The results are shown in figure 22. We observe that for
We conduct multiple experiments using different values for the batch size

---

Published as a conference paper at ICLR 2022
a)
b)
Acam
4.10-
10"
2
Loss
Adam 1e-1--
10--
10
6-10
Wall clock time [sec]
Wall clock time [sec]
Figure 3: Poisson problem: a) Example of source distribution p (bottom) and inferred potential field (top). b) Loss curves of Adam and HIG training for different learning rates 1/. c) Loss curves of Adam (17 0.0001), and HIG (17 0.02) pretrained with Adam.
Adam, all runs converge about equally quickly while HIGs and GN show improvements from larger batch sizes. This illustrates important difference between Adam and HIG: Adam uses an average of gradients of data points in the mini batch, which approaches its expectation for large Further increasing the batch size has little influence on the updates. In contrast, our method includes the individual data point gradients without averaging. As shown in equation 7 we construct updates that are optimized for the whole batch by solving a linear system. This gives our method the ability to hit target states very accurately with increasing batch size. To provide further insights into the workings of HIGs, we focus on detailed comparisons with Adam as the most popular gradient descent variant.
3.2 POISSON PROBLEM
Next consider Poisson's equation to illustrate advantages and current limitations of HIGs. Poisson problems play an important role in electrostatics, Newtonian gravity, and fluid dynamics (Ames) 2014). For source distribution p(xx), the goal is find the corresponding
potential field ((x) fulfilling the following differential equation:
A =
Classically, Poisson problems solved by solving the corresponding system of linear equations on the chosen grid resolution. Instead, we train a dense neural network with three hidden layers and 41 41408 trainable parameters solve the Poisson problem for given right hand side p. We consider a two-dimensional system with a spatial discretization of 8 x degrees of freedom. An example distribution and solution for the potential field are shown figure 3
Convergence and Runtime. Figure 3b shows learning curves for different learning rates when training the network with Adam and HIGs. As we consider a two-dimensional system, this optimization task is challenging for both methods and requires longer training runs. We find that both Adam and HIGs are able to minimize the loss by up to three orders of magnitude. The performance of Adam varies, and its two runs with larger quickly slow down. In terms of absolute convergence per time, the Adam curve with the smallest 7 shows advantages in this scenario. However, choosing a log-scale for the time axis reveals that both methods have not fully converged. In particular, while the Adam curve begins to flatten at the end, the slope of the HIG curve remains constant and decreases with steeper slope than Adam. The performance of Adam can be explained by two reasons. First, the time to compute a single Adam update is much smaller than for HIGs, which requires the SVD solve from equation 6 While these could potentially be sped up with appropriate methods (Foster et al. 2011 Allen- Zhu & Li 2016), the absolute convergence per iteration, shown in the appendix figure 7 shows how much each HIG update improves over Adam. Second, compared to the other examples, the Poisson problem is relatively simple, requiring only single matrix inversion. This represents a level of difficulty which Adam is still able to handle relatively well.
HIGs with Adam Pretraining. To further investigate the potential of HIGs, we repeat the training, this time using the best Adam model from figure 3p for network initialization. While Adam progresses slowly, HIGs able to quickly improve the state of the neural network, resulting in
7

---

Published as a conference paper at ICLR 2022
a
0.4
b)
Adam: q 2e-33
C) 0.6
I... Adam: LE
HBG:
1/0
Adem: Se-D4
Adamo HE
HIG:
0.0
10P
Adam:
HG:: gm0.7
Loss function
HIG: 0.5.5
Loss
1000
HG: n-0.33
10:4
10.4
2000
40000
0.0
500
1000
Wall clock time [sec]
Wall clock time [sec]
Figure 4: Quantum dipole: a) A transition of two quantum states in terms of probability amplitude [V, x1). controlled by a HIG-trained neural network. Its inferred output shown in blue in the background. b) Loss curves with Adam and HIGs for different n7. c) Low-energy (LE) and high- energy (HE) loss with Adam (7) = 0.0001) and HIG (17 = 0.5).
significant drop of the loss values, followed by faster descent than Adam. Interestingly, this experiment indicates that the HIG updates are able to improve aspects of the solution which Adam agnostic to. Despite outlining the potential gains from faster SVD calculations, this example also highlights the quality of the HIG updates for simpler PDEs.
3.3 QUANTUM DIPOL.
As a final example, we target the quantum dipole problem, a standard control task formulated on the Schrödinger equation and highly relevant in quantum physics (Von Neumann 2018). Given an initial and a target state, we train a neural network to compute the temporal transition function u(t) in an anfinitewwell potential V according the evolution equation of the physical state T.
a ( (-A A-A+V + u(t) 2) L
(10)
We employ a modified Crank-Nicolsen scheme (Winckel et al. 2009) for the discretization of spatial and temporal derivatives. Thus, each training iteration consists of multiple implicit time integration steps 384 in our setup - for the forward as well as the backward pass of each mini-batch. The control task consists of inferring signal that converts the ground state to given randomized linear combination of the first and the second excited state. We use a dense neural network with three hidden layers, 9484 trainable parameters and lanh activations. Similarity in quantum theories is quantified with inner products; therefore, our loss function is given by L(V., a, V.,) = 1 1-1(V,, 1,)) A visualization of this control task is shown in figure 44.
Speed and Accuracy. We observe that HIGs minimize the loss faster and reach a better final level of accuracy than Adam (figure 4b). While the Adam run with the largest learning rate drops faster initially, its final performance is worse than all other runs. In this example, the difference between the final loss values is not as large as for the previous experiments. This is due to the numerical accuracy achievable by pure physics optimization, which for choice of parameters is around 10 10-6, 6, Hence, we can not expect to improve beyond this lower bound for derived learning problems. Our results indicate that the partial inversion of the Jacobian successfully leads to the observed improvements in convergence speed and accuracy.
Low and High Energy Components. The quantum control problem also serves to highlight the weakness of gradient-based optimizers in appropriately processing different scales of the solutions. In the initial training stage, the Adam eurves stagnate a loss value of 0.5. This most pronounced for 10-4 in dark blue. To explain this effect, we recall that our learning objective targets transitions to combinations average with equal weight in the training data. Transitions to the energetically higher states are more difficult and connected to smaller scales in the physics solver, causing Adam to fit the lower- energetic component first. In contrast, Jacobian via the half-inversion more efficiently. As a consequence, the loss curves decrease faster below 0.5. We support this explanation by explicitly plotting separate loss curves in figure 4
the 1st and 2nd excited quantum states, and both states appear on
method is constructed to process small scales in the
8

---

Published as a conference paper at ICLR 2022
quantifying how well the low and high energy component of the target state was learned. Not only does Adam prefer to minimize the low-energy loss, it also increases the same loss again before it is able to minimize the high-energy loss. In contrast, we observe that HIGs minimize both losses uniformly. This is another indication for the correctness of the theory outlined above of an more even processing of different scales in joint physies and neural network objectives through our method.
4 RELATED WORK
Optimization algorithms. Optimization on continuous spaces is huge field that offers vast range of techniques (Ye et al., 2019). Famous examples are gradient descent (Curry [944)) Gauss-Newton's s method (Gill & Murray, 1978), Conjugate Gradient (Hestenes et al. 1952), or the memory BFGS algorithm (Liu & Nocedal 1989). In deep learning, the preferred methods instead rely on first order information in the form of the gradient, such as SGD (Bottou 2010) and RMSProp (Hinton et al. 2012). Several methods approximate the diagonal of the Hessian to improve scaling behavior, such as Adagrad (Duchi et al. 201T), Adadelta (Zeiler, 2012), and most prominently, Adam (Kingma & Ba, 2015). However, due to neglecting inter- dependencies of parameters, these methods are limited in their capabílities to handle physical learning objectives. Despite the computational cost, higher-order methods have also been studied in deep learning (Pascanu & Bengio, 2013) Practical methods have been suggested by using Kroenecker- faceooizziattooozitiionno of the Fisher matrix (Martens & Grosse 2013). iterative linear solvers (Martens] 2010 or by recursive approximations of the Hessian (Botev et al. 2017). To the best of our knowledge the only other technique specifically targeting optimization of neural networks with physies objectives is the inversion approach from Holl et al. (202T). However, their updates are based on inverse physics solvers, while we address the problem by treating network and solver as an entity and half-inverting its Jacobian. Thus, we work on the level of linear approximations while updates based on physics inversion are able to harness higher-order information provided that an higher-order inverse solver exists. Additionally, they compute their update by averaging gradients over differrnt data points, in line with typical gradient- neural network optimizers. HIGs instead process the feedback of different data points via collective inversion.
Incorporating physics. Many works involve differentiable formulations of physical models, e.g., for robotics [Toussaint et al., 2018), to enable deep architectures (Chen er al.. 2018). as means for scene understanding (Battaglia et al., 2013 Santoro et al. (017). or the control of rigid body environments de Avila delbute-Peres et al. (2018). Additional works have shown the advantages of
physical loss formulations (Greydanus et al. 2019 Cranmer et al.) 2020). Differentiable simulation
methods were proposed for a variety of phenomena, e.g. for fluids (Schenck & Fox 2018), PDE discretizations (Bar-Sinai et al., 2019), molecular dynamics (Wang et al. 2020), reducing numerical errors (Um et al.][222) and cloth (Liang et al., (2019; RRahhed et al. 2020). It is worth noting that none of these works question the use of standard deep learning optimizers, such as Adam. In addition, by now a variety of specialized software frameworks are available to realize efficient
implementations (Hu al. 2020 Schoenholz & Cubuk 2019 Holl al. 2020).
5 DISCUSSION AND OUTLOOK
We have considered optimization problems of neural networks in combination with physical solvers and questioned the current practice of using the standard gradient- based network optimizers for training. Derived from an analysis of smooth transitions between gradient descent and Gauss- Newton's 's method, our novel method learns physics modes more efficiently without overly straining the network through large weight updates, leading to faster and more accurate minimization of the learning objective. This was demonstrated with a range of experiments.
We believe that our work provides starting point for further research into improved learning methods for physical problems. Highly interesting avenues for future work are efficient methods for the
of the Jacobian matrix, applying HIGs to physical systems exhibiting chaotic behavior or to more sophisticated training setups (Battaglia et al. 2013 Ummenhofer et al., 2020 Pfaff et al. 2020).
9

---

Published as a conference paper at ICLR 2022
ACKNOWLEDGEMENTS
This work was supported by the ERC Consolidator Grant CoG-2019-863851 SpaTe, and by the DFG Transreero 109 DGD. We would also like to express our gratitude to the reviewers and the area chair for their helpful feedback.
REPRODUCIBILITY STATEMENT
Our code for the experiments presented in this paper is publicly available at https://githhbb com/tuu pbs/ /ha 1 inverse grwa-rraaiee ent Additionally, the chosen hyperparameters are listed in the appendix along with the hardware used to our simulations.
REFERENCES
Peter Adby. Introduction optimization methods. Springer Science and Business Media, 2013.
Zeyuan Allen-Zhu and Yuanzhi Li. Lazysvd: Even faster svd decomposition yet without agonizing pain. Advances in Neural Information Processing Systems,
2016.
William Ames. Numerical methods for partial differential equations. Academic press, 2014.
Yohai Bar-Sinai, Stephan Hoyer, Jason Hickey, and Michael Brenner. Learning data-driven discretizations for partial differential equations. Proceedings of the National Academy of Sciences, 116331):5334-115599 2019.
Peter Battaglia, Razvan Pascanu, Matthew Lai, Danilo Jimenez Rezende, et al. Interaction networks for learning about objects, relations and physics. In Advances in Neural Information Processing Systems, pp. 4502-4510 2016.
Peter W Battaglia, Jessica B Hamrick, and Joshua B Tenenbaum. Simulation as an engine of physical scene understanding. Proceedings of the National Academy of Sciences, 110(45), 2013.
Aleksandar Botev, Hippolyt Ritter, and David Barber. Practical gauss-newton optimisation for deep learning. In International Conference on Machine Learning, pp. 557-565. PMLR, 2017.
Léon Bottou. Lcage-scale machine learning with stochastic gradient descent. In Proceedings of COMPSTAT pp. 177-186. Springer. 2010.
Tian Qi Chen, Yulia Rubanova, Jesse Bettencourt, and David Duvenaud. Neural ordinary differential equations. In Advances in Neural Information Processing Systems, pp. 6571-6583, 2018.
Miles Cranmer, Sam Greydanus, Stephan Hoyer, Peter Battaglia, David Spergel, and Shirley Ho. Lagrangian neural networks, arXiv:2003.04630 2020.
Haskell B Curry. The method of steepest descent for non-linear minimization problems. Quarterly of Applied Mathematics, 2(3):258-22) 261, 1944.
Filipe de Avila Belbute-Peres, Kevin Smith, Kelsey Allen, Josh Tenenbaum, and Zico Kolter. End-to-end differentiable physics for learning and control. In Advances in Neural Information Processing Systems, 2018.
John Duchi, Elad Hazan, and Yoram Singer. Adaptive subgradient methods for online learning and stochastic optimization. Journal of machine learning research, 12(7), 2011.
Blake Foster, Sridhar Mahadevan, and Rui Wang. A gpu-based approximate svd algorithm. In International Conference on Parallel Processing and Applied Mathematics, pp. 569-578. Springer, 2011.
Philip E Gill and Walter Murray. Algorithms for the solution of the nonlinear least-squares problem.
SIAM Journal on Numerical Analysis, 15(5):977-992. 1978.
Ian Goodfellow, Yoshua Bengio, Aaron Courville, and Yoshua Bengio. Deep learning, volume 1. MIT Press, 2016.
10

---

Published as a conference paper at ICLR 2022
Samuel Greydanus, Misko Dzamba, and Jason Yosinski. Hamiltonian neural networks. In Advances in Neural Information Processing Systems, pp. 15353-15363, 2019.
Magnus R Hestenes, Eduard Stiefel, et Methods of conjugate gradients for solving linear systems. Journal of research of the National Bureau of Standards, 49(6)409-436, 1952.
Geoffrey Hinton, Nitish Srivastava, and Kevin Swersky. Lecture 6a, overview of mini-batch gradient descent. Neural networks for machine learning, 14(8):2, 2012.
Philipp Holl, Vladlen Koltun, and Nils Thuerey, Learning to control pdes with differentiable physics. In International Conference on Learning Representations (ICLR), 2020.
Philipp Holl, Vladlen Koltun, and Nils Thuerey. Physical gradients for deep learning. arXiv, 2109.15048, 2021.
Yuanming Hu, Luke Anderson, Tzu-Mao Li, Qi Sun, Nathan Carr, Jonathan Ragan-Kelley, and Frédo Durand. Difftaichi: Differentiable programming for physical simulation. International Conference on Learning Representations (ICLR), 2020.
H. Ibach and H. Lüth. Solid-State Physics. Advanced texts in physics. Springer, 2003. ISBN 9783540438700.
Diederik Kingma and Jimmy Ba. Adam: A method for stochastic optimization. In International Conference on Learning Representations (ICLR), 2015.
Dmitrii Kochkov, Jamie A Smith, Ayya Alieva, Qing Wang, Michael P Brenner, and Stephan Hoyer. Machine aeaminng-aceeeerated computational fluid dynamics. Proceedings of the National Academy of Sciences, 118(21), 2021.
Alex Krizhevsky, Ilya Sutskever, and Geoffrey E Hinton. Imagenet classification with deep convolutional neural networks. In Advances in Neural Information Processing Systems, 2012.
Peter Lax. Linear algebra and its applications, volume Hoboken Wiley, 2014.
Junhang Liang, Ming Lin, and Vladlen Koltun. Differentiable cloth simulation for inverse problems. In Advances in Neural Information Processing Systems, pp. 771-780, 2019.
Dong Liu and Jorge Nocedal. On the limited memory bfgs method for large scale optimization. Mathematical programming, 45(1-3):503-553 1989.
James Martens. Deep learning via hessian free optimization. International conference on machine learning, pp. 735-742, PMLR, 08 2010.
James Martens and Roger Grosse. Optimizing neural networks with kronecker-factored approximate curvature. In International conference machine learning, pp. 2408-2411 PMLR, 2015.
Jeremy Morton, Antony Jameson, Mykel J Kochenderfer, and Freddie Witherden. Deep dynamical modeling and control of unsteady fluid flows. In Advances in Neural Information Processing Systems, 2018
Razvan Pascanu and Yoshua Bengio. 10X15:1301.3584, 2013.
Revisiting natural gradient for deep networks.
Tobias Pfaff, Meire Fortunato, Alvaro Sanchez-Gonzalez, and Peter W Battaglia. Learning mesh- based simulation with graph networks. arXiv preprim::000..03409, 2020.
Maziar Raissi, Paris Perdikaris, and George Karniadakis.
neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. Journal of Computational Physics, 378:686-707, 2019.
Abdullah Haroon Rasheed, Victor Romero, Florence Bertails-Descouber Stefanie Wuhrer, Jean- Sébastien Franco, and Arnaud Lazarus. Learning to measure the static friction coefficient in cloth contact. In The Conference Computer Vision and Pattern Recognition, 2020.

---

Published as a conference paper at ICLR 2022
Adam Santoro, David Raposo, David GT Barrett, Mateusz Malinowski, Razvan Pascanu, Peter Battaglia, and Timothy Lillicrap. A simple neural network module for relational reasoning. 170.001706. 2017.
Connor Schenck and Dieter Fox. Spnets: Differentiable fluid dynamies for deep neural networks. In Conference on Robot Learning, pp. 317-335, 2018.
Samuel S Schoenholz and Ekin D Cubuk. Jax, md: End-to-end differentiable, hardware accelerated, molecular dynamics in pure python. arXiv:19/2.04232. 2019.
David Silver, Julian Schrittwieser, Karen Simonyan, Ioannis Antonoglou, Aja Huang, Arthur Guez, Thomas Hubert, Lucas Baker, Matthew Lai, Adrian Bolton, al. Mastering the game of Go without human knowledge. Nature, 550(7676), 2017.
Ilya Sutskever, Oriol Vinyals, and Quoc V Le. Sequence to sequence learning with neural networks. In Advances Neural Information Processing Systems, pp. 3104-3112, 12, 2014.
Jonathan Tompson, Kristofer Schlachter, Pablo Sprechmann, and Ken Perlin. Accelerating eulerian fluid simulation with convolutional networks. In Proceedings of Machine Learning Research, pp- 3434-3433, 2017.
Mare Toussaint, Kelsey Allen, Kevin Smith, and Joshua B Tenenbaum. Differentiable physics and stable modes for tool-use and manipulation planning. In Robotics: Science and Systems, 2018.
L. Trefethen and D. Bau. Numerical Linear Algebra. EngineeringPro collection. Society for Industrial and Applied Mathematics, 1997. ISBN 978089871493 14876.
Kiwon Um, Robert Brand, Yun Raymond Fei, Philipp Holl, and Nils Thuerey. Solver-in-the-lone Learning from differentiable physics to interact with iterative pde-solvers. Advances in Neural Information Processing Systems, 2020.
Benjamin Ummenhofer, Lukas Prantl, Nils Thuerey, and Vladlen Koltun.
Lagrangian fluid simulation with continuous convolutions. In International Conference on Learning Representations, 2020.
John Von Neumann. Mathematical foundations quantum mechanics. Princeton university press, 2018.
Wujie Wang. Simon Axelrod, and Rafael CombarlBo Differentiable molecular simulations for control and learning. arXiv:2003. 00868, 2020.
Greg von Winckel, Alfio Borzi, and Stefan Volkwein. globalized newton method for the accurate solution of dipole quantum control problem. SIAM Journal Scientific Computing, 31(6): 4176-4203, 2009, doi: 10.1 10.1137/097749611.
Tailin Wu and Max Tegmark. Toward an artificial intelligence physicist for unsupervised learning. Physical Review E, 100(3):03311,, 2019.
Nan Ye, Farbod
and Tiangang Cui. Optimization methods for inverse problems. In 2017 MATRIX Annals, pp. 121-140. Springer, 2019.
Matthew Zeiler. Adadelta: adaptive learning rate method. arXi:: 2012.
Olgierd Cecil Zienkiewicz, Robert Leroy Taylor, Perumal Nithiarasu, and JZ Zhu. The finite element method, volume 3. McGraw-hill London, 1977.
12

---

Published as a conference paper at ICLR 2022
APPENDIX
A FURTHER DETAILS ON OPTIMIZATION ALGORITHMS
Our work considers optimization algorithms for functions of the form f(x;0) y with 0, AO Rt, denoting weight vector and weight update vector, respectively, while æ E T and y E R denote input and output. The learning process solves the minimization problem argming 0), y) via sequence 0 and we target losses of the form (ax,, D;; 0) = Σ. (f 0), i() with being an index for multiple data points (i.c., observations). denotes the L2-loss E3 || 23 with j referencing the entries of mini batch of size b.
0 m/.. Here, are the reference solutions,
A.1 UPDATE STEP OF THE GAUSS-NEWTON ALGORITIMM
Using this notation, the update step of the Gauss-Newton algorithm (Adby 2013 for f) 1 given by:
AOGN -
((
Oy a0
T
dy 00
Oy 30
T
OL
T
(II)
8y
The size of the Jacobian matrix is given by the dimensions of y- and O-space. For a full-rank Jacobian corresponding to con-constrained optimization, the Gauss-Newton update is equivalent to:
-
18y 30
aL
A0CN
(12)
Even a constrained setting, we can reparametrize the coordinates to obtain an unconstrained optimization problem the accessible manifold and rewrite AdGN similarly. This shortened form of the update step is given in equation 3] and is the basis for our discussion in the main text.
A.2 GEOMETRIC INTERPRETATION AS STEEPEST DESCENT ALGORITHMM
It is well-known that the negative gradient of function L(0) points in the direction of steepest descent leading to the interpretation gradient descent as a steepest descent algorithm. However, the notion of steepest descent requires defining measure of distance, which is in this case the usual L2-norm in By using different metrics, we can regard Gauss-Newton and HIG steps as steepest descent algorithms well.
Causs-Newtor updates. The updates AdEN be regarded as gradient descent in y up first. order in the update step. This can be seen with a simple equation by considering how these updates change y.
(yy 00
OT. Oy
Ay
AdGN + o(AfGN)
+ o(A@GN)
(13)
In figure [ of the main paper, this property is visible the physics trajectories for the well- conditioned case, where L(y) is uniform L2-loss and hence, gradient descent in produces a straight line to the target point. The Gauss-Newton curve first shows several steps in varying directions the higher-order terms from the neural network cannot be neglected yet. However, after this initial phase the curve exhibits the expected linear motion.
The behavior of GN perform steepest descent the y- g-manifold stands in contrast gradient descent methods, which instead perform steepest descent on the m-manifold. This geometric view is the basis for alternative way derive method that is presented below.
13

---

Published as a conference paper at ICLR 2022
HIG updates. HIG updates be regarded steepest descent algorithm, again up to first order in the update step, when measuring distances of (-vectors with the following semi-norm:
| H == [|
(14)
Here || -|ll denotes the usual 12-norm and J 80 the Jacobian of network and solver. The exponentiation performed as explained the main text, with J = UAVT being the SVD, and J3/4 given by VA3/AUT Additionally, we will use the natural map between dual vector and vector aL öy
(,)) and the loss gradient g
To prove the claim above, expand the loss around an arbitrary starting point Op:
L + A0)) L(y(On)) (g J, o(20)
(15)
The first term on the right-hand side is constant and the third term is neglected according to the assumptions of the claim. Hence, we investigate for which fixed-length AO the second term decreases the most:
arg min || ||... =const. (1g - A0)) arg min (ig 22/ ,,
||0| l.......o.st.
arg min cos |1g |4/4|| -
(16)
arg min (CSS 7)
In the first step above, we split the Jacobian JT VAUT (VAVAVT J1/4
/3/4
y denotes the angle between JJ/4gT and J3/4 A0. This expression minimized for y -T, meaning the two vectors have to be antiparallel:
331=0= = -J/ig
(17)
This requirement is fulfilled by the HIG update A0HIG = T and is therefore steepest descent method, which concludes our proof.
This presents another approach view HIGs as an interpolation between gradient descent and Gauss-Newton's 's method. More precisely, gradient descent performs steepest descent in the usual L2-norm in O-space () Considering only terms up to linear order, Gauss-Newton performs steepest descent in the L2-norm y-space J[l) The HIG update (|]]3//0]]) lies between these two methods. The quarter factors in the exponents result from the additional factor of 2 that has to be compensated for when considering 1.2nnors.
A.3 STABILITY OF INVERSIONS THE CONTEXT OF PHYSICAL DEEP LEARNING.
In the following, we illustrate how the full inversion of GN can lead to instabilities at training time. Interestingly, physical solvers are not only cause of small singular values the Jacobian. They also occur when applying equation 12 to mini batch train neural network and are not caused by numerical issues. Consider the simple case of two data points (æ1, ②1) and (æ2, 92) and one-dimensional output. Let f be the neural network and J the Jacobian, which is in this case the gradient of the network output. Then equation 12 yieelss:
Jg(en) (y(2z)
91 f(z)) 322
f(
A0GN
(18)
Next, we linearly approximate the second row by using the Hessian by assuming the function be learned is f. i.e. f(x)) = yn and f(xz) = y2- Neglecting terms beyond the linear approximation, receive:
14

---

Published as a conference paper at ICLR 2022
f(a1) - y/ Jy(x1) + (2 (022 10-1) x1) ) ( [(xi)) - yl (Jz(x1) fx(x) (œ2
])(x1)
x (19)
Considering the case of two nearby data points, i.e. æ2 æ1 being small, the two row vectors in the stacked Jacobian on the left-hand side are similar, i.e. the angle between them small. This leads to small singular value of the stacked Jacobian. In the limit of æ æ both vectors are linearly dependant and hence, one singular value becomes zero.
Moreover, even if T2 is not close 11, small singular values occur if the batch size increases: for a growing number of row vectors it becomes more and more likely that the Jacobian contains similar linearly dependent vectors.
After inversion, a small singular value becomes large. This leads to a large update AdON when the right-hand side of equation 19 Oveerlaps with the corresponding singular vector.
This can easily happen if the linear approximation of the right-hand side is poor, for instance when f is a solution to an inverse physics problem. Then f can have multiple modes and can, even within a mode, exhibit highly sensitive or even singular behavior.
In turn, applying large updates to the network weights naturally can lead to the oversaturation of neurons, illustrated above. and diverging training runs general.
As illustrated in the main paper, these inherent problems of GN are alleviated by the partial inversion of the HIG. yields fundamentally different order of scaling via its square-root inversion, which likewise does not guarantee that small singular values lead to overshoots (hence the truncation), but in general strongly stabilizes the training process.
B EXPERIMENTAL DETAILS
In the following, we provide details of the physical simulations used for our experiments in section 3 of the main paper. For the different methods, we use the following abbreviations: half-inverse gradients (HIG), Gauss-Newton's 's method (GN), and stochastic gradient descent (GD). Learning rates are denoted by 7), batch sizes by b, and truncation parameters for HIG and GN by All loss results are given for the average loss over a test set with samples distinct from the training data set.
For each method, we run a hyperparameter search for every experiment, varying the learning rate by several orders of magnitude, and the batch size in factors two. Unless noted otherwise, the best runs in terms of final test loss were selected and shown in the main text. The following sections contain several examples from the hyperparameter search to illustrate how the different methods react to the changed settings.
Runtime Measurements Runtimes for the non-linear chain and quantum dipole were measured on a machine with Intel Xeon 6240 CPUs and NVIDIA GeForce RTX 2080 Ti GPUs. The Poisson experiments used an Intel Xeon W-2235 CPU with NVIDIA Quadro RTX 8000 GPU, We experimentally verified that these platforms yield an on-par performance for our implementation. As deep learning API we used TensorFlow version 2.5. If not stated otherwise, each experiment retained the default settings.
All runtime graphs in the main paper and appendix contain wall-clock measurements that include all steps of learning run, such as initialization, in addition to the evaluation time of each epoch. However, the evaluations of the test sets to determine the performance in terms of loss are not included. As optimizers such as Adam typically performs a larger number of update steps including these evaluations would have put these optimizers at an unnecessary disadvantage.
B. TOY EXAMPLE (SECTION 2.1
For the toy example, the target function is given by f(x) (sin(62), (sis)))).cos((())).. We used dense neural network consisting of one hidden layer with neurons and tanh activation, and an output layer with 2 neurons and linear activation. For training, use 1024 data points uniformly sampled
15

---

Published as a conference paper at ICLR 2022
Table 2: Hyperparameters for different optimization algorithms in figure 2
Method
Adadelta
Adagrad
Adam
GN
HIG
RMSprop
SGD
b
512
512
512
128
128
512
512
7 0.1
0.1
3-10-*
10-4
0.1
10-6
10-6
Table 3: Nonlinear oscillators: memory requirements, update duration and duration of the Jacobian computation for Adam and HIG
Optimizer
Adam
Adam
Adam
HIG
HIG
HIG
Batch size
256
512
1024
32
64
128
Memory (MB)
11.1
22.2
44.5
169
676
2640
Update duration (sec)
0.081
0.081
0.081
0.087
0.097
0.6
Jacobian duration (sec)
0.070
0.070
0.070
0.070
0.070
0.070
from the |-1, 1] interval, and a batch size of 256. For the optimizers, the following hyperparameters were used for both the well-conditioned loss and the ill-conditioned loss: Adam rj 0.3; GN has no learning rate (equivalent to 1) = 1), r = 10 10-4; HIG 1) = 1.0, T = 10-6
B.2 CONTROL OF NONLINEAR OSCILLATORS (SECTION 3.1
The Hamiltonian function given in equation leads to the following equations of motions:
x -2 + 4a(2; - x2-1)² 4a(x; 11+1)² u(t)
(20)
The simulations of the nonlinear oscillators of 12 units with a time step Al = 0.125. This results in 96 time steps via 4th order Runge- Kutta per learning iteration. We generated 4096 data points for control vector (0.0, 3.0), and an interaction strength a 1.0 with randomized conjugate variables x and The test set consists of 4096 new data points, For the neural network, we set up fully-connected network with ReLU activations passing inputs through three hidden layers with 20 neurons in each layer before being mapped 96 output layer with linear activation.
performed for two points and time interval
For the comparison with other optimizers (figure 2p) we performed broad hyperparameter search for each method, as outlined above, to determine suitable settings. The parameters for Adagrad (Duchi et al. 2011), Adadelta (Zeiler, 2012), Adam (Kingma & Ba (015)) RMSprop (Hinton et al.) 2012], Gauss-Newton (Gill & Murray 1978). HIGs, and stochastic gradient descent (Curry 1944) summarized in table For figure 2 the following hyperparameters were used: 7) 3 10 for Adam, and = 1.0, r 10 for HIG.
Further Experiments. Figure and figure contain additional hyperparameters for the method comparison of figure 2 b in the main paper. The graphs illustrate that all five method do not change their behavior significantly for the different batch sizes each plot, but become noticeably unstable for larger learning rates (plots on the right sides of each section).
with different
Details on the memory footprint and update durations can found in table Since our simulations were not limited by memory, we used an implementation for the Jacobian computation of HIGs, which scales quadratically in the batch size. Should this become a bottleneck, this scaling could potentially be made linear by exploiting that the Jacobian of the physical solver for multiple data points is blockdiagonal.
B.3 POISSON PROBLEM (SECTION 3.2
We discretize Poisson's 's equation on regular grid for two-dimensional domain Ω [0, 8] x [0.8] with a grid spacing of Az 1. Dirichlet boundary conditions of ø 0 are imposed on all four sides of Q. The Laplace operator is discretized with finite difference stencil Ames (214))
16

---

Published as a conference paper at ICLR 2022
Adadelta
b) Adagrad
10
0=0.00
10"h
6.00.1
10'
0.5
10t
0.001
10
0.00000
10'
n=
10-
10
10
10
10-2
10-2
Loss
Loss
10-4
15:00 3000
0
1500 3000
0
1500 2000
1500 2000
0
1500 3000
1500 2000
Wall clock time [sec]
Wall clock time [sec]
c) SGD
d) RMSprop
8-
100
10a
n-
105
n- 1e:05
10a
G- 0.0001
10*
G- 0.003
b-2024
b-2024
10-2
Loss
Loss
1900
1560 Wall clock time [sec]
1300
1560 Wall clock time [sec]
Figure 5: Control of nonlinear oscillators: Additional experiments with (a) Adadelta, (b) Adagrad, (c) stochastic gradient descent, and (d) RMSprop. Each showing different learning rates 1 and batch sizes b.
10a
Se-05
10a
n= 0.0001
7= 0.0003
10a
n= 0.0005
10
n= 0.001
10P
n= 0.003
b-2
10-2
10-2
10-2
10-2
10-2
Loss
10a1
10-4
10
10-1
10-4
10
L/b
0
1500
3000
0
1500
3000
1500 3000
1500
3000
1500
3000
0
1500
3000
Wall clock time [sec]
Figure 6: Control of nonlinear oscillators: Additional experiments with Adam for different learning rates 17 and batch sizes b.
For the neural network, we set up a fully-connected network with tanh activation functions. The 8x8 inputs pass through three hidden layers with 64, 256 and 64 neurons, respectively, before being mapped to 8x8 in the output layer. For training, source distributions are sampled from random frequencies in Fourier space, and transformed real space via the inverse Fourier transform. The mean value is normalized to zero. We sample data on-the-ftt fly, resulting in an effectively infinite data set. This makes a separate test set redundant all training data is previously unseen.
Further Experiments. Figure 7 shows Adam and HIG runs from figure over epochs. The HIG runs converge faster per iteration, which indicates that HIGs perform qualitatively better updates.
Additionally, the pretrained HIG run from figure 3 as starting point for further Adam training. The results are shown in Zp. We observe that the network quickly looses the progress the HIGs have made, and continues with loss value similar to the orginal Adam This again
17

---

Published as a conference paper at ICLR 2022
Table 4: Poisson problem: memory requirements, update duration and duration of the Jacobian computation for Adam and HIG
Optimizer
Adam
HIG
Batch size
64
64
Memory (MB)
1.3
3560
Update duration (sec)
0.011
13.8
Jacobian duration (sec)
0.010
0.0035
b)
Adam
10a
4.10-2
pretrained HIG
-0
Loss
Loss
2.
10a
Adam
10-2
Adam 3=
Adam 1o-04
10-2
HIG
6 10-3
10
103 Epochs
10s
103
1.2 1.3- Wall clock time [sec]
1.1-105
Figure 7: Poisson problem: a) Loss curves for Adam and HIG per epoch for different learning rates, b) Loss curves of Adam (n) =le-04) of HIG (1) = 0.02) pretrained with Adam, and Adam (1 1ee-04 pretrained with the HIGs.
supports our intuition that Adam, in contrast to HIGs, cannot harness the full potential of the physics solver.
Details on the memory footprint and update durations can be found in table
B.4 QUANTUM DIPOLE (SECTION 3.3
For the quantum dipole problem, we discretize the Schrödinger equation on a spatial domain [0, 2] with spacing of Ar 0.133 resulting in 16 discretization points. We simulate up to time of 19.2 with time step of At = 0.05, which yields 384 time steps. Spatial and temporal discretization use modified Crank-Nicolson scheme (Winckel et al., 20099) which is tailored to quantum simulations. The training data set consists of 1024 randomized superpositions of the first set of 1024 randomized superpositions.
and second excited state, while the test set contains For the neural network, we set up a fully-connected network with tanh activations passing the inputs through three hidden layers with 20 neurons in each layer before being mapped 384 neuron output layer with linear activation. Overall, the network contains 9484 trainable parameters.
Experimental details. For the training runs figure Adam used i 256, while for HIG 16, and 10 were used. For the training runs in figure 4.. Adam used 256, 17 0.0001, while HIGs used = 16. 7 10-5 5, and = 0.5. Details on the memory footprint and update
durations can be found in table 5
Figure and figure 9 show the performance both methods for a broader range of T settings for HIGs, and 77 for Adam. For Adam, trade-off between slow convergence and oscillating updates exists. The HIGs yield high accuracy in training across a wide range values for т, ranging from 10-5 to 10-3 This supports the argumentation in the main text that the truncation is not
18

---

Published as a conference paper at ICLR 2022
5e-05
10"
f1-
0.00001
11 0.0002
11- 0.0005
/1 0.001
10*
G 0.0002
10
/1 0.005
10
10-3
10
10-2
10-2
10-2
10-2
Loss
10 10
10
10
2-64
10
1-8
10.4
10
10
10
10-4
10"5
2500
5000
2500
5000
2500
5000
2500 5000
2500 5000
2500 5000
2500 5000
Wall clock time [sec]
Figure 8: Quantum dipole: Additional experiments with Adam for different learning rates and batch sizes b.
r=
1c:08
r=
T= 16:06
= 1e:05
== ==0.0001
T 0.0010
r= 0.01
10a
109
10"
10
10-2
10-2
10-2
10-2
10-2
10-2
Loss
0-0.1
10
10-4
10-4
10:4
b-8
0.0.5
15
10.
10
10
10
10
2500
5000
25005000
2500
5000
2500 5000
2500 5000
0
25005000 5000
0
25005000
Wall clock time [sec]
Figure 9: Quantum dipole: Additional experiments with HIGs for different learning rates 1), batch sizes and truncation parameters
overly critical for HIGs. As long as numerical noise is suppressed with r > 10 and the actual information about scaling of network parameters and physical variables is not cut off. The latter case is visible for an overly large
0.01 the last graph on the right.
Note that many graphs in figure 9 contain small plateau at the start of each training run. These regions with relatively small progress per wall clock time are caused by the initialization overhead of the underlying deep learning framework (TensorFlow in our case). As all graphs measure wall clock time, we include the initialization overhead of TensorFlow, which causes a noticeable slow down of the first iteration. Hence, the relatively slow convergence of the very first steps in figure 9 are not caused by conceptual issues with the HIGs themselves. Rather, they are a result of the software frameworks and could, e.g., be alleviated with pre-compilation of the training graphs. In contrast, the initial convergence plateaus of Adam with smaller y) in Figure & are a fundamentally different nature: they caused by inherent problem of non-inverting optimizers: their inability to appropriately handle the combination of large and small scale components in the physics of the quantum dipole setup (as outlined in section 3.3
Table 5: Quantum dipole: memory requirements, update duration and duration of the Jacobian computation for Adam and HIG
Optimizer
Adam
Adam
Adam
HIG
HIG
Batch size
256
512
1024
8
16
Memory (MB)
460
947
2007
1064
5039
Update duration (sec)
0.40
0.50
1.33
0.42
0.60
Jacobian duration (sec)
0.39
0.49
1.32
0.40
0.53
19

---

Published as a conference paper at ICLR 2022
5e-05
f1 0.00001
( 0.0002
11= 0.0005
/100.00
Q 1000002
0- 0.005
10c
10*
10a
1-64
b-128
10
10-3
10a3
10-2
10-2
10
10-2
Loss
10
10
10
10.
10
10
10
10
2500
5000
2500 5000
2500
5000
2500 s000
2500 5000
2500
5000
2500 5000
Wall clock time [sec]
Figure 10: Quantum dipole with Convolutional Neural Network: Experiments with Adam for different learning rates 7) and batch sizes b.
1e-00
T=1e00
1e:05
T=0.0001
0=00001
T=0.01
10P
109
109
109
107
10F
109
b-8
10
10
10
10
10
10
0.3.3
Loss
b:B
8-8
10
10-4
10-4
10-1
10
10
0-0.5
10
b-15
0-0.5
10-6
10-4
10-4
10-1
10-4
10
000.
2500
5000
2500
5000
2500 soon
2500 soon
2500
2500 5000
2500 5000
Wall clock time [sec]
Figure 11: Quantum dipole with Convolutional Neural Network: Experiments with HIGs for different learning rates 77. batch sizes b, and truncation parameters т
Loss Functions. While training evaluated terms of the regular inner product as loss function: L(W., Tb - (.. 1.) we the following modified losses to evaluate low- and high- energy states for figure 4b. Let Y be the first excited state, then we define the low-energy loss
L(V Wb) = (1(M,. W1) - ( w
Correspondingly, define the high-energy loss with the second excited state V2::
L( Mo) = (1(Va. M2) - 20 2a))
Additional Experiments with Convolutional Neural Network. Our method is agnostic to specific network architectures. To illustrate this, we conduct additional experiments with convolutional neural network. The setup is the same as before, only the fully-connected neural network is replaced by network with hidden convolutional layers each with kernel size 3, 20 features and tanh activation, followed by 384 neuron dense output layer with linear activation giving the network a total 21984 trainable parameters.
The results of these experiments plotted in figure 10 and We find that HIGs behave in line with the fully-com collecccodeeted network case (figure 9). There exists a range r-values from around 10 10-5 to 10 for which stable training is possible. Regarding optimization with Adam, likewise observe a faster and more accurate minimization the loss function for the best HIG run (n = 0.7. 16, 10-4) compared to the best Adam (1) 0.0002, b 256).
20

---

Published as a conference paper at ICLR 2022
C ABLATION STUDY
In this last section, we investigate how the HIG-hyperparameters affect the outcome. This includes ablation experiments with respect to & and т defined section 2.2 We use the nonlinear oscillator example the basis for these comparisons and consider the following HIG update step:
O0y 20
<B,K>
T
(1 B,K) -1
(21)
Here, the exponent <. B > of the Jacobian denotes the following procedure defined with the aid of the singular value decomposition J UAV
as:
J<B,>> m VA"UT,
(22)
Compared to the HIG update in the main text, update 21 has an additional scalar prefactor with an parameter 3 resulting from earlier experiments with our method. Setting 3 -1 k yields algorithms that rescale the largest singular value 1, which ensures that the resulting updates cannot produce arbitrarily large updates in yyspace. This can be thought of as weaker form of scale
invariance. Just as equation 21 defines an interpolation between gradient descent (3 = 0, A = 1) and the Gauss-Newton method (3 = 0, 11 as well.
Scalar prefactor term 3: We test B-values between 0, no scale correction, and 100--0. -0.5, which fully normalizes the effect of the largest singular value for x = -0.5. The results are shown in figure [2p. Compared to the other hyperparameters, we observe that / has only little influence the outcome, which is why we decided to present the method without this parameter in the main text.
Exponent of the diagonal singular value matrix x: We test k. for various values between 1.0, stochastic gradient descent, and 11, Gauss-Newton. The results shown in figure 12p. For positive values, curves stagnate early, while for negative k, the final loss values are several orders of magnitude better. The HIG curve corresponding to supports our argumentation that strong dependence on this parameter exists, and that a choice of -0.5 is indeed good compromise for correcinng updates of reasonable size. The strong improvement as soon as s becomes negative indicates that the collective inversion of the feedback of different data points of the mini-batch an important ingredient in our method.
-0.5 achieves the best result. This
understand the effect of this parameter, we consider the singular value decomposition (SVD) of the network-solver Jacobian, which is determined by the SVDs of the network Jacobian and the solver Jacobian. The singular values of a matrix product AB depend non-trivially the singular values of the matrices A and B. In the simplest case, the singular values of the matrix product are received by multiplication of the individual singular values of both matrix factors. In the general case, this depends each other. However, it likely that singular vectors with a small singular value of A or B overlap significantly with singular vectors with a small singular value of AB. For this reason, it is important not to truncate too much as this might remove the small-scale physics modes that we are ultimately trying preserve in order to achieve accurate results. On the other hand, less truncation leads to large updates of network weights on a scale beyond the validation of the linear approximation by first-order derivatives. These uncontrolled network modifications can lead to satursate neurons and prevent further training progress.
Truncation parameter T:
how the singular vectors of A and B overlap with
From a practical point of view, we choose according to the accuracy of the pure physics optimization problem without neural network. For the quantum dipole training, this value was set to 10-5 Trying to solve the pure physics optimization with far smaller values leads to a worse shows result or convergence at all. The network training behaves in line with this: Figure Shows that the network does not learn to control the quantum system with 7-values far smaller than 10 For the nonlinear oscillator system, the pure physics optimization is stable over large range of 7-values with similarly good results. For the network training, we chose J to be 10 10-6, We conducted further experiments for the network training with different from 10-5 10-10 presented in figure 13
21

---

Published as a conference paper at ICLR 2022
10-1
==
-0.25
b)
10a
-0.25
0.0
-0.333
1 (SGD)
-0.5 (HIG)
-0.4
0.5
-0.75
-0.166
-0.5
-1 (GN)
Loss
10-1
Loss
10
10-6
10-5
1500
3000
1500
3000
Wall clock time [sec]
Wall clock time [sec]
Figure 12: a) Ablation experiments with the 3-hyperparameter, and b) with the k-hyperparameter.
Gauss-Newton
HIG
109
10
1e-05
1e-05
t = 1e-06
= 1e-06
1e-07
1e-07
1e-08
1e-08
10-2
10
1e-09
1e-10
Loss
10-4
10-4
10-6
10-6
1500
3000
0
1500
3000
Wall clock time [sec]
Wall clock time [sec]
Figure 13: Ablation experiments with the r-hyperparameter.
which show that HIGs have similar tolerance in т. For a comparison, we also plotted Gauss- Newton curves for different т. We observe that GN curves become more unstable for smaller truncation values 7 and diverge in the case 10-9 and 10 while HIG curves achieve overall better loss values and start to converge in this parameter.
22


Reviews:
[{"review_id": "v77Tqyfeer9", "reviewer": "Reviewer_itKt", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Strengths:\nThe paper presents an efficient method for training neural networks for physical simulations. The idea is simple and is easy to understand. The written and organization of paper is clear and easy to follow.\n\nWeaknesses:\nI have some concerns on the experiments.\n1. It is better to compare the batch size of HIGs with the Gauss-Newton method.\n2. The truncation threshold \\tau is important for both accuracy and efficiency. It is better to discuss its effect in more details. Also, how to set \\tau for a problem? \n3. What does \\gamma refer to throughout our the paper?\n4. Just Below Eq.(2), \"f\" should be \"$f$\".", "summary_of_the_paper": "The paper considers the training of neural networks for physical simulations. By distributing the burden equally between network and physics, the paper presents the half-inverse gradients (HIGs) method. Experiments show its advantage for achieving a faster and more accurate minimization.", "summary_of_the_review": "The paper is well-written and organized. The idea is simple and effective. It is better to discuss more details on the settings of the hyper-parameters in the proposed method.   ", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "dgcg9IxlD5O", "reviewer": "Reviewer_6S7S", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The paper starts out with a very strong overview of past work, and cites Holl et al. (2021) for the motivation that existing optimizers do not perform well on joint optimization of NNs and physics. It would be very helpful if the paper was a little more self contained in the motivation. I see several references to \"physics solvers\" and \"ill-conditioning\" but no real motivation as to why these are related. The authors give an example of an ill-conditioned quadratic, which helps motivate why incorporating neural networks with ill-conditioned problems could lead to difficulty, but it's not clear why using \"physics solvers\" causes this.\n\nFurthermore, I would like to see \"physics solvers\" to be defined a little bit more explicitly in the beginning of the paper. I was not able to quite understand until reading the appendix and all the case studies. I now assume that you mean a time-stepping scheme that you can differentiate through end-to-end (or using adjoint methods).\n\nSection 2:\n- It would be very helpful if you made vectors more explicit in the notation.\n- e.g. In the \"Physics optimization\" subsection, you say \"the sum reduces to a single data point\". This confused me for a bit, until I realized $x$ and $y$ are vectors (but even then MSE would reduce to a sum of squared terms over the dimensions of the states).\n- Equation (3) would be helpful to emphasize that you're using a pseudoinverse, not a standard inverse.\n- It took me some time to understand that batch dimension is the number of points in the target state.\n- Did you experiment with different configurations of $\\kappa, \\beta, \\tau$? It would be interesting as ablation experiments to disentangle the benefits of these parameters. In particular, does most of the benefit come from $\\beta$? Or is there also a benefit in setting $\\kappa$ to an intermediate value?\n\nNonlinear oscillator:\n- It might be useful for people who are less familiar to write out the actual time-stepping you get with Hamilton's equations. Maybe this can go in the appendix.\n- How many Adam steps are you able to perform in one HIG step? It would also be useful to see number of steps required vs batch size.\n\nPoisson problem:\n- Correct me if I'm wrong, but it seems like this setup is different from the other ones in the sense that it's a surrogate training problem, rather than optimizing through a physics solver. Is that correct? I think it would be helpful to make this clear.\n- The pretraining is interesting. Would GN perform well when pretrained with Adam?\n- re: \"the Poisson problem is relatively simple, requiring only a single matrix inversion\": Arbitrarily ill-conditioned quadratics can also be solved by a single matrix inversion. I'm not convinced that is the reason for the observations.", "summary_of_the_paper": "The paper proposes a way to interpolate between gradient descent and Gauss-Newton's method for solving nonlinear least squares problems arising from physics informed training. The authors cite past work that physics informed training is often ill conditioned, so gradient descent often performs poorly. They give a classic example of ill-conditioning, and show that gradient descent converges slowly, while Gauss-Newton quickly saturates neuron activations.\n\nThis motivates them to introduce Half-Inverse Gradients (HIG), which interpolates between GD and GN using the SVD. The authors then try the method on several test problems in scientific computing: Control of nonlinear oscillators, Poisson, and control of a Quantum dipole. They compare Adam, HIG, and GN.", "summary_of_the_review": "The paper proposes a potentially interesting method to interpolate between GD and GN, and is certainly moving in the right direction of a very important problem in scientific machine learning, but the paper as it is right now is a little confusing to read. It would also benefit a lot from some ablation experiments showing the benefit of each individual part of their methods. I vote for marginally below threshold, but am willing to change my score if my concerns are addressed.\n\nLet me know if you have questions about my review.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "JpkAHGkH-hR", "reviewer": "Reviewer_Khvk", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "I personally feel excited about this direction. Data-driven physical learning has been a hot area recently where progress is being made in multiple parallel paths including problem formulation, network architectures, loss specifications, etc. The optimization paradigm has been an unavoidable component in deep learning, and this work fits into the place as the analogy of its counterparts in vision/NLP problems. This manuscript has shown convincing theoretical intuition and validated on valuable physical problems, nevertheless, it arouses my curiosity in several directions:\n\n1. I found this work strongly parallel with a recent pre-print [1], which I would like to refer to as a concurrent work to this manuscript. I am well aware that a comparison is not required but would like to hear more qualitative discussion.\n2. It seems to me that the loss function in the toy example should be $l(y,\\hat{y},\\lambda)=\\frac{1}{2}(y^1-\\hat{y}^1)^2+\\frac{1}{2}\\lambda(y^2-\\hat{y}^2)^2$, or did I miss something? Also the notation (superscripts and subscripts) are confusing.\n3. The computational cost seems high. The Jacobian inference, in spite of being vectorized for a subset of modules, can be a memory monster. I expect some results and discussions on the memory footprint. Also, I wonder if the method can benefit from the Krylov subspace method where an iterative method is used for speedup.\n4. My experience with Adam tells me learning rate scheduler is a good friend of first-order methods. I wonder if one is used for Adam baselines in the experiments. If not, how does it perform?\n5. A very common trick in training neural networks is adding normalization layers (in recurrent cases, layer normalization is more popular). How does this method deal with normalization layers? Does it require it? If it does, how to update it during training?\n\n[1] Holl, Philipp, Vladlen Koltun, and Nils Thuerey. \"Physical Gradients for Deep Learning.\" arXiv preprint arXiv:2109.15048 (2021).", "summary_of_the_paper": "This manuscript introduces a new optimization fashion bridging the Gauss-Newton method and the vanilla gradient descent method in the middle for physical optimization with neural networks. Traditionally, when a neural network is used for a physical problem, it will inevitably be affected by the unbalanced magnitudes dramatically. This manuscript proposes HIG, which is a middle point of two popular and individually advantageous optimization methods. It also provides a nice guideline on how to choose the hyperparameters and examine the efficacy of HIG on a concept-illustrating synthetic problem and three realistic physical problems.", "summary_of_the_review": "This manuscript fills a blank for optimization in physical deep learning with a nice and clean method. The presentation is relatively clear to me. I would recommend an accept.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}]

Author rebuttal:

Meta-review: Thank you for your submission to ICLR.  The reviewers and I are in agreement that the paper presents a substantial contribution to the field at the intersection of differentiable simulation and ML methods.  In particular, the half-inverse method is compelling, non-obvious, and hints of a nice path forward towards the goal of practical differentiable simulations within models.  Overall I'm happy to recommend the paper be accepted.
Decision: Accept (Spotlight)