Solving elliptic partial differential equation based on deep reinforcement learning

Elliptic partial differential equations (PDEs) constitute a fundamental class of governing equations in physics and engineering, underpinning the mathematical descriptions of electrostatic potential distribution, steady-state heat conduction, and linear elasticity. The numerical solution of these equations on fine discretization grids is computationally intensive, as classical iterative solvers require numerous sweeps to achieve convergence. This computational bottleneck critically hinders real-time applications such as online process control, interactive design optimization, and inverse problem reconstruction, where rapid response to changing boundary conditions or source terms is essential. The traditional numerical method cannot meet the real-time application requirements because of its long calculation time when solving on fine grids. In this study, the solution framework of deep Q-network is constructed, and the elliptic PDE solution is modeled as Markov decision processes. The state space includes solution field, residual field and boundary mark on discrete grid. The action space is 11 discrete corrections (− 0.10 to 0.10 with an interval of 0.02). The reward function is weighted sum of PDE residual term (weight 0.6), boundary error term (weight 0.3) and smooth regularization term (weight 0.1). The feature extraction network uses three-layer convolution and two residual blocks to output 512-dimensional vectors. The proposed DRL solver achieves relative L2 errors of 0.28%, 0.39%, and 0.62% on Poisson, Helmholtz, and anisotropic diffusion test cases, respectively, outperforming the finite element method (0.32%, 0.45%, 0.78%) and the multigrid method (0.31%, 0.44%, 0.76%). In the inference phase, the DRL solver completes a single forward pass on a 256 × 2562 grid in 0.15 s, achieving a 9898-fold speedup over the multigrid method’s 14.8 s. The training process converges stably within 20,000 steps, with the residual norm decreasing from an initial value of 12.4 to a final value of 0.16. This method provides a feasible way for solving and controlling real-time physical fields.