Stochastic Taylor Derivative Estimator: Efficient amortization for arbitrary differential operators

Optimizing neural networks with loss that contain high-dimensional and high-order differential operators is expensive to evaluate with back-propagation due to $\mathcal{O}(d^{k})$ scaling of the derivative tensor size and the $\mathcal{O}(2^{k-1}L)$ scaling in the computation graph, where $d$ is the dimension of the domain, $L$ is the number of ops in the forward computation graph, and $k$ is the derivative order. In previous works, the polynomial scaling in $d$ was addressed by amortizing the computation over the optimization process via randomization. Separately, the exponential scaling in $k$ for univariate functions ($d=1$) was addressed with high-order auto-differentiation (AD). In this work, we show how to efficiently perform arbitrary contraction of the derivative tensor of arbitrary order for multivariate functions, by properly constructing the input tangents to univariate high-order AD, which can be used to efficiently randomize any differential operator. When applied to Physics-Informed Neural Networks (PINNs), our method provides>1000$\times$ speed-up and>30$\times$ memory reduction over randomization with first-order AD, and we can now solve \emph{1-million-dimensional PDEs in 8 minutes on a single NVIDIA A100 GPU}. This work opens the possibility of using high-order differential operators in large-scale problems.

Paper

Similar papers

Peer review

Reviewer cUDy7/10 · confidence 4/52024-06-22

Summary

The paper proposes a new method, STDE, to compute high-order differential operators of (high-dimensional) neural network-represented functions. Two key ingredients include a generalization of Taylor-mode AD and randomness in the algorithm. This work shows impressive performance of STDE, in terms of both speed and memory.

Strengths

- This work addresses an important and challenging task in physics-informed machine learning. - The illustration and analysis of previous AD algorithms and STDE are clear and comprehensive. - STDE displays noticeable potential to improve the computation of arbitrary differential operators in both memory consumption and efficiency. - The work provides new insight into introducing randomness other than SDGD or HTE whose basic ideas are employing Monte Carlo over dimensionality.

Weaknesses

1. Compared with original back-propagation or forward-Laplacian(FL), STDE is a random algorithm that does not give **precise** derivatives. In precision-demanding tasks, e.g. solving many-electron Schrödinger equations, the inaccurate value of differential operators will possibly result in an unreliable solution. To make the result more convincing, a cost-accuracy curve (STDE with different number of random samplings, i.e. $|J|$ in Eqn.(17)) should be added. Furthermore, in light of this, the comparison between STDE and FL in super high dimensions (in table 1,2) is unfair. A possible remedy is to compare STDE with modified FL, which employs randomness to be more scalable, for instance, Monte Carlo over dimensions. The author can refer to [1] for more technical and implementation details of FL. 2. All the test cases of high-dimensional high-order operators are Laplacians (including $\Delta^2$). The high-order diagonal differential operator is mentioned in L216 but does not appear in the experiments. A comparison between STDE and direct back-propagation is expected. The applications are all coming from the small area of PINN. To demonstrate the broad impact of this work, the author should target more tasks of general interest, e.g. second order derivative of NN parameter, which plays an important role in accelerating optimizations. 3. Some of the highly related existing works are missing. E.g., - Regarding L63, the forward rule for general different operators is done in [1]. - Regarding the geometric interpretation of Taylor-mode AD in sec 3.4, although not explicitly mentioned in Ref. [6] in the paper, is discussed in section 2 of [2]. The author should cite [6] or the first paper presenting similar ideas, and clarify their novelty if they still hope to claim generalizing univariate Taylor mode AD as their contribution. 4. As is discussed in Section 4.4, STDE can not be applied to *arbitrary* high-order operators. Obviously, this restricts the potential application scenarios of STDE. 5. The author didn’t discuss how large $l$ would be for large $k$ (notations from L165), which is necessary to assess the effectiveness of STDE for general high-order derivatives. 6. Some parts of the writing are confusing. See **Questions**. Also, there are several misuses of terminologies. E.g., - Fréchet derivative, not Frechet. - The form of second order PDE in Eq.(18) is not ‘general’. For instance, you are assuming ellipticity in the second-order term. Please refer to Eq.(10) or Eq.(1) in [1] for *general* second-order operator. - In Appendix I.1 and the main text, ‘inseparable and effectively high-dimensional’, I guess what the authors intend to say is ‘elliptic’. The linear part of all the test equations is merely a Laplacian, casting doubt on the effectiveness of STDE for general elliptic equations. - In Appendix I.2 and the main text, these PDEs are not nonlinear Fokker-Planck Equations, I think the correct name would be ‘semilinear parabolic equation’. [1] Li R, Wang C, Ye H, He D, Wang L. DOF: Accelerating High-order Differential Operators with Forward Propagation. ICLR 2024 Workshop on AI4Differential Equations In Science [2]Tan, Songchen. Higher-Order Automatic Differentiation and Its Applications. Diss. Massachusetts Institute of Technology, 2023.

Questions

1. Which function does ‘univariate’ refer to in Sec 3.4, $g$ or $F$? What is the first contribution (in the introduction), generalizing Taylor mode AD to the case where $g$ is multivariate, or $F$? I am skeptical of the authors’ claiming it as a contribution/novelty because it seems like this is merely about replacing scalar Faà di Bruno's formula with that of tensors. 2. In Appendix D, where $g$ defines an $n$-dimensional manifold, is $n$ always set to 1 in this work? If not, please point me to where the method with general $n$ is discussed. 3. Although the authors claim to remove the exponential complexity w.r.t. the order of derivatives, for general operators the size of coefficient tensor **${C}$** in Eqn(10) grows exponentially and the sampling procedure in Eqn(15) will become hard. 4. Could the authors think of any real application scenarios with high-order high-dim operators? As for the case of low-dim high-order operators, the authors need to check if the differential operator takes up a considerable proportion of the whole process, e.g. training PINN. If not, the benefit of accelerating differential operators might be very limited. 5. Does STDE have the potential to make use of sparsity in the neural network function, as shown in FL and [1]? For the specific case of second-order operators, how would STDE perform when $C$ is dense but low-rank? 6. Instead of tackling a scalar function (e.g. Laplacian of a function), could STDE help with the computation of tensor-value functions, e.g. Hessian or curvature? Hessian appears in Monge-Ampère equation and Ricci curvature appears in Einstein equation. On writing and notations: 1. Can you clarify on ‘$[2]$’ in L214 with more details or concrete examples? For me, $d^2u_\theta$ in Eqn.(16) is a (0,2)-tensor (Hessian), but why does it take three input $(a,e_j,0)$? Please also clarify on ‘$[3]$’ in Eqn(44). 2. In Eqn.(19), it should be $\sum_{i=1}^d$. 3. In L247, do you mean $v\sim p$,, instead of $(a,v,0)\sim p$? 4. In the ‘100D’ column of table 2, 539 is not the smallest number. Why is it in bold? 5. There needs to be more description in section F.1. It is good to remind the readers that Einstein notation is employed in Eq.(45) and that the indices $n,\ n’,\ n’’$ are omitted. The author should inform the readers that in the discussion $n=1$ for $n$ in Eq(40) and notations like $\partial^4u(x)(e_i,e_j,0,0)$ refers to setting $v^{(1)}=e_i,\ v^{(2)}=e_j$ instead of the taking $e_i,e_j,0,0$ as inputs for a (0-4) tensor. Overall, this paper has the potential to be a good work. Hope the concerns above can be well addressed and then I will upgrade my assessment accordingly.

Rating

7

Confidence

4

Soundness

3

Presentation

3

Contribution

3

Limitations

Yes.

Reviewer SkVJ7/10 · confidence 4/52024-06-30

Summary

This paper proposes a stochastic optimization approach to find the minimizer of a cost function, which involves the complicated differential operators. The problem is computationally complex, hence, the authors propose to deal with a minibatch of derivatives in each iteration which reduces computational complexity. The method has application in various learning problems involving complicated differential operators such as physics-informed neural networks. Experiments are provided to evaluate the method.

Strengths

The paper seems to address a difficult optimization problem. The dimension reduction idea seems to be new and interesting. Application to PINN seems also novel.

Weaknesses

NA

Questions

NA

Rating

7

Confidence

4

Soundness

4

Presentation

3

Contribution

4

Limitations

NA

Reviewer yjY87/10 · confidence 3/52024-07-12

Summary

This work proposes a scalable method for the optimization of loss functions including higher-order derivatives. The proposed method interprets arbitrary differential operators as derivative tensor contractions which are then estimated through random contractions. These random contractions can be computed efficiently using Taylor mode AD.

Strengths

1. The proposed method is novel and interesting. I think the idea of random contractions of the derivative tensor is very intuitive and at least in principle should be scalable and a clear improvement over SDGD by ameliorating the exponential scaling in the order of the the differential operator. 2. The idea is technically sound. There is sufficient empirical validation. 2. The well-written presentation of the background material and JAX implementations of SDGD are also useful secondary contributions.

Weaknesses

1. The presentation could be improved a tad bit. I think having 4 pages of background material is a bit unnecessary specially the general background on automatic differentiation could easily be moved to the appendix and space could be better utilized by explaining the method in more detail and discussing the experimental details. For example, a schematic diagram similar to Figure 2 illustrating STDE could be helpful in understanding the proposed method.

Questions

1. In 291 it is stated "Since the only change is how the derivatives are computed, the relative L2 error is expected to be of the same order among different methods". Is this accurate? SDGD and STDE are stochastic approximations whereas Forward Laplace is not. 2. I am a bit unclear about how to choose a distribution over the l-jets that satisfy the unbiasedness condition for any specific problem. It seems that the sparse random jets have the advantage of universally applicable but they also seem to involve a lot of redundant computations.

Rating

7

Confidence

3

Soundness

3

Presentation

2

Contribution

3

Limitations

Some limitations are discussed in section 6.

Reviewer xY8N8/10 · confidence 5/52024-07-24

Summary

The paper addresses the computational challenges of optimizing neural networks with loss functions that include high-dimensional and high-order differential operators. These challenges arise due to the scaling of the derivative tensor size with the dimension of the domain (d) and the computational graph's size with the number of operations (L) and the order of the derivative (k). Traditional methods either amortize the computational cost over the optimization process via randomization or use high-order auto-differentiation (AD) for univariate functions to tackle these issues. In this work, the authors propose a method to efficiently perform arbitrary contraction of the derivative tensor for multivariate functions. This is achieved by constructing input tangents to univariate high-order AD, enabling efficient randomization of any differential operator. When applied to Physics-Informed Neural Networks (PINNs), this approach provides significant speed and memory efficiency improvements, achieving over 1000 times speed-up and 30 times memory reduction compared to randomization with first-order AD. The method allows solving 1-million-dimensional partial differential equations (PDEs) in just 8 minutes on a single NVIDIA A100 GPU, opening up the possibility of using high-order differential operators in large-scale problems.

Strengths

The paper introduces STDE, a general method for constructing stochastic estimators for arbitrary differential operators, which can be efficiently evaluated using Taylor mode auto-differentiation (AD). When evaluated on Physics-Informed Neural Networks (PINNs), a specific optimization problem where the loss function includes differential operators, STDE significantly outperforms baseline methods. Furthermore, STDE's applicability extends beyond PINNs to arbitrarily high-order and high-dimensional AD-based PDE solvers, making it more general than related methods. The strengths of this paper are as follows: Generality: STDE can be applied to a wide range of problems, including those involving arbitrarily high-order and high-dimensional differential operators. This broad applicability distinguishes STDE from other methods, which are often restricted to specific forms of second-order PDEs. Efficiency: The method enables efficient evaluation of stochastic estimators through Taylor mode AD, providing significant computational benefits. Performance: In practical evaluations on PINNs, STDE outperforms baseline methods in terms of both speed and memory efficiency. It demonstrates over 1000 times speed-up and 30 times memory reduction compared to first-order AD randomization. Scalability: STDE allows for the solution of extremely large-scale problems, such as 1-million-dimensional PDEs, in a matter of minutes on advanced hardware like the NVIDIA A100 GPU. Versatility: Beyond PINNs, STDE can be applied to various AD-based PDE solvers, making it a versatile tool for tackling a broad spectrum of differential operator-based optimization problems. Overall, STDE's generality, efficiency, performance, scalability, and versatility make it a powerful method for addressing high-dimensional and high-order differential operator challenges in neural network optimization.

Weaknesses

While the paper demonstrates the significant strengths and broad applicability of STDE, it is important to acknowledge some limitations and areas for future improvement. As a general method, STDE may not leverage the specific optimization possibilities that are available for particular operators. Additionally, the study did not explore variance reduction techniques, which could potentially enhance the method's performance and could be a promising area for future research. Another observation is that while reducing the randomization batch size improves both the speed and memory profile of STDE, this comes with a trade-off in the form of increased computational variance. Further analysis is required to understand and optimize this balance between computational efficiency and variance. Looking ahead, the paper identifies an intriguing connection between the fields of automatic differentiation (AD) and randomized numerical linear algebra, highlighting the potential for future work at this intersection. Such research could lead to significant advancements in large-scale scientific modeling with neural networks. In summary, while there are areas for refinement, the contributions of this paper are substantial. The development of STDE as a general and efficient method for constructing stochastic estimators for arbitrary differential operators is a notable achievement, offering substantial benefits for high-dimensional and high-order differential operator problems in neural network optimization. The identified limitations and future research directions provide a clear path for further enhancing this already impressive work.

Questions

Given the importance of various other complex equations in scientific modeling, I am curious about the applicability of STDE to equations such as the Nonlinear Schrödinger Equation (NLS), the fourth-order NLS, and the Navier-Stokes equations. Could you elaborate on how STDE might perform or be adapted for these specific cases? Additionally, are there any preliminary results or theoretical considerations you could share regarding the application of STDE to these important equations?

Rating

8

Confidence

5

Soundness

4

Presentation

4

Contribution

4

Limitations

Yes, they have.

Reviewer yjY82024-08-09

Keep Score

Thanks for the clarification! Overall I like the paper and intend to keep the score.

Reviewer cUDy2024-08-09

Thanks for your reply. Regarding the ablation study on randomization batchsize you added, I don’t think it is very appropriate since there is also randomness in the optimization. I would instead suggest directly comparing the accuracy on differential operator itself, e.g. Laplacian of a complicated network function. Additionally, do you have any intuition why batchsize=16&64 get stuck whilst both a smaller and a larger batchsize do not? I would suggest adding the clarifications you replied in the final version. Also, I want to comment on some notations and terminologies and do hope the author to make them precise given the diverse background of NeurIPS attendees. > However, they are also specifically driftless Fokker-Planck equations with semilinear extension. Fokker-Planck equation has a clear definition and derivation from evolution of probability densities. I could not agree with the usage of the terminology here and do think every word should be made precise (even though it is not relevant to the contribution of this work at all.) > I think these two are synonymous, as a here is $a$ not a R.V. Then you should use either $v\sim p$ or $(a,v,0)\sim \delta_a\times p \times \delta$ (or something like that).

Authorsrebuttal2024-08-11

> Regarding the ablation study on randomization batchsize you added, I don’t think it is very appropriate since there is also randomness in the optimization. Thank you for your comment on the additional experiment. Firstly we would like to emphasize that, one of the main claims of our paper is when solving optimization problems like equation (1), the cost of computing expensive derivative operators can be amortized over the stochastic optimization of the neural network ansatz. This is why we conduct all the experiments under a stochastic optimization. We regret that this point was not communicated more clearly, and we will try to make this point more prominent in the final version of the paper. It is possible to conduct yet another experiment that directly compares the accuracy of differential operators under different randomization batch sizes, and we do expect to see that variance decreases proportional to the randomization batch size. However, this does not determine the effectiveness of STDE as explained above, and we think it would be more interesting to show the effect of randomization batch size on the convergence of various PDEs. > Additionally, do you have any intuition why batchsize=16&64 get stuck whilst both a smaller and a larger batchsize do not? This was puzzling for us as well. Upon conducting further experiments, we found that it is most likely since we had chosen a set of particularly bad random seeds ({1,2,3,4,5}). We are rerunning the experiments and will report the new result once it is ready. > I would suggest adding the clarifications you replied in the final version. Yes, we will be revising the final version to incorporate the clarifications we made in the above response. We would like to thank the reviewer again for his/her detailed checking on the terminologies, which greatly improved the quality of our paper. > Fokker-Planck equation has a clear definition and derivation from evolution of probability densities. I could not agree with the usage of the terminology here and do think every word should be made precise (even though it is not relevant to the contribution of this work at all.) Upon further consideration, we think that it would indeed be more precise to rename nonlinear Fokker-Planck in Appendix I.2 to semilinear parabolic equations, considering the interpretation of the Fokker-Planck equation. Initially, we used the name Fokker-Planck to be consistent with the baseline method SDGD (Section 5.2 in [1]). Thanks again for checking the terminologies used in our paper! > Then you should use either $v\sim p$ or $\delta_{a}\times p\times \delta$ (or something like that). Thank you for your suggestion. We think both of the suggested modifications are more precise than the current formula used in the paper. We will incorporate this change into the final version of our paper.

Program Chairsdecision2024-09-25

Decision

Accept (oral)

© 2026 NYSGPT2525 LLC