Structure-preserving Randomized Neural Networks for Incompressible Magnetohydrodynamics Equations

Yunlong Li; Fei Wang; Lingxiao Li

physics.flu-dynMar 1, 2026

Structure-preserving Randomized Neural Networks for Incompressible Magnetohydrodynamics Equations

Authors: Yunlong Li, Fei Wang, Lingxiao Li

Organizations: School of Mathematics and Statistics, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, P.R. China. · School of Mathematics and Statistics & State Key Laboratory of Multiphase Flow in Power Engineering, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, China. · School of Mathematics and Statistics, Henan University, Kaifeng, Henan 475004, China. · Center for Applied Mathematics of Henan Province, Henan University, Zhengzhou, Henan 450046, China.

Abstract

The incompressible magnetohydrodynamic (MHD) equations are fundamental in many scientific and engineering applications. However, their strong nonlinearity and dual divergence-free constraints make them highly challenging for conventional numerical solvers. To overcome these difficulties, we propose a Structure-Preserving Randomized Neural Network (SP-RaNN) that automatically and exactly satisfies the divergence-free conditions. Unlike deep neural network (DNN) approaches that rely on expensive nonlinear and nonconvex optimization, SP-RaNN reformulates the training process into a linear least-squares system, thereby eliminating nonconvex optimization. The method linearizes the governing equations through Picard or Newton iterations, discretizes them at collocation points within the domain and on the boundaries using finite-difference schemes, and solves the resulting linear system via a linear least-squares procedure. By design, SP-RaNN preserves the intrinsic mathematical structure of the equations within a unified space-time framework, ensuring both stability and accuracy. Numerical experiments on the Navier-Stokes, Maxwell, and MHD equations demonstrate that SP-RaNN achieves higher accuracy, faster convergence, and exact enforcement of divergence-free constraints compared with both traditional numerical methods and DNN-based approaches. This structure-preserving framework provides an efficient and reliable tool for solving complex PDE systems while rigorously maintaining their underlying physical laws.

Figures & tables

Panel of Figure 1 — Figure 1: The structure of two-hidden-layers neural networks: the dotted line representing weights and bias which are tunable, and the solid line representing weights and bias which are randomly given and fixed.

Panel of Figure 2 — Figure 2: The structure of a two-hidden-layer neural network: the dotted lines represent tunable weights and biases, while the solid lines represent weights and biases that are randomly initialized and fixed.

Method	$m$ (DoFs)	$\\|\mathbf{u}-\mathbf{u}_{\rho}\\|_{1}$	$\\|\mathbf{u}-\mathbf{u}_{\rho}\\|_{0}$	$\\|p-p_{\rho}\\|_{0}$	$\\|\nabla\cdot\mathbf{u}\\|_{0}$
RaNN	$100\ (300)$	2.09E-03	1.50E-04	2.57E-03	1.83E-03
	$200\ (600)$	2.70E-05	1.65E-06	3.51E-05	2.13E-05
	$400\ (1200)$	3.91E-06	2.25E-07	4.81E-06	3.10E-06
	$800\ (2400)$	1.27E-06	7.47E-08	1.56E-06	9.63E-07
SP-RaNN	$100\ (200)$	2.09E-05	9.40E-07	2.35E-05	1.46E-12
	$200\ (400)$	8.02E-07	3.00E-08	1.00E-06	5.73E-13

Table 1: Numerical errors for RaNN and SP-RaNN with different

m

, where

N_{c}=5180,\ r_{1}=1,\ r_{2}=0.5,\ \gamma=500

, for Re

=1

in Example 5.1 ).

Method	$m$ (DoFs)	$\\|\mathbf{u}-\mathbf{u}_{\rho}\\|_{1}$	$\\|\mathbf{u}-\mathbf{u}_{\rho}\\|_{0}$	$\\|p-p_{\rho}\\|_{0}$	$\\|\nabla\cdot\mathbf{u}\\|_{0}$
RaNN	$100\ (300)$	1.61E-03	1.20E-04	1.03E-06	1.29E-05
	$200\ (600)$	4.16E-05	2.34E-06	2.28E-08	5.14E-07
	$400\ (1200)$	1.24E-05	6.43E-07	6.01E-09	1.94E-07
	$800\ (2400)$	5.86E-06	2.86E-07	3.04E-09	1.14E-07
SP-RaNN	$100\ (200)$	4.24E-06	1.80E-07	3.43E-09	1.73E-12
	$200\ (400)$	2.29E-07	9.36E-09	2.68E-10	5.84E-13

Table 2: Numerical errors for RaNN and SP-RaNN with different

m

, where

N_{c}=5180,\ r_{1}=1,\ r_{2}=0.5,\ \gamma=1

, for Re

=1000

in Example 5.1 ).

Re	DoFs	$\|\mathbf{u}-\mathbf{u}_{h}\|_{1}$	$\\|\mathbf{u}-\mathbf{u}_{h}\\|_{0}$	$\\|p-p_{h}\\|_{0}$
1	102414	1.98E-03	5.91E-06	4.0E-03
1000	102414	1.98E-03	5.91E-06	3.98E-03

Table 3: Numerical results from the mixed FEM with the modified Crouzeix–Raviart element ( [ 36 ] ) in Example 5.1 .

Method	$e$	$t=0$	$t=0.25$	$t=0.5$	$t=0.75$	$t=1$
SP-RaNN	$e_{u_{1}}$	7.39E-06	1.65E-06	1.38E-06	2.00E-06	1.31E-05
	$e_{u_{2}}$	7.01E-06	1.61E-06	1.37E-06	1.93E-06	1.21E-05
	$e_{u_{3}}$	6.08E-06	1.48E-06	1.31E-06	2.15E-06	1.17E-05
	$e_{p}$	4.83E-05	1.48E-05	2.21E-05	3.92E-05	2.25E-04
	$\\|\nabla\cdot\mathbf{u}\\|_{0}$	2.05E-11	2.08E-11	2.09E-11	2.11E-11	2.09E-11
NSFnets-VP	$e_{u_{1}}$	1.10E-04	4.08E-04	5.94E-04	9.08E-04	1.63E-03

Table 4: Comparison between NSFnets (Table 4 in [ 30 ] ) and SP-RaNN, where

m=3200,\ N_{c}=6144,\ r_{1}=r_{2}=0.2,\ \gamma=100

, for Re

=1

in Example 5.2 .

Re	$m$ (DoFs)	$e^{st}_{p}$	$e^{st}_{u}$	$e^{st}_{v}$	$e^{st}_{w}$
1	100 (200)	2.91E-01	1.17E-01	1.08E-01	1.33E-01
	200 (400)	1.61E-01	3.67E-02	3.84E-02	3.99E-02
	400 (800)	2.37E-02	5.50E-03	5.68E-03	5.41E-03
	800 (1600)	2.67E-03	4.48E-04	4.47E-04	4.40E-04
	1600 (3200)	1.49E-04	2.16E-05	2.18E-05	2.00E-05
1000	100 (200)	1.91E-01	1.10E-01	1.17E-01	1.18E-01

Table 5: Space-time relative

L^{2}

errors for SP-RaNN with different

m

, where

\ N_{c}=6144,\ r_{1}=r_{2}=0.2,\ \gamma=100

, for Re

=1

, and

\ N_{c}=6144,\ r_{1}=r_{2}=0.2,\ \gamma=10

, for

\text{Re}=1000

in Example 5.2 .

Panel of Figure 5 — Figure 5: The physical setting of lid-driven cavity flow and activation function in Example 5.3 .

Figure 6: Numerical solution

u

and its values on vertical centerline, numerical solution

v

and its values on horizontal centerline, where

r_{1}=6,\ r_{2}=3

for Re

=100

r_{1}=10,\ r_{2}=5

for Re

=400

and

r_{1}=15,\ r_{2}=10

for Re

=1000

in Example 5.3 .

Figure 7: Streamlines of the numerical solutions for different Reynolds numbers in Example 5.3 .

Methods	DoFs	$\\|\mathbf{E}-\mathbf{E}_{\rho}\\|_{0}$	$\\|\mathbf{H}-\mathbf{H}_{\rho}\\|_{0}$	$\\|\nabla\cdot\mathbf{H}_{\rho}\\|_{0}$
SP-RaNN	400	1.76E-02	6.54E-02	4.23E-14
	800	1.05E-02	2.23E-02	1.97E-13
	1600	2.15E-03	2.19E-03	1.00E-12
	3200	3.78E-04	2.35E-04	6.42E-12
	6400	5.27E-05	1.26E-05	7.46E-12
		$\\|\mathbf{E}-\mathbf{E}_{h}\\|_{0}$	$\\|\mathbf{H}-\mathbf{H}_{h}\\|_{0}$	$\\|\nabla\cdot\mathbf{H}_{h}\\|_{0}$

Table 6:

L^{2}

errors at

t=1

, where

N_{c}=6144,\,r_{1}=0.5,\ r_{2}=0.2,\ \gamma=100

in Example 5.4 .

Method	DoFs	$\\|\mathbf{A}-\mathbf{A}_{\rho}\\|_{0}$	Running time (s)
FEM	8368	1.12E-02	0.12
	62048	2.80E-03	1.02
	477376	6.98E-04	8.51
	3744128	1.74E-04	76.53
RaNN	300	2.76E-02	0.14
	600	1.90E-03	0.29

Table 7: Comparison between the preconditioned FEM and RaNN for Equation ( 5.17 ). Numerical errors and running times are reported for different DoFs, where

N_{c}=4725

r=1

, and

\gamma=100

Method	DoFs	$\\|\mathbf{u}-\mathbf{u}_{\rho}\\|_{1}$	$\\|p-p_{\rho}\\|_{0}$	$\\|\mathbf{B}-\mathbf{B}_{\rho}\\|_{\rm H(curl)}$	$\\|\nabla\cdot\mathbf{u}\\|_{0}$	$\\|\nabla\cdot\mathbf{B}\\|_{0}$	$N$	Running time (s)
FEM	4249	2.90E-03	1.84E-03	4.81E-02	-	-	8	1.72
	28749	7.08E-04	3.90E-04	2.38E-02	-	-	8	15.12
	210709	1.75E-04	9.07E-05	1.18E-02	-	-	8	135.46
	1611813	4.34E-05	2.21E-05	5.89E-03	-	-	7	981.22
SP-RaNN	300	1.39E-02	1.43E-02	6.55E-03	5.22E-16	2.10E-16	5	0.60
	600	1.56E-03	1.78E-03	5.58E-04	9.35E-16	4.12E-16	5	1.39

Table 8: Comparison between the preconditioned FEM and SP-RaNN for Equation ( 5.18 ). Numerical errors and running times are reported for different DoFs, where

N_{c}=1600

r_{1}=r_{2}=r_{3}=1

, and

\gamma=100

. Here,

N

denotes the number of Picard iterations.

$\gamma$	DoFs	$e_{\mathbf{u}}$	$e_{\nabla\mathbf{u}}$	$\\|\nabla\cdot\mathbf{u}\\|_{0}$	$e_{\mathbf{B}}$	$e_{\nabla\mathbf{B}}$	$\\|\nabla\cdot\mathbf{B}\\|_{0}$	$e_{p}$
$10$	$300$	1.03E-03	2.48E-03	6.71E-13	1.35E-02	1.02E-02	9.18E-13	5.67E-02
	$600$	5.68E-05	1.74E-04	3.39E-13	1.65E-03	1.20E-03	3.24E-13	4.93E-03
	$1200$	8.75E-06	2.95E-05	1.16E-13	3.65E-04	2.66E-04	1.28E-13	7.83E-04
	$2400$	2.84E-06	9.97E-06	9.30E-14	9.25E-05	6.85E-05	9.53E-14	2.68E-04
$10^{2}$	$300$	1.28E-04	3.87E-04	8.89E-13	3.29E-03	2.42E-03	7.46E-13	1.03E-02
	$600$	4.13E-06	1.49E-05	2.88E-13	6.16E-05	4.60E-05	2.92E-13	4.50E-04

Table 9: Numerical errors for various values of

\gamma

with the boundary condition

{\bf B}\cdot{\bf n}=0

, where

N_{c}=5180,\ r_{1}=r_{2}=1,\ r_{3}=0.5

in Example 5.6 .

Method	DoFs	$e_{\mathbf{u}}$	$e_{\nabla\mathbf{u}}$	$\\|\nabla\cdot\mathbf{u}\\|_{0}$	$e_{\mathbf{B}}$	$e_{\nabla\mathbf{B}}$	$\\|\nabla\cdot\mathbf{B}\\|_{0}$	$e_{p}$
SP-RaNN-P	$300$	9.92E-05	3.10E-04	8.36E-13	7.70E-05	2.45E-04	8.63E-13	8.81E-03
	$600$	3.53E-06	1.28E-05	2.62E-13	3.11E-06	9.07E-06	2.75E-13	3.34E-04
	$1200$	8.31E-07	3.30E-06	1.26E-13	1.31E-06	3.19E-06	1.82E-13	8.41E-05
	$2400$	3.52E-07	1.19E-06	8.44E-14	8.61E-07	1.19E-06	9.87E-14	2.72E-05
SP-RaNN-N	$300$	9.30E-05	2.75E-04	4.80E-13	8.97E-05	2.77E-04	1.12E-12	7.50E-03
	$600$	3.66E-06	1.26E-05	2.59E-13	4.11E-06	9.87E-06	2.75E-13	3.63E-04

Table 10: Comparison of WG-FEM (Table 4 in [ 54 ] ) and SP-RaNN, where

N_{c}=5180,\ r_{1}=r_{2}=1,\ r_{3}=0.5

and

\gamma=100

in Example 5.6 .

$\rm H_{a}$	DoFs	$e_{\mathbf{u}}$	$e_{\nabla\mathbf{u}}$	$\\|\nabla\cdot\mathbf{u}\\|_{0}$	$e_{\mathbf{B}}$	$e_{\nabla\mathbf{B}}$	$\\|\nabla\cdot\mathbf{B}\\|_{0}$	$e_{p}$
10	$2400$	5.72E-07	1.26E-06	2.08E-13	8.50E-07	1.14E-06	1.10E-13	1.10E-06
100	$2400$	5.08E-04	8.91E-04	4.49E-11	8.94E-07	1.37E-06	9.19E-14	8.39E-06
1000	$2400$	4.81E-03	7.89E-03	4.95E-11	1.05E-06	1.38E-06	1.01E-13	8.42E-05

Table 11: Errors for different Hartmann numbers using SP-RaNN-N, with

N_{c}=5180,\ r_{1}=r_{2}=1,\ r_{3}=0.5

and

\gamma=100

in Example 5.6 .

Method	DoFs	$e_{\mathbf{u}}$	$e_{\nabla\mathbf{u}}$	$\\|\nabla\cdot\mathbf{u}\\|_{0}$	$e_{\mathbf{B}}$	$e_{\nabla\mathbf{B}}$	$\\|\nabla\cdot\mathbf{B}\\|_{0}$	$e_{p}$
SP-RaNN-P	$300$	7.89E-05	4.60E-03	2.94E-12	1.38E-04	1.62E-04	2.88E-12	1.70E-04
	$600$	3.46E-07	6.40E-07	1.40E-12	6.83E-07	5.45E-07	1.75E-12	2.80E-05
	$1200$	6.03E-08	9.54E-08	4.56E-13	5.91E-07	2.97E-07	5.90E-13	4.27E-06
	$2400$	2.31E-08	3.78E-08	2.66E-13	6.79E-07	3.13E-07	3.41E-13	1.93E-06

Table 12: Errors for different

m

using SP-RaNN-P, where

N_{c}=5460,\ r_{1}=r_{2}=1,\ r_{3}=0.5

and

\gamma=100

in Example 5.6 .

Method	DoFs	$\\|\mathbf{u}-\mathbf{u}_{\rho}\\|_{0}$	$\|\mathbf{u}-\mathbf{u}_{\rho}\|_{1}$	$\\|p-p_{\rho}\\|_{0}$	$\\|\mathbf{B}-\mathbf{B}_{\rho}\\|_{0}$	$\\|\nabla\cdot\mathbf{u}_{\rho}\\|_{0}$	$\\|\nabla\cdot\mathbf{B}_{\rho}\\|_{0}$	time (s)
Picard	$300$	3.73E-04	3.61E-03	4.19E-03	5.07E-04	1.31E-14	2.16E-14	0.84
	$600$	1.55E-04	1.31E-03	1.79E-03	2.04E-04	2.51E-13	4.03E-13	1.48
	$1200$	4.10E-06	3.66E-05	8.35E-05	3.56E-06	7.84E-13	1.03E-12	3.39
	$2400$	8.93E-07	8.75E-06	1.96E-05	6.76E-07	1.54E-12	2.47E-12	10.72
	$4800$	5.31E-07	4.79E-06	1.12E-05	4.06E-07	1.13E-12	1.88E-12	51.24
Newton	$300$	3.88E-04	3.49E-03	4.46E-03	5.02E-04	1.55E-14	2.45E-14	1.05

Table 13: Errors of SP-RaNN at

t=0.2

and CPU time, where

N_{c}=3072,\ r_{1}=r_{2}=r_{3}=0.1

and

\gamma=100

in Example 5.7 .

	$(\tau,h)$	DoFs	$\\|\mathbf{u}-\mathbf{u}_{h}\\|_{0}$	$\|\mathbf{u}-\mathbf{u}_{h}\|_{1,h}$	$\\|p-p_{h}\\|_{0}$	$\\|\mathbf{B}-\mathbf{B}_{h}\\|_{0}$	$\\|\nabla\cdot\mathbf{u}_{h}\\|_{0}$	$\\|\nabla\cdot\mathbf{B}_{h}\\|_{0}$
k=1	$(\tau_{0},h_{0})$	604	3.87E-04	1.29E-02	1.77E-01	7.50E-02	1.71E-09	0
	$(\tau_{0},h_{0})/2$	4184	1.02E-04	6.23E-03	8.86E-02	3.69E-02	5.28E-10	0
	$(\tau_{0},h_{0})/4$	31024	2.83E-05	3.04E-03	4.43E-02	1.83E-02	7.27E-11	0
	$(\tau_{0},h_{0})/8$	238688	7.84E-06	1.50E-03	2.21E-02	9.08E-03	5.18E-10	0
k=2	$(\tau_{0},h_{0})$	1062	3.71E-04	1.28E-02	1.77E-01	5.14E-03	1.19E-09	4.43E-16
	$(\tau_{0},h_{0})/2$	7380	9.61E-05	6.19E-03	8.84E-02	1.26E-03	3.40E-10	1.33E-15

Table 14: Errors of the fully divergence-free FEM method ( [ 20 ] ) at

t=0.2

for varying time steps

\tau

and mesh sizes

h

. Here,

k

denotes the order of the piecewise polynomials used to approximate

\mathbf{A}

(

\mathbf{B}_{h}={\rm curl}\,\mathbf{A}_{h}

), while 1st-order elements and piecewise constant functions are used to approximate

\mathbf{u}

and

p

, respectively.

Method	DoFs	$\\|\mathbf{u}-\mathbf{u}_{\rho}\\|_{0}^{st}$	$\|\mathbf{u}-\mathbf{u}_{\rho}\|_{1}^{st}$	$\\|p-p_{\rho}\\|_{0}^{st}$	$\\|\mathbf{B}-\mathbf{B}_{\rho}\\|_{0}^{st}$	$\\|\nabla\cdot\mathbf{u}_{\rho}\\|_{0}^{st}$	$\\|\nabla\cdot\mathbf{B}_{\rho}\\|_{0}^{st}$	time (s)
Picard	$300$	7.51E-02	3.25E-01	2.81E-01	1.32E-01	4.23E-15	6.91E-15	0.83
	$600$	1.52E-02	9.88E-02	9.87E-02	2.94E-02	1.29E-14	1.97E-14	1.57
	$1200$	2.13E-03	1.48E-02	1.56E-02	4.17E-03	2.37E-14	3.89E-14	3.49
	$2400$	3.30E-04	2.50E-03	2.94E-03	7.42E-04	7.77E-14	1.28E-13	11.45
	$4800$	3.96E-05	3.50E-04	3.40E-04	1.39E-04	2.76E-13	4.66E-13	52.64
Newton	$300$	7.29E-02	3.39E-01	2.57E-01	1.24E-01	4.75E-15	6.58E-15	1.08

Table 15: Space-time errors over

(0,1)^{3}\times(0,5)

and CPU time for SP-RaNN, where

N_{c}=3072,\ r_{1}=r_{2}=r_{3}=0.3

and

\gamma=100

in Example 5.7 .

Structure-Preserving Neural Networks

See all 14 papers

Structure-Oriented Randomized Neural Networks for Poisson-Nernst-Planck and Poisson-Nernst-Planck-Navier-Stokes Systems
Jun 18, 2026Yunlong Li, Fei WangNeural PDE Solvers Structure-Preserving Neural Networks
School of Mathematics and Statistics, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, P.R. China. · School of Mathematics and Statistics & State Key Laboratory of Multiphase Flow in Power Engineering, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, China.
Random Neural Network Expressivity for Non-Linear Partial Differential Equations
May 24, 2026Muhammed Ali Mehmood, Lukas GononNeural PDE Solvers PDE Solving
Department of Mathematics Imperial College London UK · School of Computer Science University of St. Gallen Switzerland
Project and Generate: Divergence-Free Neural Operators for Incompressible Flows
Mar 25, 2026Xigui Li, Hongwei Zhang, Ruoxi Jiang +6Neural Operators Computational Fluid Dynamics
Artificial Intelligence Innovation and Incubation Institute, Fudan University, Shanghai, China · Shanghai Academy of Artificial Intelligence for Science, Shanghai, China

Jun 18, 2026J/K move · Enter open · S save

Structure-Oriented Randomized Neural Networks for Poisson-Nernst-Planck and Poisson-Nernst-Planck-Navier-Stokes Systems

Yunlong Li, Fei Wang

School of Mathematics and Statistics, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, P.R. China. · School of Mathematics and Statistics & State Key Laboratory of Multiphase Flow in Power Engineering, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, China.

Neural PDE Solvers Structure-Preserving Neural Networks Computational Fluid Dynamics

Open paper AskPDF arXiv

We develop a structure-oriented randomized neural network framework, termed SO-RaNN, for the Poisson-Nernst-Planck (PNP) system and the Poisson-Nernst-Planck-Navier-Stokes (PNP-NS) system. The decoupled linearized subproblems are solved iteratively by randomized neural networks in a space-time framework. For the concentration variables, a pointwise cut-off is used to enforce positivity at the value level, and discrete mass-scaling factors are computed at selected correction instants and interpolated in time, so as to ensure exact mass matching at those instants and to promote approximate mass preservation between them. To introduce an auxiliary discrete dissipation mechanism, we further employ an SAV-type post-processing correction, which yields monotonicity of the SAV auxiliary variable under the ideal SAV update. For the PNP-NS system, a structure-preserving randomized neural network (SP-RaNN) is used for the velocity field, so that the velocity approximation satisfies the incompressibility constraint pointwise by construction. On the theoretical side, we derive residual-based estimates for the raw, uncorrected RaNN solvers of the linearized subproblems, formulate a conditional local-in-time convergence result for the raw outer Picard iteration of the PNP system, and analyze the value-level positivity correction together with the mass-correction and SAV post-processing steps. For the PNP-NS system, we establish an approximation result for the SP-RaNN space and provide a conditional error statement for the corresponding linearized Oseen-type problem. Numerical experiments demonstrate approximation accuracy in the source-driven manufactured tests and illustrate the intended value-level positivity correction, selected-time mass matching, computed free-energy curves based on the final gauge-fixed potential, and divergence-free approximation in benchmark tests.

Abstract

Structure-Preserving Neural Networks

Explore similar work

Structure-Oriented Randomized Neural Networks for Poisson-Nernst-Planck and Poisson-Nernst-Planck-Navier-Stokes Systems

Random Neural Network Expressivity for Non-Linear Partial Differential Equations

Project and Generate: Divergence-Free Neural Operators for Incompressible Flows