Pub Date : 2026-06-01Epub Date: 2026-02-11DOI: 10.1016/j.jcp.2026.114756
Senwei Liang , Chunmei Wang , Xingjian Xu
Modeling stochastic differential equations (SDEs) is crucial for understanding complex dynamical systems in various scientific fields. Recent methods often employ neural network-based models, which typically represent SDEs through a combination of deterministic and stochastic terms. However, these models usually lack interpretability and have difficulty in generalizing beyond their training domain. This paper introduces the Finite Expression Method (FEX), a symbolic learning approach designed to derive interpretable mathematical representations of the deterministic component of SDEs. For the stochastic component, we integrate FEX with advanced generative modeling techniques to provide a comprehensive representation of SDEs. The numerical experiments on linear, nonlinear, and multidimensional SDEs demonstrate that FEX generalizes well beyond the training domain and delivers more accurate long-term predictions compared to neural network-based methods. The symbolic expressions identified by FEX not only improve prediction accuracy but also offer valuable scientific insights into the underlying dynamics of the systems.
{"title":"Identifying stochastic dynamics via finite expression methods","authors":"Senwei Liang , Chunmei Wang , Xingjian Xu","doi":"10.1016/j.jcp.2026.114756","DOIUrl":"10.1016/j.jcp.2026.114756","url":null,"abstract":"<div><div>Modeling stochastic differential equations (SDEs) is crucial for understanding complex dynamical systems in various scientific fields. Recent methods often employ neural network-based models, which typically represent SDEs through a combination of deterministic and stochastic terms. However, these models usually lack interpretability and have difficulty in generalizing beyond their training domain. This paper introduces the Finite Expression Method (FEX), a symbolic learning approach designed to derive interpretable mathematical representations of the deterministic component of SDEs. For the stochastic component, we integrate FEX with advanced generative modeling techniques to provide a comprehensive representation of SDEs. The numerical experiments on linear, nonlinear, and multidimensional SDEs demonstrate that FEX generalizes well beyond the training domain and delivers more accurate long-term predictions compared to neural network-based methods. The symbolic expressions identified by FEX not only improve prediction accuracy but also offer valuable scientific insights into the underlying dynamics of the systems.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114756"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146187260","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-05DOI: 10.1016/j.jcp.2026.114744
Yue Zhang , Peng Wang
Evaluation of failure probability is essential in system design and analysis. But it has become increasing prohibitive due to growing system complexity and soaring simulation costs. To address such challenge, we propose a novel active learning framework, the generalized Bayesian subset simulation (GBSS), to estimate failure probability. Based on the concept of generalized subset simulation, GBSS introduces a novel threshold-first criteria and decomposes the target failure events into an optimal sequence of intermediate events. Combined with a Bayesian surrogate model, high-quality scaled conditional samples may be generated efficiently with residual resampling in Markov chain Monte Carlo and in turn improves model accuracy. We evaluate our framework via a comprehensive set of academic examples and industrial cases of high dimensions. Comparing to five state-of-the-art methods, including Monte Carlo, subset simulation, generalized subset simulation, sequential directional importance sampling and multi-fidelity subset simulations, GBSS demonstrates substantial improvements in both numerical efficiency and stability.
{"title":"Generalized Bayesian subset simulation for the evaluation of failure probability","authors":"Yue Zhang , Peng Wang","doi":"10.1016/j.jcp.2026.114744","DOIUrl":"10.1016/j.jcp.2026.114744","url":null,"abstract":"<div><div>Evaluation of failure probability is essential in system design and analysis. But it has become increasing prohibitive due to growing system complexity and soaring simulation costs. To address such challenge, we propose a novel active learning framework, the generalized Bayesian subset simulation (GBSS), to estimate failure probability. Based on the concept of generalized subset simulation, GBSS introduces a novel threshold-first criteria and decomposes the target failure events into an optimal sequence of intermediate events. Combined with a Bayesian surrogate model, high-quality scaled conditional samples may be generated efficiently with residual resampling in Markov chain Monte Carlo and in turn improves model accuracy. We evaluate our framework via a comprehensive set of academic examples and industrial cases of high dimensions. Comparing to five state-of-the-art methods, including Monte Carlo, subset simulation, generalized subset simulation, sequential directional importance sampling and multi-fidelity subset simulations, GBSS demonstrates substantial improvements in both numerical efficiency and stability.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114744"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146187265","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-03DOI: 10.1016/j.jcp.2026.114747
Heming Bai, Xin Bian
The pressure Poisson equation, central to the fractional step method in incompressible flow simulations, incurs high computational costs due to the iterative solution of large-scale linear systems. To address this challenge, we introduce HyDEA (Hybrid Deep lEarning line-search directions and iterative methods for Accelerated solutions), a novel framework that synergizes deep learning with classical iterative solvers. It leverages the complementary strengths of a deep operator network (DeepONet) – capable of capturing large-scale features of the solution – and the conjugate gradient (CG) or a preconditioned conjugate gradient (PCG) (with Incomplete Cholesky, Jacobi, or Multigrid preconditioner) method, which efficiently resolves fine-scale errors. Specifically, within the framework of line-search methods, the DeepONet predicts search directions to accelerate convergence in solving sparse, symmetric-positive-definite linear systems, while the CG/PCG method ensures robustness through iterative refinement. The framework seamlessly extends to flows over solid structures via the decoupled immersed boundary projection method, preserving the linear system’s structure. Crucially, the DeepONet is trained on fabricated linear systems rather than flow-specific data, endowing it with inherent generalization across geometric complexities and Reynolds numbers without retraining. Benchmarks demonstrate HyDEA’s superior efficiency and accuracy over the CG/PCG methods for flows with no obstacles, single/multiple stationary obstacles, and one moving obstacle – using fixed network weights. Remarkably, HyDEA also exhibits super-resolution capability: although the DeepONet is trained on a 128 × 128 grid for Reynolds number , the hybrid solver delivers accurate solutions on a 512 × 512 grid for via interpolation, despite discretizations mismatch. In contrast, a purely data-driven DeepONet fails for complex flows, underscoring the necessity of hybridizing deep learning with iterative methods. HyDEA’s robustness, efficiency, and generalization across geometries, resolutions, and Reynolds numbers highlight its potential as a transformative solver for real-world fluid dynamics problems.
压力泊松方程是不可压缩流动模拟中分步法的核心，由于大规模线性系统的迭代求解，其计算成本很高。为了应对这一挑战，我们引入了HyDEA (Hybrid Deep lEarning line-search direction and iterative methods for Accelerated solutions)，这是一个将深度学习与经典迭代求解器相结合的新框架。它利用了深度算子网络（DeepONet）的互补优势——能够捕获解决方案的大规模特征——和共轭梯度（CG）或预条件共轭梯度（PCG）（具有不完全Cholesky， Jacobi或多网格预调节器）方法，有效地解决了精细尺度误差。具体而言，在线搜索方法的框架内，DeepONet预测搜索方向以加速求解稀疏、对称-正定线性系统的收敛，而CG/PCG方法通过迭代细化确保鲁棒性。框架通过解耦浸入式边界投影法无缝扩展到实体结构上的流动，保留了线性系统的结构。最重要的是，DeepONet是在合成的线性系统上进行训练的，而不是在特定流的数据上进行训练，这使其具有跨越几何复杂性和雷诺数的固有泛化能力，而无需重新训练。基准测试表明，对于无障碍物、单个/多个固定障碍物和一个移动障碍物（使用固定网络权重）的流动，HyDEA比CG/PCG方法具有更高的效率和准确性。值得注意的是，HyDEA还展示了超分辨率能力：尽管DeepONet是在雷诺数Re=1000的128 × 128网格上训练的，但混合求解器通过插值在Re=10,000的512 × 512网格上提供准确的解，尽管离散化不匹配。相比之下，纯数据驱动的DeepONet无法处理复杂的流，这强调了将深度学习与迭代方法相结合的必要性。HyDEA在几何、分辨率和雷诺数方面的稳健性、高效性和通用性突出了其作为现实世界流体动力学问题变革性求解器的潜力。
{"title":"Hybrid deep learning and iterative methods for accelerated solutions of viscous incompressible flow","authors":"Heming Bai, Xin Bian","doi":"10.1016/j.jcp.2026.114747","DOIUrl":"10.1016/j.jcp.2026.114747","url":null,"abstract":"<div><div>The pressure Poisson equation, central to the fractional step method in incompressible flow simulations, incurs high computational costs due to the iterative solution of large-scale linear systems. To address this challenge, we introduce HyDEA (Hybrid Deep lEarning line-search directions and iterative methods for Accelerated solutions), a novel framework that synergizes deep learning with classical iterative solvers. It leverages the complementary strengths of a deep operator network (DeepONet) – capable of capturing large-scale features of the solution – and the conjugate gradient (CG) or a preconditioned conjugate gradient (PCG) (with Incomplete Cholesky, Jacobi, or Multigrid preconditioner) method, which efficiently resolves fine-scale errors. Specifically, within the framework of line-search methods, the DeepONet predicts search directions to accelerate convergence in solving sparse, symmetric-positive-definite linear systems, while the CG/PCG method ensures robustness through iterative refinement. The framework seamlessly extends to flows over solid structures via the decoupled immersed boundary projection method, preserving the linear system’s structure. Crucially, the DeepONet is trained on <em>fabricated</em> linear systems rather than flow-specific data, endowing it with inherent generalization across geometric complexities and Reynolds numbers without retraining. Benchmarks demonstrate HyDEA’s superior efficiency and accuracy over the CG/PCG methods for flows with no obstacles, single/multiple stationary obstacles, and one moving obstacle – using <em>fixed network weights</em>. Remarkably, HyDEA also exhibits super-resolution capability: although the DeepONet is trained on a 128 × 128 grid for Reynolds number <span><math><mrow><mi>R</mi><mi>e</mi><mo>=</mo><mn>1000</mn></mrow></math></span>, the hybrid solver delivers accurate solutions on a 512 × 512 grid for <span><math><mrow><mi>R</mi><mi>e</mi><mo>=</mo><mn>10</mn><mo>,</mo><mn>000</mn></mrow></math></span> via interpolation, despite discretizations mismatch. In contrast, a purely data-driven DeepONet fails for complex flows, underscoring the necessity of hybridizing deep learning with iterative methods. HyDEA’s robustness, efficiency, and generalization across geometries, resolutions, and Reynolds numbers highlight its potential as a transformative solver for real-world fluid dynamics problems.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114747"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146147369","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-04DOI: 10.1016/j.jcp.2026.114753
Xin Liu , Yong Zhang
The convolution potential arises in a wide variety of application areas, and its efficient and accurate evaluation encounters three challenges: singularity, nonlocality and anisotropy. We introduce a fast algorithm based on a far-field smooth approximation of the kernel, where the bounded domain Fourier transform, one of the most essential difficulties, is well approximated by the whole space Fourier transform which usually admits explicit formula. The convolution is split into a regular and singular integral, and they are well resolved by trapezoidal rule and Fourier spectral method respectively. The scheme is simplified to a discrete convolution and is implemented efficiently with Fast Fourier Transform (FFT). Importantly, the tensor generation procedure is quite simple, highly efficient and independent of the anisotropy strength. It is easy to implement and achieves spectral accuracy with nearly optimal efficiency and minimum memory requirement. Rigorous error estimates and extensive numerical investigations, together with a comprehensive comparison, showcase its superiorities for different kernels.
{"title":"Fast convolution solver based on far-field smooth approximation","authors":"Xin Liu , Yong Zhang","doi":"10.1016/j.jcp.2026.114753","DOIUrl":"10.1016/j.jcp.2026.114753","url":null,"abstract":"<div><div>The convolution potential arises in a wide variety of application areas, and its efficient and accurate evaluation encounters three challenges: singularity, nonlocality and anisotropy. We introduce a fast algorithm based on a far-field smooth approximation of the kernel, where the bounded domain Fourier transform, one of the most essential difficulties, is well approximated by the whole space Fourier transform which usually admits explicit formula. The convolution is split into a regular and singular integral, and they are well resolved by trapezoidal rule and Fourier spectral method respectively. The scheme is simplified to a discrete convolution and is implemented efficiently with Fast Fourier Transform (FFT). Importantly, the tensor generation procedure is quite simple, highly efficient and independent of the anisotropy strength. It is easy to implement and achieves spectral accuracy with nearly optimal efficiency and minimum memory requirement. Rigorous error estimates and extensive numerical investigations, together with a comprehensive comparison, showcase its superiorities for different kernels.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114753"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146147536","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-04DOI: 10.1016/j.jcp.2026.114749
Daniel S. Finn , Joseph V. Pusztay , Matthew G. Knepley , Mark F. Adams
We present a novel structure-preserving framework for solving the Vlasov-Poisson-Landau system of equations using a particle in cell (PIC) discretization combined with discrete gradient time integrators. The Vlasov-Poisson-Landau system is an accurate model for studying hot plasma dynamics at a kinetic scale where small-angle Coulomb collisions dominate. Our scheme guarantees conservation of mass, momentum and energy as well as preservation of the monotonicity of entropy production in both the time-continuous and discrete systems. We employ the conservative integrator for both the Hamiltonian Vlasov-Poisson equations and the dissipative Landau equation using the PETSc library (www.mcs.anl.gov/petsc) to showcase structure-preserving properties.
{"title":"Structure preservation using discrete gradients in the Vlasov-Poisson-Landau system","authors":"Daniel S. Finn , Joseph V. Pusztay , Matthew G. Knepley , Mark F. Adams","doi":"10.1016/j.jcp.2026.114749","DOIUrl":"10.1016/j.jcp.2026.114749","url":null,"abstract":"<div><div>We present a novel structure-preserving framework for solving the Vlasov-Poisson-Landau system of equations using a particle in cell (PIC) discretization combined with discrete gradient time integrators. The Vlasov-Poisson-Landau system is an accurate model for studying hot plasma dynamics at a kinetic scale where small-angle Coulomb collisions dominate. Our scheme guarantees conservation of mass, momentum and energy as well as preservation of the monotonicity of entropy production in both the time-continuous and discrete systems. We employ the conservative integrator for both the Hamiltonian Vlasov-Poisson equations and the dissipative Landau equation using the <span>PETSc</span> library (<span><span>www.mcs.anl.gov/petsc</span><svg><path></path></svg></span>) to showcase structure-preserving properties.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114749"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146187261","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-10DOI: 10.1016/j.jcp.2026.114755
Huihui Cao , Yunqing Huang , Zhuoyun Li , Kailiang Wu
This paper proposes an innovative, compact, characteristic-decomposition-free, and non-intrusive convex-oscillation-suppressing (COS) framework, termed COS(DG), for high-order discontinuous Galerkin (DG) discretizations of hyperbolic systems on general meshes. The method performs an entropy-guided convex scaling per time step, blending each DG polynomial with its cell average via coefficients determined by locally scale- and evolution-invariant entropy-induced distances. We prove that COS(DG) preserves the optimal convergence rate of the underlying DG scheme for smooth solutions, is L2-nonexpansive, and inherits entropy stability whenever the base DG formulation is entropy stable. The analysis bounds entropy-induced discrepancies of DG polynomials through L2 estimates and constructs a finite, compact, convex covering inside the admissible state set. We also provide a theoretical justification that the problem-independent COS procedure guarantees uniform oscillation suppression across diverse problems, thanks to the combined local scale- and evolution-invariant structures.
The convex scaling strategy of COS(DG) offers a natural framework for integrating oscillation-suppressing and bound-preserving techniques. The COS approach is derivative-free, highly compact (using only immediate neighbors), and mode-independent (unified for modal and nodal DG), which facilitates parallelization and deployment on general meshes. A simple symmetry-preserving, scale-invariant entropy criterion detects interfaces near potential shocks and enables a local COS variant that further enhances resolution and reduces cost. Overall, COS(DG) offers nine salient advantages: (1) no characteristic decomposition; (2) provably optimal high-order convergence; (3) preservation of entropy stability (if present in the base DG method); (4) seamless integration with the Zhang–Shu bound-preserving limiter; (5) excellent compactness and easy parallelization using only neighboring data; (6) mode independence (modal/nodal unified); (7) simplicity via convex blending of solution and cell averages; (8) low overhead–applied only once per time step; and (9) local scale invariance (together with evolution invariance), enabling robust, problem-independent parameter choices across multi-scale discontinuities. Extensive benchmarks demonstrate the robustness, very high resolution, and superior performance of COS(DG), including comparisons with nine existing oscillation-control techniques, such as the classical TVB/WENO/MP-type limiters, as well as the recently developed OFDG and OEDG methods in [J. Lu, Y. Liu, and C.-W. Shu, SIAM J. Numer. Anal., 59 (2021)] and [M. Peng, Z. Sun, and K. Wu, Math. Comp., 94 (2025)].
{"title":"COS(DG): A convex oscillation-suppressing framework for high-order discontinuous Galerkin methods","authors":"Huihui Cao , Yunqing Huang , Zhuoyun Li , Kailiang Wu","doi":"10.1016/j.jcp.2026.114755","DOIUrl":"10.1016/j.jcp.2026.114755","url":null,"abstract":"<div><div>This paper proposes an innovative, compact, characteristic-decomposition-free, and non-intrusive convex-oscillation-suppressing (COS) framework, termed <span>COS(DG)</span>, for high-order discontinuous Galerkin (DG) discretizations of hyperbolic systems on general meshes. The method performs an entropy-guided convex scaling per time step, blending each DG polynomial with its cell average via coefficients determined by locally scale- and evolution-invariant entropy-induced distances. We prove that <span>COS(DG)</span> preserves the optimal convergence rate of the underlying DG scheme for smooth solutions, is <em>L</em><sup>2</sup>-nonexpansive, and inherits entropy stability whenever the base DG formulation is entropy stable. The analysis bounds entropy-induced discrepancies of DG polynomials through <em>L</em><sup>2</sup> estimates and constructs a finite, compact, convex covering inside the admissible state set. We also provide a theoretical justification that the problem-independent COS procedure guarantees uniform oscillation suppression across diverse problems, thanks to the combined local scale- and evolution-invariant structures.</div><div>The convex scaling strategy of <span>COS(DG)</span> offers a natural framework for integrating oscillation-suppressing and bound-preserving techniques. The COS approach is derivative-free, highly compact (using only immediate neighbors), and mode-independent (unified for modal and nodal DG), which facilitates parallelization and deployment on general meshes. A simple symmetry-preserving, scale-invariant entropy criterion detects interfaces near potential shocks and enables a local COS variant that further enhances resolution and reduces cost. Overall, <span>COS(DG)</span> offers <em>nine salient advantages</em>: (1) no characteristic decomposition; (2) provably optimal high-order convergence; (3) preservation of entropy stability (if present in the base DG method); (4) seamless integration with the Zhang–Shu bound-preserving limiter; (5) excellent compactness and easy parallelization using only neighboring data; (6) mode independence (modal/nodal unified); (7) simplicity via convex blending of solution and cell averages; (8) low overhead–applied only once per time step; and (9) <em>local</em> scale invariance (together with evolution invariance), enabling robust, problem-independent parameter choices across multi-scale discontinuities. Extensive benchmarks demonstrate the robustness, very high resolution, and superior performance of <span>COS(DG)</span>, including <em>comparisons with nine existing oscillation-control techniques</em>, such as the classical TVB/WENO/MP-type limiters, as well as the recently developed OFDG and OEDG methods in [J. Lu, Y. Liu, and C.-W. Shu, <em>SIAM J. Numer. Anal.</em>, 59 (2021)] and [M. Peng, Z. Sun, and K. Wu, <em>Math. Comp.</em>, 94 (2025)].</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114755"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146187264","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-11DOI: 10.1016/j.jcp.2026.114750
Jianghua Liu , Jin Zhao , Xuan Zhao , Xiaoli Li
Recently, a novel scalar auxiliary variable (SAV)-based optimization approach has been proposed in [J. Comput. Phys. 505 (2024): 112911]. The main idea is to view the optimization process of the loss function as an energy descent process for the gradient flow model. By leveraging the SAV method, which ensures linearity and monotonic modified energy decrease, an efficient optimization solver can be constructed. Compared with other gradient-based optimization algorithms such as gradient descent (GD), adaptive moment estimation (Adam), and limited memory quasi-Newton (L-BFGS) methods, the SAV based optimization approach allows the use of large step sizes, which may accelerate the optimization process. However, SAV-type schemes for the gradient flows decrease only the modified energy, not the original energy. At the same time, although many machine learning methods currently perform well in solving PDEs, most of them are highly sensitive to hyperparameters, such as the initial learning rate. To address the above limitations, we develop two effective and stable methods for solving PDEs by introducing the new Lagrange multiplier-based subspace approach, which provides a general framework for constructing linear schemes that dissipate the modified energy, as well as for developing schemes that dissipate the original energy through nonlinear algebraic equation. The constructed schemes first utilize the Lagrange multiplier method for preliminary optimization to identify subspaces that effectively capture the solution structure. Subsequently, the approximate solutions are obtained within these subspaces using the least squares approach. Numerical experiments demonstrate that our methods not only achieve high prediction accuracy but also exhibit robustness with respect to the initial learning rate.
{"title":"A machine learning method for solving PDEs with a robust initial learning rate based on Lagrange multiplier optimization","authors":"Jianghua Liu , Jin Zhao , Xuan Zhao , Xiaoli Li","doi":"10.1016/j.jcp.2026.114750","DOIUrl":"10.1016/j.jcp.2026.114750","url":null,"abstract":"<div><div>Recently, a novel scalar auxiliary variable (SAV)-based optimization approach has been proposed in <em>[J. Comput. Phys. 505 (2024): 112911]</em>. The main idea is to view the optimization process of the loss function as an energy descent process for the gradient flow model. By leveraging the SAV method, which ensures linearity and monotonic modified energy decrease, an efficient optimization solver can be constructed. Compared with other gradient-based optimization algorithms such as gradient descent (GD), adaptive moment estimation (Adam), and limited memory quasi-Newton (L-BFGS) methods, the SAV based optimization approach allows the use of large step sizes, which may accelerate the optimization process. However, SAV-type schemes for the gradient flows decrease only the modified energy, not the original energy. At the same time, although many machine learning methods currently perform well in solving PDEs, most of them are highly sensitive to hyperparameters, such as the initial learning rate. To address the above limitations, we develop two effective and stable methods for solving PDEs by introducing the new Lagrange multiplier-based subspace approach, which provides a general framework for constructing linear schemes that dissipate the modified energy, as well as for developing schemes that dissipate the original energy through nonlinear algebraic equation. The constructed schemes first utilize the Lagrange multiplier method for preliminary optimization to identify subspaces that effectively capture the solution structure. Subsequently, the approximate solutions are obtained within these subspaces using the least squares approach. Numerical experiments demonstrate that our methods not only achieve high prediction accuracy but also exhibit robustness with respect to the initial learning rate.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114750"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146187266","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-04DOI: 10.1016/j.jcp.2026.114757
Tariq D. Aslam, Eduardo Lozano
Here, we present an efficient level set method to track an arbitrary number of materials. The algorithm is optimal in the sense that it only needs to store a single unsigned distance-like function and a single integer indicator function, independent of the number of materials or distinct regions being tracked. For smooth velocity fields and smooth interface shape, arbitrarily high order solutions can be demonstrated. For interfaces that are or become kinked, the solution is limited to second-order convergence rates in the L1 norm and first-order in the L∞ norm.
{"title":"An efficient level set method for tracking many materials","authors":"Tariq D. Aslam, Eduardo Lozano","doi":"10.1016/j.jcp.2026.114757","DOIUrl":"10.1016/j.jcp.2026.114757","url":null,"abstract":"<div><div>Here, we present an efficient level set method to track an arbitrary number of materials. The algorithm is optimal in the sense that it only needs to store a single unsigned distance-like function and a single integer indicator function, independent of the number of materials or distinct regions being tracked. For smooth velocity fields and smooth interface shape, arbitrarily high order solutions can be demonstrated. For interfaces that are or become kinked, the solution is limited to second-order convergence rates in the <em>L</em><sub>1</sub> norm and first-order in the <em>L</em><sub>∞</sub> norm.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114757"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146187263","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-04DOI: 10.1016/j.jcp.2026.114746
Yiming Ren , Chuan Li , Guangqing Long , Shan Zhao
In this paper, we propose the first fourth-order Cartesian grid method built upon the augmented matched interface and boundary (AMIB) framework for solving two-dimensional elliptic interface problems with Lipschitz continuous interfaces involving sharp-edged corners. By using non-body-fitted meshes, traditional Cartesian grid methods are limited to be second-order accurate near sharp-edged corners, because there are insufficient grid nodes to support high-order approximations. To improve the efficiency and stability of interface methods, an augmented formulation is commonly used in the literature, which introduces auxiliary variables, such as Cartesian derivative jumps, so that the discrete Laplacian could be inverted by fast Poisson solvers, including fast Fourier transform and multigrid. In this work, a new correction to the fourth-order central difference approximation of the Laplacian is proposed, by treating the fictitious jumps, which measure the differences between fictitious values (FVs) and function values at irregular grid points near the interface, as the auxiliary variables in the AMIB method. Besides being simpler and involving less approximations, this new correction allows us to generate FVs in an entirely different manner, i.e., one FV could potentially depend on all other unknown FVs. This makes a lot of FVs available near sharp-edged corners to discretize the jump conditions, so that fourth-order approximations could be accomplished for problems with piecewise C6 continuous solutions. These unknown fictitious jumps will be solved together with function values in an enlarged linear system by using a Schur complement procedure combined with fast Poisson solvers. Numerical experiments demonstrate that the AMIB-FV scheme achieves fourth-order accuracy for both solutions and their gradients, while maintaining the overall efficiency as O(n2log n) or O(n2) for an n × n uniform grid, in solving complex interface problems with mixed boundary conditions. Moreover, the condition numbers of the AMIB-FV scheme are significantly smaller than those of the existing AMIB methods.
{"title":"A fourth-order Cartesian grid method for elliptic problems with sharp-edged interfaces","authors":"Yiming Ren , Chuan Li , Guangqing Long , Shan Zhao","doi":"10.1016/j.jcp.2026.114746","DOIUrl":"10.1016/j.jcp.2026.114746","url":null,"abstract":"<div><div>In this paper, we propose the first fourth-order Cartesian grid method built upon the augmented matched interface and boundary (AMIB) framework for solving two-dimensional elliptic interface problems with Lipschitz continuous interfaces involving sharp-edged corners. By using non-body-fitted meshes, traditional Cartesian grid methods are limited to be second-order accurate near sharp-edged corners, because there are insufficient grid nodes to support high-order approximations. To improve the efficiency and stability of interface methods, an augmented formulation is commonly used in the literature, which introduces auxiliary variables, such as Cartesian derivative jumps, so that the discrete Laplacian could be inverted by fast Poisson solvers, including fast Fourier transform and multigrid. In this work, a new correction to the fourth-order central difference approximation of the Laplacian is proposed, by treating the fictitious jumps, which measure the differences between fictitious values (FVs) and function values at irregular grid points near the interface, as the auxiliary variables in the AMIB method. Besides being simpler and involving less approximations, this new correction allows us to generate FVs in an entirely different manner, i.e., one FV could potentially depend on all other unknown FVs. This makes a lot of FVs available near sharp-edged corners to discretize the jump conditions, so that fourth-order approximations could be accomplished for problems with piecewise <em>C</em><sup>6</sup> continuous solutions. These unknown fictitious jumps will be solved together with function values in an enlarged linear system by using a Schur complement procedure combined with fast Poisson solvers. Numerical experiments demonstrate that the AMIB-FV scheme achieves fourth-order accuracy for both solutions and their gradients, while maintaining the overall efficiency as <em>O</em>(<em>n</em><sup>2</sup>log <em>n</em>) or <em>O</em>(<em>n</em><sup>2</sup>) for an <em>n</em> × <em>n</em> uniform grid, in solving complex interface problems with mixed boundary conditions. Moreover, the condition numbers of the AMIB-FV scheme are significantly smaller than those of the existing AMIB methods.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114746"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146147535","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-06-01Epub Date: 2026-02-06DOI: 10.1016/j.jcp.2026.114759
Zhengshou Lai , Shuai Huang , Yong Kong , Shiwei Zhao , Jidong Zhao , Linchong Huang
This study presents a hybrid resolved and unresolved computational fluid dynamics-discrete element method (CFD-DEM) coupling framework for modeling fluid-particle systems involving irregularly shaped and polydisperse particles. The framework integrates high-fidelity signed distance field (SDF) representations for coarse particles with efficient coarse-fine contact algorithms and unresolved or semi-resolved fluid-particle interaction schemes for fine particles. Key features such as accurate representation of irregularly shaped particles, consideration of varying particle sizes, and evaluation of fluid-particle interactions in complex granular settings are highlighted. The hybrid CFD-DEM solver is implemented and verified through a series of standard benchmark tests. A detailed filtration example is provided to demonstrate its capability in modeling multiscale fluid-particle interactions and in analyzing the influence of particle size and shape on filtration behavior. Additional case studies, including channelized sorting, jet-induced destabilization, and dam-break flows, further illustrate the flexibility and effectiveness of the proposed approach in practical applications.
{"title":"Hybrid resolved-unresolved CFD-DEM framework for multiscale fluid-particle systems with irregular-shaped and polydisperse particles","authors":"Zhengshou Lai , Shuai Huang , Yong Kong , Shiwei Zhao , Jidong Zhao , Linchong Huang","doi":"10.1016/j.jcp.2026.114759","DOIUrl":"10.1016/j.jcp.2026.114759","url":null,"abstract":"<div><div>This study presents a hybrid resolved and unresolved computational fluid dynamics-discrete element method (CFD-DEM) coupling framework for modeling fluid-particle systems involving irregularly shaped and polydisperse particles. The framework integrates high-fidelity signed distance field (SDF) representations for coarse particles with efficient coarse-fine contact algorithms and unresolved or semi-resolved fluid-particle interaction schemes for fine particles. Key features such as accurate representation of irregularly shaped particles, consideration of varying particle sizes, and evaluation of fluid-particle interactions in complex granular settings are highlighted. The hybrid CFD-DEM solver is implemented and verified through a series of standard benchmark tests. A detailed filtration example is provided to demonstrate its capability in modeling multiscale fluid-particle interactions and in analyzing the influence of particle size and shape on filtration behavior. Additional case studies, including channelized sorting, jet-induced destabilization, and dam-break flows, further illustrate the flexibility and effectiveness of the proposed approach in practical applications.</div></div>","PeriodicalId":352,"journal":{"name":"Journal of Computational Physics","volume":"554 ","pages":"Article 114759"},"PeriodicalIF":3.8,"publicationDate":"2026-06-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"146187262","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
扫码关注我们
求助内容：
应助结果提醒方式：
京公网安备 11010802042870号