Solving stiff dark matter equations via Jacobian Normalization with Physics-Informed Neural Networks

M. P. Bento; H. B. Câmara; J. R. Rocha; J. F. Seabra

hep-phFeb 25, 2026

Solving stiff dark matter equations via Jacobian Normalization with Physics-Informed Neural Networks

Authors: M. P. Bento, H. B. Câmara, J. R. Rocha, J. F. Seabra

Organizations: Departamento de Física and CFTP, Instituto Superior Técnico, Universidade de Lisboa, Av. Rovisco Pais 1, 1049-001 Lisboa, Portugal · Institute of Experimental and Applied Physics, Czech Technical University in Prague, Prague 160 00, Czech Republic

Abstract

Stiff differential equations pose a major challenge for Physics-Informed Neural Networks (PINNs), often causing poor convergence. We propose a simple, hyperparameter-free method to address stiffness by normalizing loss residuals with the Jacobian. We provide theoretical indications that Jacobian-based normalization can improve gradient descent and validate it on benchmark stiff ordinary differential equations. We then apply it to a realistic system: the stiff Boltzmann equations (BEs) governing weakly interacting massive particle (WIMP) dark matter (DM). Our approach achieves higher accuracy than attention mechanisms previously proposed for handling stiffness, recovering the full solution where prior methods fail. This is further demonstrated in an inverse problem with a single experimental data point - the observed DM relic density - where our inverse PINNs correctly infer the cross section that solves the BEs in both Standard and alternative cosmologies.

Figures & tables

Panel of Figure 1 — Figure 1: PINN results for ODE ( 8 ), with normalized (unnormalized) residuals in green (red). Left: Maximum Hessian eigenvalue $\lambda_{\text{max}}$ ( 9 ) in terms of the Jacobian $C$ . The circles (crosses) indicate normalized (unnormalized) residuals. The linear fit slope gives the power-law scaling exponent of $\lambda_{\text{max}}$ with $C$ . Right: Mean squared difference between PINN and FEM solutions ( 26 ) in terms of $C$ .

Panel of Figure 2 — Figure 2: PINN results for ODEs Eq. ( 27 ) (left) and ( 28 ) (right), with normalized [unnormalized] residuals in magenta [blue]. Top panels: Maximum Hessian eigenvalue $\lambda_{\text{max}}$ ( 9 ) in terms of $C$ . Bottom panels: Mean squared difference between PINN and FEM solutions $\overline{\varepsilon}$ ( 26 ) in terms of $C$ .

Panel of Figure 3 — Figure 3: PINN $95\%$ confidence intervals for $\varepsilon_{\text{rel}}(z)$ ( 43 ) with (magenta) and without (blue) Jacobian normalization, using linear (left) or negative sigmoid (right) output layers.

Panel of Figure 4 — Figure 4: Solution of the freeze-out DM BE obtained with FEM (black dashed line) and various PINN approaches: vanilla (blue), Jacobian normalization (magenta), residual-based attention (green), and soft attention (brown). Left (Right): Results without (with) a negative sigmoid activation in the output layer. Top: Evolution of the WIMP DM particle yield. The horizontal dotted line satisfies the CDM abundance ( 1 ). Middle: Mean relative squared difference ( 45 ) between PINN and FEM during training. Bottom: Comparison of residual weights.

Figure 5: Diagrammatic representation of the inverse PINN structure modeling the WIMP DM BE for a given cosmological scenario (see text for details).

Figure 6: Inverse PINN solution: particle yield evolution $W_{n}(z)$ for Standard (magenta), GB (green) and RS (brown) cosmology. We set $z_{t}=\ln 500$ . Black-dashed curves were obtained via FEM. The horizontal dotted line corresponds to the $W_{n}(z)$ value reproducing the observed DM abundance ( 1 ).

$\gamma$	$C_{\text{PINN}}$	$\langle\sigma v\rangle$ $\left(\text{GeV}^{-2}\right)$	$W(z_{f})_{\text{FEM}}$	$\Omega h^{2}_{\text{FEM}}$
$-2/3$	$27.6359$	$3.1174\times 10^{-10}$	$-26.1585$	$0.1200$
$0$	$29.2037$	$1.4951\times 10^{-9}$	$-26.1585$	$0.1200$
$2$	$31.7022$	$1.8187\times 10^{-8}$	$-26.1585$	$0.1200$

Table 1: Inverse PINN solution: finding

C

[see Eq. ( 37 )] for GB (

\gamma=-2/3

), Standard (

\gamma=0

) and RS (

\gamma=2

) cosmology –see Fig. 6 . Last two columns cross-check the PINN results with FEM.

Explore similar work

Variational Boosting for Physics-Informed Neural Networks
Jul 27, 2026Kaylee Vo, Pavlos ProtopapasEnsemble Learning PDE Solving
Harvard University
Curvature-Aware Optimization for High-Accuracy Physics-Informed Neural Networks
Apr 6, 2026Anas Jnini, Elham Kiyani, Khemraj Shukla +5Physics-Informed ML Neural PDE Solvers
Department of Information Engineering and Computer Science, University of Trento, Trento, Italy · Division of Applied Mathematics, Brown University, Providence, RI, USA · Departament de Física, Universitat d’Alacant, Spain +3
Gradient Scaling Effects in Adaptive Spectral PINNs for Stiff Nonlinear ODEs
May 6, 2026Isabela M. Yepes, Pavlos ProtopapasDynamical Systems Neural Network Optimization
Harvard University

Jul 27, 2026J/K move · Enter open · S save

Variational Boosting for Physics-Informed Neural Networks

Kaylee Vo, Pavlos Protopapas

Harvard University

Ensemble Learning PDE Solving Second-Order Optimization Physics-Informed ML

Open paper AskPDF arXiv

Physics-Informed Neural Networks (PINNs) solve differential equations by minimizing the residual of a nonlinear operator over a neural parameterization of the solution. However, monolithic PINNs often suffer from ill-conditioning, spectral bias, and optimization instability. We introduce a variational boosting framework in which solutions are constructed additively in function space. Each stage trains a weak learner whose converged correction satisfies a local orthogonality condition, equivalent to a projected functional gradient descent step onto the tangent space of the network's function manifold. Because each correction network is deliberately small, the restricted minimization admits full Newton or conjugate gradient updates, which are typically infeasible in large PINNs. The resulting method separates global nonlinear refinement into a sequence of well-conditioned subproblems while preserving the full variational structure of the operator. This framework provides a geometric interpretation of multi-stage PINNs as projected functional gradient descent and enables stable second-order optimization for nonlinear differential equations.

Abstract

Explore similar work

Variational Boosting for Physics-Informed Neural Networks

Curvature-Aware Optimization for High-Accuracy Physics-Informed Neural Networks

Gradient Scaling Effects in Adaptive Spectral PINNs for Stiff Nonlinear ODEs