CodePDE: An Inference Framework for LLM-driven PDE Solver Generation
Authors: Shanda Li
Tanya Marwah
Junhong Shen
Weiwei Sun
Andrej Risteski
Yiming Yang
Ameet Talwalkar
Venue: Accepted by TMLR ()

Abstract: Partial differential equations (PDEs) are fundamental to modeling physical systems, yet solving them remains a complex challenge. Traditional numerical solvers rely on expert knowledge to implement and are computationally expensive, while neural-network-based solvers require large training datasets and often lack interpretability. In this work, we frame PDE solving as a code generation task and introduce CodePDE, the first inference framework for generating PDE solvers using large language models (LLMs). With CodePDE, we present a thorough evaluation on critical capacities of LLM for PDE solving: reasoning, debugging, self-refinement, and test-time scaling. CodePDE shows that, with advanced inference-time algorithms and scaling strategies, LLMs can achieve strong performance across a range of representative PDE problems. We also identify novel insights into LLM-driven solver generation, such as trade-offs between solver reliability and sophistication, design principles for LLM-powered PDE solving agents, and failure modes for LLM on hard tasks. These insights offer guidance for building more capable and reliable LLM-based scientific engines.

Published in Transactions on Machine Learning Research (02/2026)
CodePDE: An Inference Framework for LLM-driven PDE Solver Generation
Shanda Li Cornegie Mellon University
shandal@es.co www.:l/wwwwc.cm..c.mmmmmmm. edu
Tanya Marwah Flatiron Institute; Polymathic AT
Smarwaalfffairoonnnssttttto.y
Junhong Shen
junhongs@@s.ccmu..l
Cornegie Mellon University
Weiwei Sun Carnegic Mellon University
weiweis@cs. ctttteiisssscmuu.com cmnu, edu
Andrej Risteski Carnegie Mellon University
aristesk@cs..
cesstes@@cs.cmm.oom cma, edu
Yiming Yang
yiming@os.
comiing@esccu.cooo cmu. edu
Carnegie Mellon University
Ameet Talwalkar
atalwalk@ccs.m cmn, edu
Carnegie Mellon University; Datadog
Reviewed OpenReview: https: openreview. net/ formm? id- eG3Qy50uz6
Abstract
Partial differential equations (PDEs) are fundamental to modeling physical systems, yet solving them remains a complex challenge. Traditional numerical solvers rely on expert knowledge to implement and are computationally expensive, while neural-network- based solvers require large training datasets and often lack interpretability. In this work, we frame PDE solving as code generation task and introduce CodePDE, the first inference framework for generating PDE solvers using large language models (LLMs). With CodePDE, present thorough evaluation critical capacities of LLM for PDE solving; reasoning, debugging, self-refinement, and test-time scaling. CodePDE shows that, with advanced inference-time algorithms and scaling strategies, LLMs can achieve strong performance across a range of representative PDE problems. We also identify novel insights into LLM-driven solver generation, such trade-offs between solver reliability and sophistication, design principles for LLM-powered PDE solving agents, and failure modes for LLM on hard tasks. These insights offer guidance for building more capable and reliable L.LM-based scientific engines.
Github
HuggingFace
1 Introduction
Partial differential equations (PDEs) are foundational for modcling complex phenomena science and engineering, yet developing effective numerical solvers remains a substantial challenge. Traditional numerical methods demand significant computational resources and deep domain expertise for implementation. Crafting

---

Published in Transactions on Machine Learning Research (02/2026)
robust and efficient solvers typically involves meticulous manual tuning, extensive debugging, and specialized numerical analysis knowledge. The success of deep learning has motivated the development of neural PDE solvers (e.g,, Raissi et al., 2019; Lu et al., 2021; Li et al., 2022b) and multiphysics foundation models (e.g., McCabe et al., 2024; Hao et al., 2024; Shen et al., 2024b; Herde et al., 2024) to automate PDE solving. However, training these models requires large amounts of data, and their outputs lack transparency, making it difficult to understand how and why a particular solution is produced.
Meanwhile, large language models (LLMs) have been increasingly applied to complex mathematical and scientific challenges (e.g,, Romera-Paredee et al., 2024; Shen et al., 2024a; 2025; Chan et al., 2025; Sun et al,, 2025; Zhou et al., 2025; Li et al., 2026). Central to these advancements is the observation that code acts as a versatile intermediary between natural language and structured symbolic scientific computations. Inspired by this insight, a compelling alternative paradigm for PDE solving emerges: automated generation of executable solver code directly from high-level natural language descriptions of PDEs. This emerging paradigm raises a critical and open rescarch question:
Can LLMs generate effective and efficient solver code for PDEs?
Despite the conceptual simplicity, naively applying LLMs to PDE solving may underutilize their capabilities, neglecting recent progress in inference-time algorithms (Welleck et al., 2024) and scaling strategies (Snell et al., 2025). We introduce CodePDE, inference framework for LLM-driven automated PDE solver generation. Given a problem description, CodePDE instructs LLMs to produce diverse solver implementations (e.g., finite difference, spectral methods, or other novel approaches) CodePDE further incorporates mechanisms for code repair, iterative refinement, test-time scaling, and rigorous evaluation of accuracy, efficiency, and convergence, as depicted in Figure 1. Designed with forward compatibility and modularity, CodePDE supports both local deployments and APl-based interfaces, and readily integrates with a variety of inference strategies and agentic workflows.
Deploying CodePDE with 16 strong LLMs, we conduct systematic evaluations across four core capabilities for PDE solver generation:
Generating numerically sound solvers through thain-of-thought. LLMs can explore the combination of wide range of methods, such as finite difference, finite volume, and spectral methods for spatial domains, well as Explicit Euler, Runge-Kutta, and IMEX for time integration, yielding diverse solver strategies.
Reasoning:
Identifying and autonomously correcting code errors based on runtime feedback. The average rate of bug-free solver generation increases from 41% to 84% after applying iterative self-debugging, demonstrating the effectiveness of automated error repair.
Debugging:
Improving solver accuracy through feedback-based improvement. We find that enabling self-refinement consistently enhances solution quality across all tested PDEs, underscoring the value of feedback-driven optimization.
Refinement:
Test-Time Scaling:
Boosting solution quality by scaling up inference compule. Using best-of-n sampling, we observe that solution accuracy improves with increased inference budget (Figure 3), highlighting practical scaling law for LLM-generated solvers.
These capabilities constitute essential building blocks for LLM-powered scientific computing. Our experiments demonstrate that LLMs, when paired with inference-time techniques such as automated debugging (Chen et al., 2024), self-refinement (Madaan et al., 2023), and test-time scaling (Snell et al., 2025) -cchivee performance comparable to tailored numerical solvers and specialized PDE solving softwares. Further ablation study shows that these techniques are critical ingredients for building highly accurate PDE-solving agentic workflows.
CodePDE also provides a platform for systematically evaluating LLMs on solver generation, uncovering several key insights. For example, our analysis identifies a trado-off between solver reliability and sophistication, where some models conservatively produce simple, robust solvers while others explore more diverse complex mothods. In addition, we find that code generation and refinement may be distinct skills for LLMs: while
2

---

Published in Transactions on Machine Learning Research (02/2026)
du ar
a² 0
dx²
def solver(...
Task S Specification
Code Generation
Debugging
Evaluation
Solver Refinement
An overview of CodePDE framework. Critical LLM capabilities unlocked in the steps: code generation leverages chain-of-thone thoughh reasoning, code repair enables autonomous debugging, evaluation implements test-time scaling through best-of-n sampling, and solver refinement optimizes performance through feedback-driven improvement.
Figure 1:
advanced reasoning models like o3 and DeepSeek-RI excel at generating high-quality solvers from scratch, they do not consistently outperform standard models like GPT-4o and DeepSeek-V V3 in the refinement. stage.
In summary, our contributions threefold: (1) We frame PDE solving code generation task and develop CodePDE, the first inference framework to combine PDE domain knowledge with the celf-improvement and test-time scaling capabilities of LLMs. (2) We show that CodePDE unlocks the capabilities of LLM to generate high-quality PDE solvers. (3) We present the first comprehensive study of LLM--enerated PDE solvers, analyzing their accuracy, efficiency, numerical scheme diversity, and failure modes-ofe offering new insights into the potential and limitations of L.LSs in scientific computing.
2 Related Works: PDE Solving and the Emerging Role of LLMs
Numerical Methods. Traditional numerical solvers approximate PDE solutions via domain discretization. Finite difference method (FDM) (LeVeque, 2007) uses grid-based differences to estimate derivatives; finite element method (FEM) (Zienkiewizz & Taylor, 2013) approximates solutions over mesh elements; and spectral methods (Canuto et al., 2007) represent solutions as sums of global basis functions. These solvers offer convergence guarantees and error bounds, but they require expert knowledge to implement and are computationally intensive, particularly for high-dimensional problems. In recent years, a few softwares have emerged to streamline the use of these numerical methods (Zwicker, 2020; Burns et al., 2020; Scroggs et al., 2022).
Neural PDE Solvers. Neural network methods such as PINNs (Raissi et al., 2019; Wang et al., 2022; He et al., 2023) and neural operators (Lu et al., 2021; Li et al., 2022b; 2025) provide data-driven alternatives. While bypassing some traditional bottlenccks, they are often tailored to specific equations. Recent works explore diverse architectures on PDE problems, such as Transformers (Cao, 2021; Li al., 2022a; Shen et al., 2023), GNNs (Li et al., 2020; Pfaff et al., 2021; Brandstetter al., 2022), state-space models (Zheng et al., 2024; Buitrago et al., 2025), and other specialized networks (Long et al., 2018; Marwah al., 2023; Shen et al., 2022; Liang et al., 2025). Pretrained multiphysics foundation models (Shen et al., 2024b; Subramanian et al., 2023; McCabe et al., 2024) are also developed to improve generalization ability. However, these models require expensive offline training. Their solution process is not interpretable, and they do not take advantage of well-cstablished mumerical techniques.
LLMs for Scientific Problem Solving. LLMs have shown strong performance generating exccutable Guo et al., 2023; Shen et al., 2024c), physics (e.g., Arlt et al., 2024), mathematics (e.g., Hong et al, 2024; Wang et al., 2023; Feng et al., 2026), and computational biology (e.g., Tang et al, 2024; Haase et al., 2024; Cheng et al, 2024). Agentic workflows further enhance capabilities through planning and itcrative refinement (c.g., Romera-Paredes et al., 2024; Ma et al., 2024). For example, FunSearch KRomera-Paredee et al., 2024) uses sell-prompting with an evolutionary algorithm for mathematical discovery. AIDE (Jiang et al., 2025) performs tree search over solution strategies for machine learning engincering. Another very recent work, PDE-Contron
code for tasks in chemistry (e.g., Bran al., 2023;
(Soroco et al., 2025), uses LLMs to control systems governed by PDEs. In contrast, our work is the first study PDE solver code generation with LLM inference algorithms and scaling strategies and systematically evaluate L.M capabilities in this domain.
3

---

Published in Transactions on Machine Learning Research (02/2026)
3 Method: CodePDE
In this section, we present CodePDE, our inference framework for LLM-driven PDE solver generation. We begin by introducing the problem background, followed by an overview of our five-step framework design as visualized in Figure
Preliminaries.
We focus on PDEs of the general form:
Lu( ((() Bu(() = (()
EENCH
(1)
Ø E 0n
where is the peatiilleempora//ppaaaia variable for
PDE, Q denotes the spatial-tomporal/spatiiii domain, and L and B are partial differential operators defining the governing equation and boundary condition, respectively. PDEs in this form encompass Navier-Stokes, Reaction-Diffusio, and Burgers' Equations, as well as numerous others that describe phenomena in all sorts of domains. These equations constitute the majority PDE benchmarks machine learning research (TTu et al., 2022; Roberts et al., 2021; Takamoto et al., 2022).
A Analytical solutions can only be derived for highly restricted classes of PDEs, highlighting the need for computing numerical solutions. Numerical solvers discretize the domain into grids and output solution values at each grid point over specified time steps.
Framework Design. CodePDE implemented in Python, leveraging its widespread use in the ML community and the availability of modern numerical solver libraries in Python (Takamoto et al., 2022; Kochkov et al., 2021; Zwicker, 2020). Our full codebase released at https://githhn: com/LLthhiuDD//Codeeee Our evaluation data is released https: /tt.ggggggaacceeccm/////// httge//ttttt////gggggggggaa hwwpt//ffae...oaataaaeeeeeeeeeee /tt/aaaattt/LLAA220/// cttps//DDDDD02///oeeeeeeeeeeeeee wotppe/ddttatteeeeeeeeeeeeeeeeee CodePDE has also been intergrated into Terminal Bench (Merrill et al., 2026) via Harbor Adapters. Below, detail each step in CodePDE.
Step 1: Task Specification. Given a PDE problem of the form specified in Equation 1, we format it into natural language descriptions so that it can be processed by LLMs. Specifically, we include the governing equations, domain specifications, boundary conditions, and initial conditions. For example, a Burgers Equation task specified as follows:
Your task is to solve partial differential equation (PDE) using Python. The PDE the Burgers equation with the periodic boundary condition, given by
Oeu(x, 1) 0. (="(a)) uI,, f), x (0, 1), LE (0, u(x,0) uz(),,
(0,1)
where a constant representing the viscosity.
Step 2: Code Generation. After specifying the input, we prompt models generate complete solver implementations alongside any necessary auxiliary functions. To facilitate evaluation and enhance code readability, we instruct language models to implement solutions based on predefined function signature, which receives initial conditions and time grids and returns the predicted solution trajectory:
daf solver (uo, batch, c_coordnnate, nu): Solves the Burgers' equation for all times in t_coordinate.
Arga:
w0_batch (ap adarray); Initial condition [batch_eize, X]. t_coordinate (ap adarray): Time coordinates of shape [T+1]. (fleat): Viscosity coafficient.
Returns:
solutions ndarray): Shape (batch_size, T+1, N].
Implement solver for Burgers' equation.
# Bint: Consider PyTorch/JAX with GPU acceleration for officiency. solutions
4

---

Published in Transactions on Machine Learning Research (02/2026)
We use chain-of-thought prompting (Wei et al., 2022) to instruct the model to navigate complex numerical algorithms. The full prompt used can be found in Appendix B.1.
Step 3: Debugging. After the LLM returns a solver, we execute it to check its validity. If the solver crashes, we feed the error traces, the original specification, and the failed code to the LLM for it to revise its solution. This process instructs the model to diagnose root causes and produce corrected implementations without human intervention. Debugging proceeds iteratively until the issue is resolved or predefined iteration limit is reached. The debugging prompt is detailed in Appendix B.2. We quantify effectiveness by recording the debug success rate for each model on each problem.
Step 4: Evaluation. For executable solvers, we evaluate their performance by calling the solver, obtaining the predicted solution, and comparing it against reference solutions. We investigate three metrics. First, we compute the error with respect to the ground truth solution. We follow prior works (Takamoto et al., 2022; Shen et al., 2024b) and use the scale-independent normalizod root mcan squared error (nRMSE), defined as:
S [1(()) t) a (x,t)ll (x,
RRMSE
(2)
S
where S denotes the number of examples in PDE family. Second, we measure the quality of the solver using a convergence test (LeVeque, 2007), which assesses how the solution error decreases as the grid is refined. Mathematically, a solver is considered convergent the difference between solutions at successive resolutions decreases with finer discrotization. That is, for a grid spacing h, we test whether [|u,, - ua/22l2 + 0 as h O. This test ensures that the numerical solution remains stable and consistent as resolution increases, even in the absence of an exact solution. The test result also serve as solver selection criteria. Finally, we record code execution time measure of computational efficiency.
Step 5: Solver Refinement. susplp the nRMSE obtained in step 4 along with the solver implementation to the LLM for further refining the solution. Models analyze execution results to identify bottlenecks and numerical instabilities, then generate improved implementations accordingly. The refinement prompt can be found in Appendix B.2.
We highlight that our framework is general and forward-compatibie The entire CodePDE framework accepts both local hosting of LLMs as well as standardized APIs. When new specialized models agents are developed later, they can be seamlessly incorporated into our framework. The framework's modular design facilitates systematic analysis of LLM performance across atomic operations (e.g., implementation, debugging, and refinement), decomposing performance along these dimensions to identify each model's strengths and limitations.
4 Evaluation Setup
This section briefly introduces our experiment setup, with more details presented in Appendix B.3.
Datasets. We focus on 5 representative PDE families that are commonly employed as testbeds for ML methods: Advection, Burgers, Reaction-Diffusion Compressible Navier-Stokes (CNS), and Darcy Flow. These span a range of spatial dimensions, numerical stiffness, and boundary conditions (Saad et al., 2023). The datasets are drawn from PDEBench (Takamoto et al., 2022) and FNO paper (Li et al., 2022b), both with MIT license. We sample 100 instances for cach family.
Language Models and Baseline Methods. We benchmark 16 LLMs, spanning both proprictary models and open-weights models. The full model list can be found in Appendix Table 3. We compare LLM-generated solvers against standard numerical solvers, including implementations based on established PDE solving softwares (Zwicker, 2020; Burns et al., 2020) as well as manually crafted solvers using accurate schemes. For example, we use the analytical solution for Advection, a spectral method for Burgers, and Strang splitting for the Reaction-Diffusion Equation. The code can be found in the supplementary materials. We also include a variety of standard neural network solvers and foundation models (Ronneberger et al., 2015; Li al., 2022b; Raissi et al., 2019; Ye et, al., 2024; Shen et al., 2023; 2024b). Additionally, we evaluate 2 agentic
5

---

Published in Transactions on Machine Learning Research (02/2026)
workflows which feature search and refinement in the code space: FunSearch (Romera-Paredee et al,, 2024) and ATDE (Jiang et al., 2025).
Evaluation Metrics. We use normalized root mean squared error (nRMSE) as the primary measure of solution quality. We also report (1) debug success rate: the proportion of generations successfully fixed after error-driven feodback; (2) convergence rate: empirical convergence rate with incrcased spatial resolution; and (3) exccution time: runtime of generated solvers.
5 Results and Analysis
In this section, we present our main evaluation results and discuss the key insights from few axes.
5.1 How Well Do LLMs Generate PDE Solvers?
We begin by examining the comparative performance of CodePDE solvers against baselines. Table 1 presents the normalized RMSE (nRMSE) across five representative PDE families. Only subset of LLMs are included due to space constraint. Full results are provided in Table 3 in Appendix A.
The results demonstrate that LLMs under the CodePDE framework can generate high-quality PDE solvers that are competitive with the reference solvers. Specifically, we find that CodePDE with refinement outperforms the hand-crafted reference solvers on 4 out of evaluated tasks. Specifically, on Burgers Equation, Gemini 2.0 Flash achieves an nRMSE of 1.06 10-1 via refinement, compared to the reference's 's 3.55 x 10-4, representing 70% improvement in accuracy. Neural PDE solvers, including both task-specifie models (e.g., FNO) and foundation models (e.g., UPS), generally perform worse than the best CodePDE implementations, although they lypically require training from scratch or additional finetuning.
However, we also note limitations with LLM-generated solvers. The Reaction- Reaction-Diffusion Equation remains challenging for all LLM-generated solvers, which yield higher nRMSE values than the reference. This difficulty arises because the task requires nuanced insights into the problem structure, which we elaborate on in Section 5.7.
Takeaways. LLM-generated solvers, particularly those enhanced with debugging and refinement mechanisms, hold strong promise. They demonstrate strong performance on many representative PDE problems, but can falter when faced with subtle challenges.
5.2 Can LLMs Generate Executable Code or Debug Their Own Code?
The results in Table are with self-debugging for up to rounds, which is essential for obtaining reliable performance from LLM solvers. Without it, many generated programs fail to run, leading significantly worse results. To quantify this, we measure the solver success rate, defined as the percentage of bug-free solver implementations out of all samples. We compare two settings: single-shot generation without debugging vs. multi-round interactions with up to 4 rounds of debugging.
Figure 2 shows the performance difference resulting from debugging. The detailed values can be found Appendix Table and Table 6. For singlo-shot generation, only 41% of solvers run successfully on average across all PDEs. With debugging, success rates rise to 84% on average across PDEs. Notably, frontier models such Gemini 2.5 Pro, Claude 3.7 Somnet, DeepSeek- R1, GPT-4.1, and o3 achieve 90% bug- free rates after self-debugging. We also note that self-debuggin is particularly critical for challenging PDEs such as Compressible Navier-Stokes Equation (CNS), on which only 16.6% of the solvers are bug-free in the single-shot generation. We further report the average number of debug iterations in Appendix Figure 7.
Takcaways. While LLMs frequently produce non-executable code a single attempt, they are remarkably effective at self-debugging. This capability dramatically boosts the rate of generating functional, bug-free solvers, making them reliable for complex scientific coding tasks.

---

Published in Transactions on Machine Learning Research (02/2026)
Table 1: Normalized RMSE comparisons. "Geom. Mean" refers to the geometric mean of aRMSE values. Neural Network baseline results are from (Takamoto et al., 2022; el al., 2022b; al., 2024; Shen et al., 2023; 2024b). The best results of the "Reasoning + Debugging" and the "Reasoning + Debugging + Relinement" settings for each problem are highlighted in gray cells
nRMSE ())
Advection
Burgers
React- Diff
CNS
Darcy
Geom. Mean
Reference Numerical Solvers
1.03 x 10
3.55
10
2.29 10
1.89 x 10 2 4.80 10
2.38 x 10
PDE Solving Softwares
Dedalus PyPDE
6.20 x 10
2.78 x 10 3 5.20 10 2 2.01 x 10
2.03 x 10
5.92 x 10 44
2.21 x 10
2.52 x 10
Neural Network & Foundation Model baselines
U-Net
5.00 x 10
2.20 x 10-1
6.00 10
3.60 x 10-1
FNO
7.70 x 10
7.80 x 10
1.40 x 10
9.50 x 10
9.80 x 10
9.52 x 10
PINN
7.80 x 10-3
8.50 10-1
8.00 x 10-2
ORCA
9.80 x 10-3
1.20 10-2
3.00 10-3
6.20 x 10-2
PDEformer
4.30 x 10
1.46 x 10-2
UPS
2.20 x 10-3
3.73 x 10 5.57 x 10 4.50 x 10
CodePDE: Reasoning + Debugging (best of 32)
Gemini 2.0 Flash
1.14 x 10
2.97 10
2.19 10
4.50 10
4.80 10
6.92 10
Gemini 2.0 Thinking
5.54 x 10-3
3.21 x 10 4
2.19 x 10
7.56 x 10
4.80 x 10
1.07 x 10
Gemini 2.5 Pro
1.01 x 10-3
1.23 x 10
1.49 x 10-1
7.39 x 10
4.89 10 -3
5.82 x 10
Qwwen-2.5-Max QwQ-Plus
x
1.35 10-3
x x
2.36 10
6.76 10
x 10
x x
3.05 x 10-4
2.20 x 10-1
7.43 x 10-2
4.80 x 10-3
7.55 x 10-3
Claude-3.55 3aauh--
3.70 x 10-3
3.23 x 10-4
2.20 x 10-1
2.34 x 10-1
1.00 x 10"
3.61 x 10-2
Claude-3.77 Claudee3.--sonnet
x 10
2.59 10
x 10
4.26 10
4.80 10
5.15 10
DeepSeek-Vi V3
x.37 x 10-3
6.02 X 10
x
2.41 X 10
1.00 10
2.98 X 10
DeepSeek-R1
1.05 x 10-3
2.76 10
2.20 x 10-1
7.46 X 10
4.80 x 10
7.44 X 10-3
GPT-4o
1.55 x 10-3
3.69 10
1.99
10-2
1.81 X 10
7.67 10-1
x 10
GPT--411
1.50 x 10-3
3.63 10-4
1.44 10-1
5.20 X 10-2
4.88
10-3
7.23 x 10-3
o3
9.74 x 10-4
3.69 x 10-4
1.98 X 10-1
2.80 x 10-2
4.92 10-3
6.28 x 10-3
CodePDE: Reasoning + Debugging + Refinement (best of 12)
Gemini 2.0 Flash
8.24 x
1.06 10
1.76 10
1.72 x 10
4.78 10
2.63 x
Gemini 2.0 Thinking
9.74 x 10
2.70 10
1.74
10
2.41 x 10
4.80
10
3.51 x 10
Gemini 2.5 Pro
9.74 x 10
3.15 x 10
9.94 x 10 2
2.77 x 10
4.92 x 10
5.29 x 10
Qwen-2.5-MMx
9.74 x 10
2.60 10
9.13 x 10
x 10
4.80 10
4.40 x 10
QwQ-Plus
9.74 x 10-4
1.07 10-4
5.19
10-2
1.50 x 10
4.80 10-3
3.29 x 10
Claude-3.5hhaike
1.01 x 10-3
2.60 x 10
5.19 x 10
1.06 x 10
1.92 10
x x 10
Claude-3.7-son 3.7-nnne
9.74 x 10
2.60 x 10
5.19 x 10
2.72 x 10
4.83 x 10 3
4.44 x 10
DeepSeek- V3
9.74 10
2.59
10
1.74 10
1.56 10
4.80 10
3.18 10
DeepSeek-BR R1
9.74 x 10
2.59 x 10
1.43 X 10
1.52 x 10
4.80 10
4.83 x 10
GPT-40
9.74 x 10-4
2.59 x 10
1.76
10
2,41 x 10-2
4.80 x 10-3
3.48 x 10-3
GPT-4.1
1.01 x 10-3
3.15
10
1.44
10-1
1.53 x 10
4.88
10-3
5.08 x 10-3
o3
1.01 x 10-3
3.15
10
2.01
10-1
1.31 X 10-2
4.88 10-3
5.27 x 10-3
5.3 Can LLMs Improve Solvers via Celf-Refinement?
Next, we evaluate whether models can refine their solver implementations using nRMSE as the leedback signal. For each PDE, we select the five most accurate solvers generated in the "reasoning + debugging" stage "scod" programs for refinement. Each L.M generates 12 refined programs and we report the best nRMSE, As shown in Table refinement yields substantial improvements. In particular, without refinement, the best LIM generated solver CNS Equation obtains nRMSE of 2.41 x 10-2, underperforming the

---

Published in Transactions on Machine Learning Research (02/2020)
Bug-Free Rates across LLMs
Bug-Free Rates across PDEs
Bug-free Rate (%)
100
93.1
57.5
88.4
93.3
79.9
92.7
762
62.0
52.2
02.0
69.6
66.6
Initial rate Self-debug
45.3
30.5
43.2
40
33.5
20
16.5
Gemini 2.0 Flash
Gemini 2.0 2.00hinking Thinking k
Gemini 2.5 Pro
Qwen-2.5 Max Max
QwQ QwQ-Plus Plus
Claude-3 3.7 3.7-sonnet
3.5-haik haiku
DeepSeek DeepSeek-R1
GPT-40
GPT-4.1
03
Advection
Burgers
CNS
Darcy React-Diff
Claude 3,5
Figure 2: Bug-free rate before and after introducing automated debugging. Left: Averaged across all five PDE datasets for each model. Right: Averaged across LLMs for each PDE family.
Advection
Burgers
React-Diff
CNS
Darcy
Best log nRMSF
3.0
0.6
0.6
-066
-2.2
3.2
-1.0
-2.4
3.44
0.9
1.4
-2.6
3.6
1.0
-1.2
-1.8
-2.8
-388
1.2
-3.0
-1.5
-2.2
20
30
10 20 25
10 15
25
30
30
20 25 Sample size
10
Sample size
Sample size
Sample size
Sample
Gemini 2.5 Pro
QwQ-Plus
Claude-3.7
DeepSeek-R1
63
Figure 3: Smoothed test-time scaling curves for representative models each PDE family. We increase in the best-of-n sampling to study test-time scaling performance.
reference solver, which has nRMSE of 1.89 x 10-2 After refinement, more than half of the LLMs are able to surpass reference performance. Overall, our results show that LLMs effectively coarse feedback signals to optimize solvers.
Interestingly, while advanced reasoning models (e.g,, o3 and DeepSeek-RI) typically lead to better solvers in the "reasoning + debugging" stage, they are not necessarily better than standard ones (e.g., GPT-4o and DeepSeek-V33 V3) in the refinement stage. This divergence suggests that the current post-training strategies for reasoning models may lack coverage on the skills required for refinement.
Takeaways. LLMs can improve code for better accuracy using simple performance feedback. Interestingly, the best models at generating code are not always the best at refining it, suggesting these are two different. skills.
5.4 Does Solution Quality Improve with Scaling Test- Time Compute?
Finally, study the impact of scaling test-time compute. Similar to prior works (Brown et al., 2024; Wu et al., 2025; Snell et al., 2025), we vary E [4, 32] for the best-of-n sampling strategy. The candidate with the lowest nRMSE is selected. The test-time scaling curves for the top-performing models in each family are presented in Figure 3. In general, solution quality generally improves with increasing sample count n, with the most significant gains obscrved between n = 4 and n = 16. Beyond this point, returns diminish, suggesting that moderate sampling budgets often suffice to reach near-optimal performance. Note that because the absolute scale of nRMSE varies across models and tasks, some curves may appear flatter despite meaningful error reductions. To help interprotation, include per-model plots in Appendix Figure 12.
Takeaways. Test-time scaling is critical for generating high-quality numerical solvers from LLMs.
5.5 Convergence Analysis
To assess numerical correctness, we conduct convergence tests that measure how solution error diminishes with increasing grid resolution. A solver is considered convergent if the error decreases at expected rate
8

---

Published in Transactions on Machine Learning Research (02/2026)
Convergence Rate Comparisons across LLMs/PDES
Fail 2.. 19.5% 28.5%%88.144 17..%..7.9...9.....
17.3% 23.6% 8.6%
23.6% 20.3% 13.9%
0.0%
13.8%
1.5%
60.5%
22.5%
75%
Order 1 43.8%%4.%%.629%%%... 58. 44.9% 67.9%
69.4%%%9%%..%%....... 43.3% 43. 81.9%
52.4%% 41.8%
88.0%
53.6%
63.5%
36.7%
48.2%
50%
Order 2 20..%% 9.1% 28.7%
2.7% 14.9% 10.7% 30.6%
6.5%
29.2%
31.0%
34.1%
25%
Order = 3 2.9% 3.5%
1.2%
5.7% 6.8%
1.3%
3.8%
2.6%
2.5%
3.1%
1.9%
1.4%
15.1%
8.6%
1.6%
0.8%
2.5%
7.7%
0%
Geminiz. OFlash Geminiz. OThinking
Gemini2... .5Pro
Qwen Max
QWQ Plus Claude-3, 3, S-Haiku
Claude-3.7-Somi --onnet
3.7-Sonnet
DeepSeek-V3
DeepSeek.R)
GPT-4o mini
GPT-40
GPT-4.1
Advection
Burgers
React-Diff Diff
CNS
Darcy
Claude-3. 3.7dde3377-a
Figure 4: Empirical convergence rate comparisons across LLMs and PDEs. "Fail" represents the failure rate on
the convergence test, Higher orders imply faster convergence.
Table 2: Ablation study on CodePDE and evaluation of existing agentic workflows. results are based Claude-3.7ssonmet for consistent comparison. Full results on 16 LLMs are provided in Appendix Table 4.
nRMSE (4)
Advection
Burgers
React-Diff
CNS
Darcy
Geom. Mean
CodePDE Ablation
CodePDE
9.74 10
2.60 10
5.19 10
2.72
10
4.83 x 10
4.44 x 10
w/o refine
1.32 x 10-3
2.59 x 10-1
5.19 x 10-2
4.26 x 10-2
4.80 x 10
5.15 x 10-3
w/o refine-scale
1.37 x 10-3
4.60 x 10-3
4.26 x 10-1
3.02 x 10-1
6.06 x 10-1
8.68 x 10-2
w/o delin---cale-deluug
2.16 10-3
1.00
10-2
7.12 10-1
6.62 10-1
7.16 x 10-1
1.49 x 10-1
Evaluation of Existing Agentic Workflows
FunSearch AIDE
9.84 x 10 4
2.59 x 10 4
5.12 x 10 10-2 3.79 x 10 2
4.80 x 10 3
4.73 x 10 3
1.03 10-3
1.05
10-1 5.07
10-2 5.77 x 10-2
4.78 x 10-3
4.33 x 10-3
as the grid is refined. Figure 4 reports both the failure rate of the convergence test and the distribution of convergence orders across different LLMs and PDE families, where higher orders imply faster convergence.
The results reveal trade-off between solver reliability and sophistication. For example, GPT-4o mini achieves very low failure rate (8.6%) by conservatively defaulting to simple first-order methods. In contrast, Gemini 2.5 Pro and o3 demonstrate a more advanced capability, generating a significantly higher proportion of higher order methods while maintaining reasonably low failure rates.
From the problem perspective, the CNS Equation proves challenging with a high average failure rate of 60.5%. The result underscores that single-shot generation can be insufficient. for complex problems, reinforeing the need for test-time compute scaling (Section 5.4).
Takeaways. A low failure rate alone is not a sufficient measure of a model's capability; the ability to generate diverse, high-order numerical methods is also critical. For challenging PDFs where single-shot generation generations are prone to failure, test-time scaling essential for obtaining diverse, correct, and robust solvers.
5.6 Ingredients for Accurate PDE Solving: Insigts into Agent Designs
To quantitatively understand the key drivers of performance in LLM-driven PDE solver generation framework, we conduct an ablation study on the core components of CodePDF. The results are presented Table 2.
The ablation study clearly demonstrates that each component refinement, scaling, and debugging is critical for achieving high accuracy. Removing these components progressively degrades performance of the full CodePDE. A naive prompting approach that lacks all three components yields an overall nRMSE of 1.49 x 10-1 and leads to unusable solvers, especially for challenging PDEs. This result confirms the value of a well-designed inference framework.
We further investigate whether these components are universally beneficial by adapting two recent agentic workflows to our tasks, FunSearch (Romera-Paredee et al., 2024) and AIDE (Jiang al., 2025). FunSearch employs evolutionary search that naturally incorporates refinement and scaling; augment it with
9

---

Published in Transactions on Machine Learning Research (02/2026)
Spatial Discretization Distribution
100%
19%
12%
80%
34%
Proportion
50%
91%
100%
90%
40%
06%
53%
20%
40%
6.5
31%
12%
0%
Gemini 2.0 Flash Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.55 Owen-2.5.Max M.x
QwQ Plus
Claude Claude-3.5-haiku h.5-h
Claude-3. Claude-3.7-sae 3.7-sonnet
DeepSeek-V3
DeepSeek-R1
GPT-4o mini
GPT-4o
GPT-4.1
03
o1-minl
03 03-minl
04 o4-mini
Finite Difference
Finite Volume
Spectral
Figure 5: Distributions of the spatial discretization schemes employed by each LLM.
Time Integration Distribution
100%
6%
80%
38%
2.3%
29%
Proportion
0%
62%6
60%
100%
62%
9%
075
37%
100%
40%
88%
84%
5.9
62%
20%
25%
37%
22%
Gemini 2.0 Flash
Gemini 2.0 Thinking
0%
Gemini 2.5 Pro
Owen-2.5-Mar 5-xx
Claude Claude-3.5-haid h.5-u
QwQ-Plus
Claude-3. 7-sonnet
DeepSeekV3
DeepSeek-R1
GPT-4o mini
GPT-40
GPT-4.1
ol ol-mini
03 03-mini
04 04-mini
Explicit Euler
RK2
RK3-SSP
RK4
IMEX
Others
Figure 6: Distributions of the time integration schemes employed by each LLM. "RK2", "RK3-SSP", "RK4", and "IMEX" refers to Runge-Kon Kutta 2rd order, Kunge-Kutta 3rd order (strong stability preserving), Runge- Kutta 4th order, and Implicit-Explicr schemes, respectively.
a debugging mechanism. AIDE implements a tree search and already features all three capabilities. As shown in Table 2, once equipped with these essential components, both FunSearch (4.73 x 10-3) and AIDE (4.33 10-3) achieve performance comparablc CodePDE. This suggests that the presence of these core functionalities is more important. than the specific high-level search strategy.
Takeaways. There is substantial performance gap between naive prompting and using a structured inference framework. The combination of self-refinement, test-time scaling, and automated debugging is essential for LLM-driven PDE solving workflow /agent designs.
5.7 Avenues for Interpretability: Insights into the Generated Solvers
In this subsection, investigate the generated solvers and make few observations.
Numerical scheme diversity. We conduct detailed analysis of the numerical schemes employed in the generated solvers for Burgers Equation. Figures 5 and 6 illustrate the distribution of spatial discretization and time integration methods. Overall, the results show a reasonably broad coverage of different schemes, although most LLMs predominantly adopt simple finite differenne and explicit Euler methods. Notably, DeepSeek-R1 diverges from this trend, lavoring spectral discretization and IMEX integration. We also find that o3 exhibits the most diversity in its scheme choices, indicating better exploration ability. More detailed discussions, well the distributions of convective flux schemes and timc-stepping strategics, are provided in Appendix C.
Understanding the failure on Reaction-Diffusion Equation. Despite strong overall performance, all LLMs and agentic methods consistently struggle with the Reaction-DDilfusionn Equation (Tables and 2). Upon inspection of the LLM generated code and the reference implementation, we identify a key distinction: LLM-generatod solvers typically apply finitoddifeerence schemes the reaction term, whereas the reference solver exploits the existence of an analytical solution for this component. We provide detailed comparisons and include relevant solver code in Appendix D. (Table 9).
10

---

Published in Transactions on Machine Learning Research (02/2026)
This analysis highlights advantage of LLM-driven PDE solver generation: solver transparency and interpretability. Unlike neural solvers, LLMs produce human-readable code that exposes the model's s reasoning. This interpretability enables us to diagnose errors such as incorrect discretizations, improper boundary conditions, or unstable time integration- and opens the door to human-in-the-loop corrections or targeted model refinement.
6 Conclusions
We introduce CodePDE, an inference framework for LLM-driven PDE solving. The investigation reveals key insights into the capabilitics of LLMs for scientific computing. We find that a structured inference framework with the combination of automated debugging, refinement, and test-time scaling is essential. We also uncover trade-off between solver reliability and sophistication and show that frontier models can balance solver correctness and diversity in their exploration, desired property for effective test-time scaling and even scientific discovery. We believe these insights pave the way for developing more capable LLMs and agents for scientific computing.
Acknowledgments
The authors would like to thank Nikola B. Kovachki (NVIDIA), Yiping Lu (Northwestern University), Wuyang Chen (Simon Fraser University), Chuwei Wang (Caltech), and Haowei Lin (Peking University) for the helpful discussions.
This work was supported in part by the National Science Foundation grants IIS1705121, IIS1838017, IIS2046613, IIS2112471. and funding from Datadog. Any opinions, findings and conclusions recommenda- tions expressed in this material are those of the authors and do not necessarily reflect the views of any of these funding agencies.
11

---

Published in Transactions on Machine Learning Research (02/2026)
References
Anthropic. Introducing the next generation of claude, 2024. URL claude-3- f-mmily..
Sören Arlt, Haonan Duan, Felix Li, Sang Michael Xie, Yuhuai Wu, and Mario Krenn. Meta-designing quantum experiments with language models, 2024. URL https://arxiv. wttp///w....bb.....//oo 02470.
Andres M Bran, Sam Cox, Oliver Schilter, Carlo Baldassari, Andrew D White, and Philippe Schwaller. Chemcrow: Augmenting langaale models with chemistry tools. arXiv preprint arXiv:2304.05376 2023.
Johannes Brandstetter, Daniel E. Worrall, and Max Welling. Message passing neural PDE solvers. In International Conference on Learning Representations, 2022. URL https://openrreviex. http:://ppparreeeeeeeeeeetuuuut net/forum?id vSix3HPYKSU.
Bradley Brown, Jordan Juravsky, Ryan Ehrlich, Ronald Clark, Quoc V Le, Christopher Ré, and Azalia Mirhoseini. Large language monkeys: Scaling inference compute with repeated sampling. arXiv preprint arXiv:2407.21887 2024.
Ricardo Buitrago, Tanya Marwah, Albert Gu, and Andrej Risteski. On the benefits of memory for modeling time-dependent PDEs. In The Thirleenth International Conference on Learning Representations, 2025. URL attpss//openreviiee wttp:/errrreeeeetttefeooomiaoooe ret/ffrrum2dd==00g55833
Keaton J Burns, Geoffrey M Vasil, Jeffrey S Oishi, Daniel Lecoanet, and Benjamin P Brown. Dedalus: A flexible framework for numerical simulations with spectral methods. Physicaa Review Research, 2(2):023068, 2020.
Claudio Canuto, M Yousuff Hussaini, Alfio Quarteroni, and Thomas A Zang. Spectral methods: evolution to complez geometries and applications to fluid dynamics. Springer Science & Business Media, 2007.
Shuhao Cao. Choose transformer: Fourier or galerkin. Advances in neural information processing systems, 34:24924-2490 2021.
Jun Shern Chan, Neil Chowdhury, Oliver Jafse, James Aung, Dane Sherburn, Evan Mays, Giulio Starace, Kevin Liu, Loon Maksin, Tejal Patwardhan, Alcksander Madry, and Lilian Weng. MLFh-benhh Evaluating machine learning agents on machine learning engineering, In The Thirteenth International Conference on Learning Representations, 2025. URL https://openrevierr wwt. /to.f/orrrmd---5500W11
Xinyun Chen, Maxwell Lin, Nathansel Schärli, and Denny Zhou. Teaching large language models to self-debug. In The Twelfth International Conference on Learning Representations, 2024. URL https: //openreview.c wttp:///oooormiddd.i........
Wenduo Cheng, Junhong Shen, Mikhail Khodak, Jian Ma, and Ameet Talwalkar. L2g: Repurposing language models for genomics tasks. bioRrie, 2024. doi: 10.1101 1111012202..12.0.6222 ///pw.biorrk/ wwwlvvooggoooooeooo http:/oonn////r////22222 11///////1111 1//11///11///////// 12. 29..99.622. 627422.
URL https:
AI, Aixin Liu, Bei Feng, Bing Xue, Bing-Li Wang, Bochao Wu, Chengda Lu, Chenggang Zhao, Chengqi Deng, Chenyu Zhang, Chong Ruan, Damai Dai, Daya Guo, Dejian Yang, Deli Chen, Dong-Li Li Ji, Erhang Li, Fangyun Lin, Fucong Dai, Fuli Luo, Guangbo Hao, Guanting Chen, Guowei Li, H. Zhang, Han Bao, Hanwei Xu, Haocheng Wang, Haowei Zhang, Honghui Ding, Huajian Xin, Huazuo Gao, Hui Li, Hui Qu, J. L. Cai, Jian Liang, Jianzhong Guo, Jiaqi Ni, Jiashi Li, Jiawei Wang, Jin Chen, Jingchang Chen, Jingyang Yuan, Junjie Qiu, Junlong Li, Jun-Mei Song, Kai Dong, Kai Hu, Kaige Gao, Kang Guan, Kexin Huang, Kuai Yu, Lean Wang, Lecong Zhang, Lei Xu, Leyi Xia, Liang Zhao, Litong Wang, Liyue Zhang, Meng Li, Miaojun Wang, Mingchuan Zhang, Minghua Zhang, Minghui Tang, Mingming Li, Ning Tian, Panpan Huang, Peiyi Wang, Peng Zhang, Qiancheng Wang, Qihao Zhu, Qinyu Chen, Qiushi Du, R. J. Chen, Jin, Ruiqi Ge, Ruisong Zhang, Ruizhe Pan, Runji Wang, Runxin Xu, Ruoyu Zhang, Ruyi Chen, S. S. Li, Shanghao Lu, Shangyan Zhou, Shanhuang Chen, Shao-Ping Wu, Shengfeng Ye, Shirong Ma, Shiyu Wang, Shuang Zhou, Shuiping Yu, Shunfeng Zhou, Shuting Pan, T. Wang, Tao Yun, Tian Pei, Tianyu Sun,
12

---

Published in Transactions on Machine Learning Research (02/2026)
W. L. Xiao, Wangding Zeng, Wanjia Zhao, Wei An, Wen Liu, Wenfeng Liang, Wenjum Gao, Wen-Xuan Yu, Wentao Zhang, X. Q. Li, Xiangyu Jin, Xianzu Wang, Xiaoling Bi, Xiaodong Liu, Xiaohan Wang, Xi-Cheng Shen, Xiaokang Chen, Xiaokang Zhang, Xiaosha Chen, Xiaotao Nie, Xiaowen Sun, Xiaoxiang Wang, Xin Chong, Xin Liu, Xin Xie, Xingchao Liu, Xingkai Yu, Xinnan Song, Xinxia Shan, Xinyi Zhou, Xinyu Yang, Xinyuan Li, Xuecheng Su, Xuheng Lin, Y. K. Li, Y. Q Wang, Y. X. Wei, Y. X. Zhu, Yang Zhang, Yanhong Xu, Yanping Huang, Yao Li, Yao Zhao, Yaofeng Sun, Yao Li, Yaohui Wang. Yi Yu, Yi Zheng, Yichao Zhang, Yifan Shi, Yi Xiong, Ying He, Ying Tang, Yishi Piao, Yisong Wang, Yixuan Tan, Yi-Bing Ma, Yiyuan Liu, Yongqiang Guo, Yu Wu, Yuan Ou, Yuchen Zhu, Yuduan Wang, Yue Gong, Yuheng Zou, Yujia He, Yukun Zha, Yunfan Xiong, Yunxiang Ma, Yuting Yan, Yu-Wei Luo, Yu mei You, Yuxuan Liu, Yuyang Zhou, Z. F. Wu, Zehui Ren, Zehui Ren, Zhangli Sha, Zhe Fu, Zhean Xu, Zhen Huang, Zhen Zhang, Zhenda Xie, Zhen guo Zhang, Zhewen Hao, Zhibin Gou, Zhicheng Ma, Zhigang Yan, Zhihong Shao, Zhipeng Xu, Zhiyu Wu, Zhongyu Zhang, Zhuoshu Li, Zihui Gu, Zijia Zhu, Zijum Liu, Zi-An Li, Ziwei Xie, Ziyang Song, Ziyi Gao, and Zizheng Pan. Deepscck-v3 technical report. ArXiv, abs/2412.19433 2024.
DeepSock-AT, Daya Guo, Dejian Yang, Haowei Zhang, Jun-Mei Song, Ruoyu Zhang, Runxin Xu, Qihao Zhu, Shirong Ma, Peiyi Wang, Xiaoling Bi, Xiaokang Zhang, Xingkai Yu, Yu Wu, Z. F. Wu, Zhibin Gou, Zhihong Shao, Zhuoshu Li, Ziyi Gao, Aixin Liu, Bing Xue, Bing-Li Li Wang, Bochao Wu, Bei Feng, Chengda Lu, Chenggang Zhao, Chengqi Deng, Chenyu Zhang, Chong Ruan, Damai Dai, Deli Chen, Dong-Li Ji, Erhang Li, Fangyun Lin, Fucong Dai, Fuli Luo, Guangbo Hao, Guanting Chen, Guowei Li, H. Zhang, Han Bao, Hanwei Xu, Haocheng Wang, Honghui Ding, Huajian Xin, Huazuo Gao, Hui Qu, Hui Li, Jianzhong Guo, Jiashi Li, Jiawei Wang. Jingchang Chen, Jingyang Yuan, Junjie Qiu, Junlong Li, Jiong Cai, Jiaqi Ni, Jian Liang, Jin Chen, Kai Dong, Kai Hu, Kaige Gao, Kang Guan, Kexin Huang, Kuai Yu, Lean Wang, Lecong Zhang, Liang Zhao, Litong Wang, Liyue Zhang, Lei Xu, Leyi Xia, Mingchuan Zhang, Minghua Zhang, M. Tang, Meng Li, Miaojun Wang, Mingming Li, Ning Tian, panpan Huang, Peng Zhang, Qiancheng Wang, Qinyu Chen, Qiushi Du, Ruiqi Ge, Ruisong Zhang, Ruizhe Pan, Runji Wang, R. J. Chen, R. L. Jin, Ruyi Chen, Shanghao Lu, Shangyan Zhou, Shanhuang Chen, Shengfeng Ye, Shiyu Wang. Shuiping Yu, Shunfeng Zhou, Shuting Pan, S. S. Li, Shuang Zhou, Shao-Kang Wu, Tao Yun, Tian Pei, Tianyu Sun, T. Wang, Wangding Zeng, Wanjia Zhao, Wen Liu, Wenfeng Liang, Wenjum Gao, Wen-Xia Yu, Wentao Zhang, W. L. Xiao, Wei An, Xiaodong Liu, Xiaohan Wang, Xiaokang Chen, Xiaotao Nie, Xin Cheng, Xin Liu, Xin Xie, Xingchao Liu, Xinyu Yang, Xinyuan Li, Xuecheng Su, Xuheng Lin, X. Q. Li, Xiangyu Jin, Xi-Cheng Shen, Xiaosha Chen, Xiaowen Sun, Xiaoxiang Wang, Xinnan Song, Xinyi Zhou, Xianzu Wang, Xinxia Shan, Y. K. Li, Y. Q Wang, Y. X. Wei, Yang Zhang, Yanhong Xu, Yao Li, Yao Zhao, Yaofeng Sun, Yaohui Wang, Yi Yu, Yichao Zhang, Yifan Shi, Yi Xiong, Ying He, Yishi Piao, Yisong Wang, Yixuan Tan, Yiyang Ma, Yiyuan Liu, Yongqiang Guo, Yuan Ou, Yuduan Wang, Yue Gong, Yu-Jing Zou, Yujia He, Yunfan Xiong, Wei Luo, Yu mci You, Yuxuan Liu, Yuyang Zhou, Y. X. Zhu, Yanping Huang, Yao Li, Yi Zheng, Yuchen Zhu, Yunxiang Ma, Ying Tang, Yukun Zha Yuting Yan, Zehui Ren, Zehui Ren, Zhangli Sha, Zhe Fu, Zhean Xu, Zhenda Xie, Zhen guo Zhang, Zhewen Hao, Zhicheng Ma, Zhigang Yan, Zhiyu WU, Zihui Gu, Zijia Zhu, Zijun Liu, Zi-An Li, Ziwei Xic, Ziyang Song, Zizhong Pan, Zhen Huang, Zhipeng Xu, Zhongyu Zhang, and Zhen Zhang. Deepseek-r1: Incentivizing reasoning capability in Iims via reinforcement learning. ArXin, abs/2501.1 12948, 2025.
Shengyu Feng, Weiwei Sun, Shanda Li, Ameet Talwalkar, and Yiming Yang. FrontierCO: Real-world and large- scale evaluation of machine learning solvers for combinatorial optimization. In The Fourteenth International Conference on Learning Representations, 2026. URL https:
Google. Gemini 2.0 flash thinking, 2025. URL https: http////wmm..gggggge//////////// flash-thinking/.
Taicheng Guo, Kehan Guo, Bozhao Nan, Zhenwen Liang, Zhichun Guo, Nitesh V Chawla, Olaf Wiest, and Xiangliang Zhang. What can large language models do in chemistry? a comprehensive benchmark on eight. tasks. In Thirth-seventh Conference on Neural Information Processing Systems Datascts and Benchmarks Track, 2023. URL
Robert Haase, Christian Tischer, Jcan-Karim Hériché, and Nico Scherf. Benchmarking large language models for bio-image analysis code generation. bioRriv, 2024. doi: 1111101/2224.44119100 URL https://wwe https:///ww.bborxiv www..iiioog/oooneeeeeeeeeeeeeeee /ttp:///teee////aa////////////// 19900009
13

---

Published in Transactions on Machine Learning Research (02/2026)
Zhongkai Hao, Chang Su, Songming Liu, Julius Berner, Chengyang Ying, Hang Su, Anima Anandkumar, Jian Song, and Jun Zhu. DPOT: Auto-regressive denoising operator transformer for large-scale PDE pre- training. In Forty-first International Conference on Machine Learning, 2024. URL https://openrevial net/forum? iid X7UnDevHOM.
Di He, Shanda Li, Wenlei Shi, Xiaotian Gao, Jia Zhang, Jiang Bian, Liwei Wang, and Tie- Yan Liu. Learning physics-infornn neural networks without stacked pacpapropp.ation. In Francisco Ruiz, Jennifer Dy, and Jan- Willem van de Meent (eds.), Proceedings of The 26th International Conference on Artificial Intelligence and Statistics, volume 206 of Proceedings of Machine Learning Research, pp. 3034-3044 PMLR, 25-27 Apr 2023. URL pettpa:/prrceediiee
mlr. wwt.pprees/vv0 wttp///222/e....//////// html.
Maximilian Herde, Bogdan Raonic, Tobias Rohner, Roger Käppeli, Roberto Molinaro, Emmanuel de Bézenac, and Siddhartha Mishra. Poseidon: Elficient foundation models for pdes. Advances in Neural Information Processing Systems, 37:72525 72624, 2024.
Pengfei Hong, Navonil Majumder, Dcepanway Ghosal, Somak Aditya, Rada Mihalcea, and Soujanya Poria. Evaluating 1lms' mathematical and coding compelency through onlology-guided interventions, 2024. URL
https://arxiv. http://axxiv.../abbbba owg/abs/////eeeeoooeooooooee ab//a/aa2...... 2ttp//2............
Zhengyao Jiang, Dominík Schmidt, Dhruv Srikanth, Dixing Xu, Ian Kaplan, Deniss Jacenko, and Yuxiang Wu. Aide: Ai-driven exploration in the space of code. arXiv preprint ar Xit::5522.1333, 2025.
Dmitrii Kochkov, Jamie A. Smith, Ayya Alieva, Qing Wang, Michael P. Bremner, and Stephan Hoyer. Machine learning-acceeeraaeel accelerated computational fluid dynamics. Proceedings of the National Academy of Sciences, 118 (21), 2021. ISSN 0027-8424. doi: 10.1073/pmae 21017888.. URL https: //www///wwwww/// pnas. org/content/118, 21/e2101784188.
Randall J LeVeque. Finite difference methods for ordinary and partial differential equations: steady-state and time-dependent problems. SIAM, 2007.
Shanda Li, Shinjae Yoo, and Yiming Yang. Maximal update parametrization and zero-shot hyperparameter transfer for fourier neural operators. In Forly-second International Conference on Machine Learning, 2025.
URL
Sijie Li, Weiwei Sun, Shanda Li, Ameet Talwalkar, and Yiming Yang. Comind: Towards community- driven agents for machine learning engineering. In The Fourteenth International Conference on Learning Representations, 2026. URL https: /ttp::////pereeeeeeeee http://wwwveeeeeeeeeoor//eeeeeee nett/foruum/dd=4777000000
Zijie Li, Kazem Meidani, and Amir Barati Farimani. Transformer for partial differential equations' operator
learning. Xriio preprint arXiv:2205. 116,11 2022a.
Zongyi Li, Nikola Kovachki, Kamyar Aziwadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Neural operator: Graph kernel network for partial differential equations. arXiv preprint XrXi:22033.00 03485, 2020.
Zongyi Li, Daniel Zhengyu Huang, Burigede Liu, and Anima Anandkumar. Fourier neural operator with learned deformations for pdes on general geometries. arXie preprint wXXi:2207.05200 2022b.
Weixin Liang, Junhong Shen, Genghan Zhang, Ning Dong, Luke Zettlemoyer, and Lili Yu. Mixture-of- mamba: Enhancing multi-modal state-space models with modality sparsity, 2025. URL https:
//arxivvorg/abs//5000 16295.
Zichao Long, Yiping Lu, Xianzhong Ma, and Bin Dong. PDE-net: Learning PDEs from data. In Jennifer Dy and Andreas Krause (eds.), Proceedings of the 35th International Conference on Machine Learning, volume 80 of Proceedings of Machine Learning Research, pp. 3208-322 PMLR, 10-15 Jul 2018. URL
https://proceedings mlr. press/v80/longg html.
Lu Lu, Pengzhan Jin, Guofei Pang. Zhongqiang Zhang, and George Em Karniadakis. Learning nonlincar operators via deeponet based on the universal approximation theorem of operators. Nature machine intelligence, 3(3):218 2211822292222 2021.
14

---

Published in Transactions on Machine Learning Research (02/2026)
Pingchuan Ma, Tsun-Hsuan Wang, Minghao Guo, Zhiqing Sun, Joshua B. Tenenbaum, Daniela Rus, Chuang Gan, and Wojciech Matusik. Llm and simulation as bilevel optimizers: A new paradigm to advance physical scientific discovery, ArXiv, abs/2405.097883 2024. URL https://apt. cttps://ap..emmatiiccoolllr semantccacholar I:: 269790814.
Aman Madaan, Niket Tandon, Prakhar Gupta, Skyler Hallinan, Luyu Gao, Sarah Wiegrefle, Uri Alon, Nouha Dziri, Shrimai Prabhumoye, Yiming Yang, Shashank Gupta, Bodhisattwa Prasad Majumder, Katherine Hermann, Scan Welleck, Amir Yazdanbakhsh, and Peter Clark. Solf-refine: Iterative refinement with self-feedback. In Thirby-seventh Conference on Newral Information Processing Systems, 2023. URL https://openreviiee http://ppeeeeriee...eoooorrorddd
Tanya Marwah, Ashwini Pokle, Zico Kolter, Zachary Chase Lipton, Jianfeng Lu, and Andrej Risteski. Deep equilibrium based neural operators for steady-state PDEs. In seirty- Conference on Neural Information Processing Systems, 2023. URL https: https://openreviex
Michael McCabe, Bruno Régaldo-Saint Blancard, Liam Holden Parker, Ruben Ohana, Miles Cranmer, Alberto Bietti, Michael Eickenberg, Siavash Golkar, Geraud Krawezik, Francois Lanusse, Mariel Pettee, Tiberiu Tesileanu, Kyunghyun Cho, and Shirley Ho. Multiple physics pretraining for spatiotemporal surrogate models. In The Thirty-eighth Annual Conference on Neural Information Processing Systems, 2024. URL http///eeperreeeeeeeeee/ee rett/forummitttxSSSBBU...
Mike A Merrill, Alexander Glenn Shaw, Nicholas Carlini, Boxuan Li, Harsh Raj, Ivan Bercovich, Lin Shi, Jeong Yeon Shin, Thomas Walshe, E. Kelly Buchanan, Junhong Shen, Guanghao Ye, Haowei Lin, Jason Poulos, Maoyu Wang, Jenia Jitsev, Marianna Nezhurina, Di Lu, Orfeas Menis Mastromichalakis, Zhiwei Xu, Zizhao Chen, Yue Liu, Robert Zhang, Leon Liangyu Chen, Anurag Kashyap, Jan-Lucaa Uslu, Jeffrey Li, Jianbo Wu, Minghao Yan, Song Bian, Vedang Sharma, Ke Sun, Steven Dillmann, Akshay Anand, Andrew Tanpouthakoun, Bardia Koopah, Changran Hu, Etash Kumar Guha, Gabriel H. S. Dreiman, Jiacheng Zhu, Karl Krauth, Li Zhong, Niklas Muennighoff, Robert Kwesi Amanfu, Shangyin Tan, Shreyas Pimpalgaonkar, Tushar Aggarwal, Xiangning Lin, Xin Lan, Xuandong Zhao, Yiqing Liang, Yuanli Wang, Zilong Wang, Changzhi Zhou, David Heineman, Hange Liu, Harsh Trivedi, John Yang, Junhong Lin, Manish Shetty, Michael Yang, Nabil Omi, Negin Raoof, Shanda Li, Terry Yue Zhuo, Wuwei Lin, Yiwei Dai, Yuxin Wang, Wenhao Chai, Shang Zhou, Dariush Wahdany, Ziyu She, Jiaming Hu, Zhikang Dong, Yuxuan Zhu, Sasha Cui, Ahson Saiyed, Arinbjörn Kolbeinsson, Christopher Michael Rytting, Ryan Marten, Yixin Wang, Alex Dimakis, Andy Konwinski, and Ludwig Schmidt. Terminal-bench: Benchmarking agents on hard, realistic tasks in command line interfaces. In The Fourteenth International Conference on Learning Representations, 2026. URL
AL. Gpt-4o, 2024a. URL http:://pppaaaa.ccc//dddd///lllp
OpenAL. Openai o3-mini, 2024b. URL
Tobias Pfaff, Meire Fortunato, Alvaro Gonzaez- and Peter Battaglia. Learning mesh-based simulation with graph networks. In International Conference on Leorning Representations, 2021. URL https: //openreview. boopearrvvii.nne//fruu/itii net//forum?.ddroogg.))) _PP
Maziar Raissi, Paris Perdikaris, and George Karniadakis. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlincar partial differential equations. Journal of Compulational Physics, 378:686-700) 2019.
Nicholas Roberts, Samuel Guo, Cong Xu, Ameet Talwalkar, David Lander, Lvfang Tao, Linhang Cai, Shuaicheng Niu, Jianyu Heng, Hongyang Qin, Minwen Deng, Johannes Hog. Alexander Pfefferle, Sushil Ammanaghatte (a Shivakumar, Arjun Krishnakumarr Yubo Wang, Rhea Sanjay Sukthanker, Frank Hutter, Euxhen Hasanaj, Tien-Dung Le, Mikhail Khodak, Yuriy Nevmyvaka, Kashif Rasul, Frederic Sala, Anderson Schneider, Junhong Shen, and Evan R. Sparks. Automl decathlon: Diverse tasks, modern methods, and efficiency at scale. In Neural Information Processing Systems, 2021. URL https: //api. wesmamticccholaarooo///oooaooeee :2255366.
15

---

Published in Transactions on Machine Learning Research (02/2026)
Bernardino Romera-Paredes, Mohammadamin Barekatain, Alexander Novikov, Matej Balog, M Pawan Kumar, Emillie Dupont, Francisco JR Ruiz, Jordan S Ellenberg, Pengming Wang, Omar Fawzi, et al. Mathematical discoveries from program search with large language models. Nature, 625(7995):468 475, 2024.
Olaf Ronneberger, Philipp Fischer, and Thomas Brox. U-net: Convolutional networks for biomedical image segmentation. In International Conference Medical image compuling and compuler-ussisled intervenlion, pp. 234 241. Springer, 2015.
Nadim Saad, Gaurav Gupta, Shima Alizadeh, and Danielle C. Maddix. Guiding continuous operator learning through physics-based boundary constraints. In The Eleventh International Conference on Learning Representations, 2023. URL http:p:://aaameeiooooee////aaaaa ner//oorun/idrggfM0G000.
Matthew W Scroggs, Igor A Baratta, Chris N Richardson, and Garth N Wells. Basix: a runtime finite element. basis evaluation library. Journal of Open Source Software, 7(73):3982 2022.
Junhong Shen, Mikhail Khodak, and Ameet Talwalkar. Efficient architecture search for diverse tasks. In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho (eds.), Advances in Neural Information Processing Systems, 2022. URL https: http://vpeeeeeeeeeeeem htte://fommmiddd/m/aaa
Junhong Shen, Liam Li, Lucio M. Dery, Corey Staten, Mikhail Khodak, Graham Neubig, and Ameet Talwalkar. Cross-modal fine-tuning: align then refine. In Proceedings of the 40th International Conference Machine Learning, 2023.
Junhong Shen, Atishay Jain, Zedian Xiao, Ishan Amlekar, Mouad Hadji, Aaron Podolny, and Ameet Talwalkar. Scribeagent: Towards specialized web agents using production-scale workflow data, 2024a. URL https://arxiv http::///aaav..org/aa//22111 arg/abs/2411. 15004.
Junhong Shen, Tanya Marwah, and Ameet Talwalkar. UPS: Efficiently building foundation models for PDE solving via modll adaptation. Transactions on Machine Learning Research, 2024b. ISSN 2835-8856. URL httpp:://ppneeeee
Junhong Shen, Neil Tenenholtz, James Brian Hall, David Alvarez-Melis, and Nicolo Fusi. Tag-LLM: Repurposing general-purpose LLMs for specialized domains. In Forty-first International Conference on Machine Learning, 2024c. URL nett/forumiiddOrnnit.
Junhong Shen, Kushal Tirumala, Michihiro Yasunaga Ishan Misra, Luke Zettlemoyer, Lili Yu, and Chunting Zhou. Cal: Content-adaplive image lokenization, 2025. URL
Charlie Vietor Snell, Jaehoon Lee, Kelvin Xu, and Aviral Kumar. Scaling LLM test-time compute optimally can be more effective than scaling parameters for reasoning. In The Thirteenth International Conference on Learning Representations, 2025. URL https://openreview http://wwwweeeeefrereeeee
net./torumPldd44WAAZZ..0
Mauricio Soroco, Jialin Song, Mengzhou Xia, Kye Emond, Weiran Sun, and Wuyang Chen. PDE-controller: LLMs for autoformalization and reasoning of PDEs. In Forty-second International Conference on Machine Learning, 2025. URL
Shashank Subramanian, Peter Harrington, Kurt Keutzer, Wahid Bhimji, Dmitriy Morozov, Michael W. Mahoney, and Amir Gholami. Towards foundation models for scientific machine learning: Characterizing scaling and transfer behavior. In Thirly-seventth Conference on Neural Information Processing Systems, 2023. URL. btetps://operreee wwtpeerreieeeneetffeeeeeeee net:/forum?:d=zzAlavzZiil
Weiwei Sun, Shengyu Feng, Shanda Li, and Yiming Yang. Co-bench: Benchmarking language model agents in algorithm search for combinatorial optimization. arXiv preprint arXiv:254.003310, 2025.
Makoto Takamoto, Timothy Praditia, Raphael Leiteritz, Daniel MacKinlay, Francesco Alesiani, Dirk Pflüger, and Mathias Niepert. Pdebench: An extensive benchmark for scientific machine learning. Advonces in Neural Information Processing Systems, 35:1596-1611 2022.
16

---

Published in Transactions on Machine Learning Research (02/2026)
Xiangru Tang, Bill Qian, Rick Gao, Jiakang Chen, Xinyun Chen, and Mark Gerstein. Biocoder: A benchmark for bioinformaties code generation with large language models, 2024. URL settps://arxiv.oor Btrgg/bbs//n 2308. 16458.
Gemini Team, Rohan Anil, Sebastian Borgeaud, Jean-Bapt iste Alayrac, Jiahui Yu, Radu Soricut, Johan Schalkwyk, Andrew M Dai, Anja Hauth, Katie Millican, et al. Gemini: a family of highly capable multimodal models. arXiv preprint arXiv:2312. 11805, 2023.
Qwen Team. Qwen2.5: A party of foundation models, September 2024. URL https:///uenlnnn wea.llm.iittuu.... / blog/quen2.5/.
Qwen Team. Qwq-32b: Embracing the power of reinforcement learning, March 2025. URL https:://q// httbub..o/blog//rgg//// 32b/.
Renbo Tu, Nicholas Roberts, Mikhail Khodak, Junhong Shen, Frederic Sala, and Ameet Talwalkar. NAS- bench-360: Benchmarking neural architecture search on diverse tasks. Thirty-siath Conference on Neural Information Processing Systems Dalasels and Benchmarks Track, 2022. URL https://openreview..rt formm?id=xUXTooggWB..
Chuwei Wang, Shanda Li, Di He, and Liwei Wang. Is S1 physics informed loss always suitable for training physics informed neural network? In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho (cds.), Advances in Newral Information Processing Systems, 2022. URL https: fttt::://openevieee /ttpp//rrrree..t///orrnrmr net/forum? id=cy1TKLRAEM.
Ke Wang, Houxing Ren, Aojun Zhou, Zimu Lu, Sichun Luo, Weikang Shi, Renrui Zhang, Lingi Song, Mingjie Zhan, and Hongsheng Li. Mathcoder: Seamless code integration in llms for euhanced mathematical reasoning, 2023. URL https: /ttp:////axivoooo org/ aarg/bbs/3100.... 03731.
Jason Wei, Xuezhi Wang, Dale Schuurmans, Maarten Bosma, brian ichter, Fei Xia, Ed H. Chi, Quoc V Le, and Denny Zhou. Chain of thought prompting elicits reasoning in large language models. In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho (eds.), Advances in Neural Information Processing Systems, 2022. URL https: / /ttpr/////p/eeeeeeeeeeee wtwperrwwwie.eeeffrunnn... net/forummid= lttrfoorumnnddiiiiedeead _ijQMMeBB__..
Sean Welleck, Amanda Bertsch, Matthew Finlayson, Hailey Schoelkopf, Alex Xie, Graham Neubig, Tlia Kulikov, and Zaid Harchaoui. From decoding to meta-generation: Inference-time algorithms for large language models. Transactions on Machine Learning Research, 2024. ISSN 2835-8856. URL https:
//openreview. ///ppe/reviei..nt///ttrt/t/////e net./EoruumiddssskMcc..aa
Survey Certification.
Yangzhen Wu, Zhiqing Sun, Shanda Li, Sean Welleck, and Yiming Yang. Inference scaling laws: An empirical analysis of compute-optimal inference for L.M problem-solving. In The Thirteenth International Conference on Learning Representations, 2025. URL https://openrevview www. settforrumidd=WWikpp JEHn.
Zhanhong Ye, Xiang Huang, Leheng Chen, Hongsheng Liu, Zidong Wang, and Bin Dong. PDEformer:
Towards foundation model for one-dimensional partial differential equations. In ICLR 2024 Workshop A14DiffernttalEguttion In Science, 2024. URL https: /ttps://openreviee m
Jianwei Zheng, LiweiNo, Ni Xu, Junwei Zhu, XiaoxuLin, and Xiaoqin Zhang. Alias-free mamba neural operator. In The Thirty-eight eight Annual Conference on Neural Information Processing Systems, 2024. URL https://openreview net//forumm (ig.uuuu?dddUUEEEx......
Lianhao Zhou, Hongyi Ling, Keqiang Yan, Kaiji Zhao, Xiaoning Qian, Raymundo Arróyave, Xiaofeng Qian, and Shuiwang Ji. Toward greater autonomy in materials discovery agents: Unifying planning, physics, and scientists. arXiv preprint 2rx: 05666066666 2025.
O.C. Zienkiewiez and R.T. Taylor. The Finile Flement Method: Its Basis and Fundamentals. The Finite Element Method, Butterworth- Hnntervorth-HHiimmaan. 2013. ISBN 9780080951355. URL /ttps://bbodte hwwgg/o...o..........e com//bokssFid
David Zwicker. py-pde: A python package for solving partial differential equations. Journal of Open Source Software, 5(48):2158, 2020.
17

---

Published in Transactions on Machine Learning Research (02/2026)
A Additional Experimental Results
A.1 Detailed results on test error
In Table 3, we present the full results measured in normalized root mean squared error (nRMSE) for all the inference strategies and agentic approaches that we sludy in our evaluation.
A.2 Detailed results on code correctness and debugging success rate
In Tables & 6, we present detailed bug-free rates for each LLM on each PDE. For each LLM and each PDE problem, we generate 32 ii.d. samples (i.e., solvers), and the rates are computed based on these samples.
In Figure 7, we present the average numbers of debugging iterations required for the language models generate executable code free from numerical issues." 1 For the solvers that are executable upon the initial generation, treat their numbers of debugging iterations O.
A.3 Detailed results on scaling test-time compute
In Figure 3, we compare the test-time scaling curves via repeated sampling among different models. Here, we additionally present the test-time scaling curve for each model on each problem individually in Figure 12.
A.4 Additional visualization on convergence rates
In Figures 8, we present visualizations on the convergence rales and solver libraries, in addition to Figures in the main paper.
A.5 Analysis on execution time
Section 5.5 presents analysis on solver convergence rates. To further the analysis, we measure average code runtime as another proxy for computational efficiency. Figure 11 shows average code exccution time across PDEs/LLMs as an aggregated metric. Its left panel shows that execution times vary significantly across models, with Gemini 2.5 Pro generating the most efficient code average. However, LLMss can occasionally introduce redundant operations inefficient looping structures, which leads to slower execution. According the right panel of Figure 11, CNS solvers tend to run slower among the tested PDE families, likely due to its complexity. The finer-grained average code execution time for each LLM on each PDE family can be found in Figure 13.
A.6 Analysis on solver library usage
CodePDE allows LLMs to flexibly select solver libraries (e.g., PyTorch, JAX). We analyze the Python libraries utilized in the generated solvers, with summary statistics shown in Figures 9 & 10. We observe that Gemini 2.5 Pro most frequently employs PyTorch and JAX with GPU acceleration, which may partially account for superior solver efficiency.
A.7 Analysis on LLM token usage
To quantify the inference costs CodePDE, we report output token counts across our PDE benchmarks. Figure details the tokens required for initial generation; Figure 15 accounts for the total tokens including up to four self-debugging rounds; Figure 16 presents the token counts in the refinement stage. As expected, challenging PDEs like the CNS equation exhibit higher computational demands due to increased token length. Refinement costs more tokens than the solver generation stage. While pricing varies across different LLM providers, solving a single PDE using the full CodePDE framework (with repeated sampling and debugging) generally incurs a total cost of less than $5 USD when using APl-based services.
This analysis considers only solvers that successfully execute within debugging iterations. Solvers that fail after iterations are excluded from the average.
18

---

Published in Transactions on Machine Learning Research (02/2026)
Table 3: Normalized RMSE comparisons. Geom. Mean refers the geometric mean of nRMSE values. Gemini 2.0 Thinking short for Gemini 2.0 Flash Thnnking. Neural Network baseline results are from (Takamoto et al., 2022; Li et al., 2022b).
nRMSE (J)
Advection
Burgers
React- Diff
CNS
Darcy
Geom. Mean
Reference Numerical Solvers
1.03 10-3 3.55
10-4
2.29 10-3 1.89 x 10-2
4.80 10-3
2.38 x 10-3
PDE Solving Softwares
Dedalus PyPDE
6.20 x 10
2.78 x 10-3
5.20
10-2 2.01 x 10-2
2.03 x 10-3
5.92 x 10
2.21 x 10-1
2.52 x 10-1
Neural Network & Foundation Model baselines
U-Net
5.00 x 10
2.20 10
6.00 10
3.60 x 10
6.98
10
FNO
7.70 x 10
7.80 x 10
1.40 x 10 3
9.50 x 10
9.80 x 10
9.52 x 10 3
PINN
7.80 x
8.50 x 10
8.00 x 10
8.09 x
PDEformer
4.30 10
4.60
10-3
4.45 10
UPS
2.20 x 10-3
3.73 x 10-2
5.57 x 10-2
4.50 x 10-3
1.20 x 10
ORCA
9.80 x 10-3
1.20 x 10-2
3.00 x 10-3
6,20 x 10-2
1.22 x 10-2
CodePDE: Reasoning + Debugging (best of 32)
Gemini 2.0 Flash
1.14 x 10-3
2.97 10-4
2.19 10-1
4.50 10-2
4.80 10-3
6.92 10-3
Gemini 2.0 Thinking
5.54 x 10-3
3.21 x 10-4
2.19 x 10-1
7.56 x 10-2
4.80 x 10
1.07 X 10-2
Gemini 2.5 Pro
1.01 x 10-3
1.23 x 10
1.49 x 10-1
7.39 x 10-2
4.89 10
5.82 x 10-3
Qwenn-2.5-Max
4.97 x 10
1.35 10
9.57 10
2.36 10
6.76 10
4.00 10
QwQ-Plus
1.03 x 10-3
3.05 x 10-4
2.20 x 10-1
7.43 x 10-2
4.80 10
x 10
Clande-3.5-haiie
3.70 x 10-3
3.23 x 10-4
2.20 x 10-1
2.34 x 10-1
1.00 x 10"
3.61 x 10
Claudo-3.7-somae
1.32 x 10
2.59 10
5.19 10
4.26 x 10
4.80 10
5.15 x 10
DeepSeek-V3
7.37 x 10
6.02 10
2.18 10
2.41 x 10
1.00 10
2.98 X 10
DeepSeek-R1
1.05 x 10
2.76 10
2.20
10
7.46 x 10
4.80 x 10 3
7.44 X 10
GPT-4o mini
5.41 x
9.93 10
3.07
2.36 X 10
1.00 10"
2.08 X 10
GPT-4o
1.55 x 10
3.69 x 10
1.99
10
1.81 10
7.67 10-1
1.74 X 10
GPT-4.1
1.50 x 10
3.63 x 10
1.44
x 10
5.20 x 10-2
4.88 10
7.23 x 10
o3
9.74 x 10
3.69 x 10
1.98
x 10
2.80 x
4.92 x 10
6.28 x 10
ol-mini
6.45 x
3.40 10
1.55
10
2.36 x
8.20 10
3.66 x 10
03-mini
1.05 x 10
2.76
10
2.10
X 10
4.48 x 10
4.83 10
6.67 X 10
04-mini
1.21 x 10-2
3.37 X 10
1.4
10
6.11 X 10
4.92 x 10
1.12 x 10
CodePDE: Reasoning + Debugging + Refinement (best of 12)
Gemini 2.0 Flash
8.24 x 10
1.06 10-4
1.76 -0-2 2
1.72 10
4.78
10-3
2.63 10
Gemini 2.0 Thinking
9.74 x 10-1
2.70 x 10-1
1.74 x 10-2
2.41 x 10-2
4.80 10
3.51 x 10
Gemini 2.5 Pro
9.74 x 10-4
3.15 x 10-4
9.94 x 10-2
2.77 x 10-2
4.92 x 10-3
5.29 x 10-3
Qwenn22.5-Max
9.74 10-4
2.60 10-4
9.13 10-2
1.49 10-2
4.80 10-3
4.40 10-3
QwQ-Plus
9.74 x 10
1.07 10
5.19 x 10-2
1.50 x 10-2
4.80 10-3
x x 10-3
Claude-3.5-haiie
1.01 x 10-3
2.60 x 10-4
5.19
10-2
1.06 x 10-1
1.92 x 10-1
1.22 10
Clande-3.7-somme
9.74 x 10-1
2.60 10
5.19
10-2
2.72 x 10-2
4.83 10-3
4.44 10
DeepSeck-V3
9.74 x 10-4
2.59 10
1.74
10-2
1.56 x 10-2
4.80 10-3
x 10-3
DeepSeck-R1
9.74 x 10
2.59 x 10
1.43 x 10
1.52 X 10-2
4.80 10
4.83 x x 10
GPT-4o mini
1.32 x 10
2.60 x 10
1.76 x 10
4.41 X 10-2
4.80 X 10
4.18 x 10
GPT-4o
9.74 x 10
2.59 10
1.76 10
2.41 10
4.80 10
3.48 x 10
GPT-4.1
1.01 x 10-3
3.15 10
1.44 10
1.53 x 10
4.88 10
5.08 x 10
o3
1.01 x 10-3
3.15 x 10
2.01 x 10-1
1.31 x 10-2
4.88 x 10-3
5.27 x 10-3
ol-mini
9.74 x 9.74 x 10
3.15
10
9.07 10
8.33 x
4.92 10-3
6.48 x 10
03-mini
2.59
10
5.19
10
1.57 x 10
4.80
10
3.97 x 10
04-mini
9.94 x 10
3.15 x 10
2.01 X 10
2.35 X 10
4.92 10
5.92 x 10
19

---

Published in Transactions on Machine Learning Research (02/2026)
Table 4: Normalized RMSE comparisons. Geom. Mean refers to the geometric mean of nRMSE values. Gemini 2.0 Thinking short for Gemini 2.0 Flash Thinking. Neural Network baseline results are from (Takamoto et al., 2022; Li al., 2022b).
Reference Numerical Solvers
1.03 10-3
3.55 x 10-1 2.29 x 10-3 1.89 x 10-2
4.80 x 10-3 2.38 x 10-3
PDE Solving Softwares
Dedalus PyPDE
6.20 10 10-3
2.78 10 5.20 10 2.01 10 10-2
2.03
10
5.92 x 10
2.21
10
2.52 10-1
Neural Network Foundation Model basclines
U-Net
5.00 x 10-2
2.20 x 10-1
6.00 x 10-3
3.60 x 10-1
6.98 x 10-2
FNO
7.70 10-3
7.80 x 10-3
1.40 x 10-3 3
9.50 10-2
9.80 x 10-3
9.52 x 10
PINN
7.80 10-3
8.50 x 10-1
8.00 10
8.09 10
PDEformer
4.30 10-3
4.60 x 10-3
4.45 x 10
UPS
2.20 10-3
3.73 x 10-2
5.57 x 10-2
4.50 x 10-3
1.20 x 10
ORCA
9.80 10-3
1.20 x 10-2
3.00 x 10-3
6.20
10-2
1.22 x 10
FunSearch- PDE
Gemini 2.0 Flash Gemini 2.0 Thinking
1.06 10
2.97 x 10-4
2.06 10-1
4.51 10
2.21 x 10-1
1.45 10
1.17 10
1.00 x 10"
2.19 x 10-1
1.00 x 10
4.80 x 10-3
6.57 x 10 10
Gemini 2.5 Pro
1.05 10
1.13 x 10
3.72 x 10
6.86 10
4.78 x 10
4.29 x 10
Qwen-2.5, Swenn22.5-Mas 2.5-Max
1.17 10-3
3.32 x 10-4
x.47
2.52 10-1
3.01 x 10-1
1.59 10-2
1.05 10
2.75 x 10-4
2.20 10-1
7.31 10-2
4.80 x 10-3
7.41 10
-laude-3.5-haile
4.30 x 10
1.39 x 10-3
x x
2.49 x 10-1
1.00 x 10
3.78 x 10
7--onntt
9.84 10
2.59 x 10
x x
3.79 x 10
4.80 x 10
4.73 x 10
DeepSeek-VS
2.27 10-3
1.30 x 10-3
3.79 x 10-2
1.00 10
1.00 x 10
4.07 X 10
DeepSeck-R1
1.05 10
2.76 x 10
3.77 x 10
7.36 10
4.83 x 10-3
5.22 10
GPT-4o mini
5.41 x 10
9.93 x 10
5.92 x 10-1
2.52 x 10
1.00 x 10
2.40 X 10-1
GPT-4o
1.55 10-3
2.99 x 10-1
1.97 x 10-2
9.54 x 10-1
1.00 x 10"
2.44 X 10-2
GPT-4.1
1.05 10-3
2.22 10
1.68 x 10-1
1.00 10
4.78 x 10-3
1.13 10
o3
1.03 10-3
2.74 10
2.20 10-1
4.40 10-2
4.65 x 10-3
6.61 10
ol-mini
1.05 10
2.76 x 10
1.14 x 10-1
9.36 x 10
6.76 X 10
1.84 X 10
03-mini
1.05 x 10
2.58 x 10-1
2.20 x 10-1
4.57 x 10-2
4.78 x 10-3
6.65 x 10-3
04-mini
9.69 10-1
2.66 x 10-1
2.09 x 10-1
10-2
4.80 X 10-3
6.91 10-3
AIDE-PDE
Gemini 2.0 Flash
6.24 10
2.90 x 10
1.00 10
1.00 10
9.97 x 10
1.78 10
Gemini 2.0 Thinking
6.24 x 10-1
2.90 x 10-4
1.00 x 109
1.00 x 100
1.00 x 100
1.78 x 10
Gemini 2.5 Pro
1.05 x 10
2.76 x 10
1.73 x 10
7.47 x 10
4.78 x 10-3
7.09 x 10
Owen-225-Mas.
1.03 10
3.92 x 10
2.21 10
1.00 10
1.00 x 10
3.89 10
QwQ-PPus
1.05 10-3
2.75 x 10-4
2.21 x 10-1
2.52 x 10-1
4.91 x 10-3
9.53 x 10-3
Claude-3.5-haide
2.95 10
8.56 x 10
1.03 x 10-1
2.55 x 10
1.00 x 10
9.21 x 10 10
Claude-3.7 7.-de-t.
1.03 10
1.05 x 10
5.07 x 10
5,77 x 10
4.78 x 1.00 x 100
4.33 x 10
DeepSeek- V3
1.06 10
1.05 10
2.21 10
1.00 10
3.01 10
DeepSeek-RI R1
1.03 10-3
4.08 x 10-4
9.68 x 10-2
1.24 x 10-1
4.80 x 10
7.53 X 10
GPT-4o mini
5.41 10
1.00 x 10
5.92 x 10
2.52 x 10
1.00 x 10
2.41 X 10
GPT-4o
5.41 x 10
7.29 x 10
5.92 x 10
2.51 10
1.00 x 10
5.67 x 10
GPT-4.1
1.17 10 10-3
3.26 x 10
1.59 x 10-1
4.65 10
4.78 10
6.69 10
o3
1.05 10-3
2.60 x 10-4
5.19 x -
7.53
10-2
3.73 x 10
5.25 x 10
ol-mini
9.27 10
3.61 10
1.87 x 10
9.95 x 10-1
3.61 x 10-1
4.68 x 10
03-mini
1.34 10
3.23 x 10
2.18 x 10 10
1.00 x 10
4.79 x 10
2.14 x 10
o4-mini
1.05 10
2.69
10
2.00
10
7.24
4.67 10
7.18 10
20

---

Published in Transactions on Machine Learning Research (02/2026)
Table 5: Bug-Free Rates in the initial round of code generation.
Advection
Burgers
React-Doe
CNS
Darcy Flow
Average
Gemini 2.0 Flash
87.50%
62.50%
18.75%
0.00%
21.88%
38%
Gemini 2.0 Thinking
93.75%
43.75%
43.75%
6.25%
25.00%
42.30%
Qwen-2.5-Man
56.25%
9.38%
43.75%
12.50%
6.25%
25.63%
QwQ-Plus
46.88%
59.38%
43.75%
28.13%
40.63%
43.75%
Claude-3.5-haiiee
65.63%
28.13%
9.38%
43.75%
31.25%
35.63%
Claude-3. -
90.63%
78.13%
37.50%
3.13%
59.38%
53.75%
DeepSeck-V3.
87.50%
15.63%
21.88%
0.00%
25.00%
30.00%
DeepSeek-RI
93.75%
71.88%
59.38%
34.38%
59.38%
63.75%
GPT-4o mini
84.38%
15.63%
31.25%
25.00%
12.50%
33.75%
GPT-4o
78.13%
9.38%
28.13%
0.00%
0.00%
23.13%
03-mini
96.88%
34.38%
90.63%
3.13%
75.00%
60.00%
Average
80.11%
38.92%
38.92%
14.20%
32.39%
40.91%
Table 6: Bug-Free Rates after up to 4 iterations of self-debugging.
Advection
Burgers
React-diff
CNS
Darcy Flow
Average
Gemini 2.0 Flash
96.88%
93.75%
96.88%
50.00%
62.50%
80.00%
Gemini 2.0 Thinking
100.00%
96.88%
93.75%
34.38%
84.38%
81.88%
Qwen-2.5-Maa
93.75%
68.75%
96.88%
75.00%
34.38%
73.75%
QwQ-Plus
100.00%
96.88%
96.88%
78.13%
100.00%
94.38%
Claude-3.5-hait
100.00%
93.75%
90.63%
93.75%
96.88%
95.00%
Claude-3. 3--somnt
100.00%
93.75%
96.88%
68.75%
100.00%
91.88%
DeepSeck-V. V3
100.00%
87.50%
81.25%
62.50%
68.75%
80.00%
DeepSeek-R1
100.00%
96.88%
100.00%
50.00%
100.00%
89.38%
GPT-4o mini
100.00%
100.00%
90.63%
87.50%
62.50%
88.13%
GPT-4o
100.00%
68.75%
78.13%
40.63%
37.50%
65.00%
03-mini
100.00%
87.50%
100.00%
34.38%
90.63%
82.50%
Average
99.15%
89.49%
92.90%
61.36%
76.14%
83.81%
Number of Iterations
Average Number Debugging Iterations Comparisons across LLMs/PDES
2.05
1.76
1.5
1.35
1.16
1.07
1.0
0.78
0.62
0.89
0.5
0.33
0.37
0.46
0.0
Gemini2 Gemini2.0Flash 0Flash
Gemini2.00 Gemini2.0Thinkinng @Thinking
Gemini2.5Pro
Qwen Max
QWQ Plus
Claude-3. H-Halku
Claude-3.7 3.noon..
Sannet
DeepSeek-V3
DeepSeek-R1
GPT-40
GPT-4.1
03
Advection
Burgers
React-Diff
CNS
Darcy
Figure 7: We report the average number of debug iterations. Left:: Averaged across all five PDE datasets for each model. Right: Averaged across LLMs for each PDE family.
21

---

Published in Transactions on Machine Learning Research (02/2020)
Solver Convergence Comparisons across LLMs/PDEs
100%
80%
27%
23%
7%
15%
11%
31%
15%
11%
9%
22%
8%
Proportion
29%
60%
60%
58%
65%
40%
56%
44%
63%
45%
69%
43%
60%
8.8%
54%
48%
6.%%
61%
20%
18%
17%
24%
9%
24%
20%
14%
23%
22%
22%
Geminiz OFlash Geminiz OThinking Geminiz .pro
Qwen Max QWQ plus Claude-3 Haiku Claude-3 Deepseek V3
DeepSeek GPT-40 4o mini
R1
GPT-40
GPT.4.1
o3
ol-mini
03 mini
oi-mini
Advection
Burgers
React-Diff
CNS
Darcy
Fail
Order
Order
Order 2
Figure 8: The distribution of solver convergence rates for each model and PDE dataset.
Library Usage Comparisons across LLMs/PDES
PyTarch 34.7%
2.7%
95.0%
35.6%
37.6%
31.7%
38.6%
57.7%
85.1%
76.1%
48.4%
4.0%
50.0%
64.1%
46.7%
55.0%
37.8%
JAX- 0.0%
0.0%
5.0%
10.5%
49.8%
6.3%
2.5%
0.0%
9.6%
6.0%
0.0%
0.0%
6.3%
5.8%
6.2%
7.4%
3.5%
50%
SciPy
NumPy 46.6%
46.6%%
84.4%
0.0%
40.0%
10.7%
31.2%
36.9%
24.0%
2.1%
9.6%
38.2%
96.0%
44.6%
29.0%
38.8%
36.6% 16.4%
Geminiz. OFlash
Geminiz. OThinking
Gemini2. 5Pro
Owen Max
Claude-3.5 3.5-ea-.5-Haiku Haiku
QWQ Plus
Claude-3.7 Moadee3..7-onne -Snneet
DeepSeek-VD 3.7-onnet
DeepSeek-R1
GPT-40
GPT-4.1
Advection Burgers
React-Diff
CNS
Darcy
Figure 9: The solver libraries picked by each model or used on each PDE problem.
Solver Library Usage Distribution Comparisons across LLMs/PDEs
100%
80%
47%
40%%
13%
31.1%
24%
10%
0%
39%
37%
10%
Proportion
40%
60%
19%
844
95%
144%
1.0%
50%
31%
22%
13%
96%
7%
4215
65%
58%
50%
64%
20%
35%
36%
38%
38%
0%
Gemini2. Gemini2.0Flash 0Flash Geminiz. GeminiZ.07hinking OThinking
Gemini2.5Pro .Prr
Owen Max
Claude 3.5ude3.5-Haiku 5-Haikuu
QWQ Plus
Caude-3.7 1.ree-3.---a 3.7--onneti Sonnet
DeepSeek-V3
DeepSeek-R1
GPT-40
GPT-4.1
03
Advection
Burgers
React-Diff
CNS
Darcy
Gemini2.
PyToch
JAX
SciPy
NumPy
Figure 10: Distributions of the solver libraries used by selected LLMs.
Average Code Execution Time Comparisons across LLMs/PDES
Average Time (s)
378.41
360
275.97
298.33
260.86
240
221.50
240.78
219.85
17584
144.89
137.73
178.75
14R03
120
101.81
73, 24.34
Gemini2 OFlash
Geminiz.oom Gemini2..Thinkkinn oThinking
5Pro 5pro
Qwen Max
QWQ Plus
Claude-3 Haiku
Claude-3.7-8 Sonnet
SeepSeek.v3
DeepSeek-RI R1
GPT-40
GPT-4.1
Advection
Burgers
React-Diff
CNS
Darcy
Gemini2.55
Figure 11: Code execution time comparisons. Left: Average code execution time across PDEs for each model. Right: Average code execution time across LLMs for each PDE family.
22

---

Published in Transactions on Machine Learning Research (02/2026)
Test-time scaling of Gemini 2.5 Pro React-Diff
Advection
Burgers
CNS
Darcy
Best log nRMSE
-1.6
-2.985
-3.6
-1.11
-1.8
-2.990
-1.12
-2.0
3.88
-0.00
-2.995
-2.2
-1.13
1011502230 Sample size
1015202533 30 Sample size
101520253 Sample
10155202330 Sample size
1015220530
Sample size
Test-time scaling of QwQ-Plus React-Diff
Advection
Burgers
CNS
Darcy
Best log nRMSE
-0.6
-1.5
-2.94
-3.2
-0.455
00.8
-2996
0.650
-2.88
-3.4
0.655
-1.0
-2.0
5 10 25 30 Sample size
1015 20 30
20 10 10155225 30
10 15 20
5 10 15 20
Sample size
Sample size
Sample size
Sample size
Test-time scaling of Claude 3.7
Advection
Burgers
React-Diff
CNS
Darcy
Best log ARMSE
-2.6
-3.4
-1.1
-1.0
-0.8
-2.7
-3.5
-1.0
-1.2
-1.5
-1.3
-2.0
-1.2
11115 20
10 15
10 15 20
10 30
1015 Sample size
Sample size
Sample size
Sample size
Sample size
Test-time scaling of DeepSeek-R1 React-Diff
Advection
Burgers
CNS
Darcy
Best log NRMSE
-2.
-3.2
0.6565
-1.122
0.6570
-1.124
-2.5
-3.4
0.6575
-1.126
0.6580
-3.0
20152025
1015 20 25
10 20 Sample size
10 2015 Sample size
1015 25
Sample size
Sample size
Sample size size
rest-time scaling of o3 React-Diff
Advection
Burgers
CNS
Darcy
Best log nRMSE
-3.010
33.0
-06998
-1.2
-0.5
-1.0
-3.011
-3.2
-1.4
-1.5
0.702
-2.0
-3.4
-3.012
101520255
1.5 25
1015 Sample size
20 Sample size
25 Sample size
Sample size
Sample size
Figure 12: Smoothed test-time scaling curves for each representative models on each PDE family. We increase in the best-of-n sampling to study test-time scaling performance.
23

---

Published in Transactions on Machine Learning Research (02/2026)
Average Code Execution Time for Advection Equation
32.69
Time (s)
21.37
22.97
20
12:55
0.74
2.47
1.13
0.40
0.63
0.39
0
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5-Max
QwQ-Plus Plus
Claude-3. 3l 5-halku
Claude 7-sonnet
DeepSeek-V3
DeepSeek-RI
GPT-40 mini
GPT-40
GPT-4.1
0l-mini
03 03-mini
04-mini
o3
QWQ
Average Code Execution Time for Burgers Equation
221.
200
Time (s)
123:94
93.42
72.61
66.71
$2.65
10.40
5.53
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.55 awen-225.Max M.Max
QwQ Plus
Claude Gaude-3.5-haiky 5-hhalku
Claude-3. --sonnet
DeepSeek-RI DeepSeek-V3
-1
GPT-4o mini
GPT-40
GPT-4.1
ol-mini
03-mini
04-mini
o3
Average Code Execution Time for CNS Equation
1000
Time (s)
696.84
825.31
500
377.29
176.55
61.76
17.30
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5-Ma
QwQ
QwQ-Plus Plus
Claude-3 -laude-3.--halku h-hhha
Claude-3.7 3,dde-377-so --sonnet
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-40
GPT-4.1
0l-mini
03 03-mini
04 04-mini
03
Average Code Execution Time for Darcy Flow Equation
400
366.76
299.62
273.77
304.29
200
60.96
46.95
84.74
118:78
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2 5.Max
QwQ-Plus
Claude 5-haiku
3.-haku
Claude-3.1 7 sonnett
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-40
GPT-4.1
01-mini
03-mini
04-mini
03
Average Code Execution Time for Reaction Diffusion Equation
771.39
Time
500
388.96
354.48
213.91
260.22
217.46
238.75
103.97
113,2
5.43
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5.Ma Max
Claude-3.5 3.ude-3.5-halku 5-halku
DwQ-Plus
Claude-3.7 sonneet
DeepSeek.V3
DeepSeek-Ri
GPT-40 mini
GPT-40
GPT-4.1
01-mini
03-mini
04-mini
o3
Figure 13: The average code execution time for each LLM each PDE family.
24

---

Published in Transactions on Machine Learning Research (02/2026)
Token Count for Solver Generation (w/o debugging) for Advection Equation
Average Tokens
2297.22
2000
1000
230.31
1227.50
71341
942.84
1161.78
1206.78
604.06
1010.09
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5-Max
QwQ-Plus Plus
Claude-3. 3laude-3.5-halku 5-halku
Claude 7 --soneet
DeepSeek-R1 DeepSeek-V3
GPT-40 mini
GPT-40
GPT-4.1
ol-mini
o3-mini
04-mini
03
QWQ
Token Count for Solver Generation (w/o debugging) for Burgers Equation
Average Tokens
2736.12
2000
1c90.56
1242.44
337.97
1069.12
1466.80
963.09
974.28
1415.50
1450.53
1568.71
1286.84
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.55 Dwen-2.5.Max M.x
QwQ Plus
Claude-3.5 daude-3.5-haiku 5-hhaku
Claude-3. --snnett
DeepSeek-V3
DeepSeek-RI R1
GPT-4o mini
GPT-40
GPT-4.1
ol-mini
03-mini
04-mini
03
Token Count for Solver Generation (w/o debugging) for CNS Equation
Average Tokens
4000
2566.88
2000
1572.34
1510-47
1304.97
1785.25
1564.16
1160.53
2100.95
1830.03
2187.53
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5-Mai
QwQ-Plus
Claude-3 3.5-aa3 hahhu
Claude-3.7 M.ude-3,7-sonee 7-sonnet
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-40
GPT-4.1
0l-mini
03-mini
04 04-mini
03
Token Count for Solver Generation (w/o debugging) for Darcy Flow Equation
Average Tokens
4000
3474.11
2000
1043.66
821.20
1444-09
1032.6 M
1844.53
1027.94
942.25
957.22
932.03
1568.06
1231,16
1499.34
1471.66
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5A4 M.xxx
QwQ-Plus Plus
Claude-3.5 hahiu
Claude-3.7 3.7 sonnet
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-4o
GPT-4.1
01-mini
03-mini
04-mini
03
Token Count for Solver Generation (w/o debugging) for Reaction Diffusion Equation
Average Tokens
2895.1
2000
1486.41
1435.00
1671.00
1487.47
930.69
1176.75
887.26
$96.75
729.76
992.31
972.41
1175.25
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5.Max
QWQ-PIn Plus
Claude-3 Moude-3.5-halle h-hakku
3.. sonnet
DeepSeek-V3
DeepSeek-Ri
GPT-40 mini
GPT-40
GPT-4.1
01-mini
03-mini
04-mini
03
Figure 14: Average output token count (without debugging) for each LLM each PDE family.
25

---

Published in Transactions on Machine Learning Research (02/2026)
Token Count for Solver Generation & Debugging for Advection Equation
Average Tokens
2391.82
2252.06
2401.12
2000
1216.88
1263.28
1272,88
1353.09
1396.78
1480.66
1410.69
1131.97
1332.62
1000
792.84
654.53
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5-Max
QwQ QwQ-Plus Plus
Claude-3. 3laude-3.5-halku 5-halku
Claude 7 --sonnet
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-40
GPT-4.1
ol-mini
o3-mini
04-mini
03
Token Count for Solver Generation & Debugging for Burgers Equation
Average Tokens
4498
4799.82
4000
3715.72
3034.97
X14 06
3137.10
2308.13
2325.64
2311-03
2000
1741.72
1751.47
1668.91
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2 Dwen-2.5.Max 5.Max
QwQ Plus
Claude Gaude-3.5-haiku 5-hhaku
Claude-3. 7 --soneet
DeepSeek-R) DeepSeek-V3
GPT-4o mini
GPT-40
GPT-4.1
ol-mini
03-mini
o4-mini
03
Token Count for Solver Generation & Debugging for CNS Equation
Average Tokens
11479.65
10422.25
10000
6712.06
122
7961.91
6525.62
5338.02
7494.88
6099.81
7639.03
3331.34
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5-Mai
QwQ Plus
QwQ-Plus
Claude-3 3.5-ha3. 5-hhakku
Claude-3.7 3aude-3,7-somee 7-sonnet
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-40
GPT-4.1
0l-mini
03-mini
04 04-mini
03
Token Count for Solver Generation & Debugging for Darcy Flow Equation
Average Tokens
471547
5000
3798.31
3000.03
4867.00
2989.97
3392.62
3802.91
4131.53
298 3.56
4583,51
3004.23
1541.34
2015.81
1743.88
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2 Maxxx
QwQ-Plus Plus
Claude 3.ude-3.5-haiki) 5-haiku
Claude-3.7 3.7 sonnet
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-4.1 GPT-40
01-mini
03-mini
04-mini
03
Taken Count for Solver Generation & Debugging for Reaction Diffusion Equation
Average Tokens
4642.51
4000
3658, 64
2373.22
2028.68
2743.72
2804.19
2761.34
2564112
2000
3450.31
1680,78
1082.50
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5.M 2.5Maxx Max
QwQ-PIn Plus
Claude-3 Saudee3.5-hallee h-halku
3.aude-3.7 7 sonnet
DeepSeek-V3
DeepSeek-Ri
GPT-40 mini
GPT-40
GPT-4.1
01-mini
03-mini
04-mini
o3
Figure 15:
Average output token count (with debugging) for each LLM on each PDE family.
26

---

Published in Transactions on Machine Learning Research (02/2026)
Token Count for Solver Refinement & Debugging for Advection Equation
Average Tokens
5993.38
5000
3145.00
2957.33
3084.00
2993.50
2500
2121.4
3720.42
1977.50
2333.83
1894.33
1716.25
1121.50
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
owen-2.5-Mal
QwQ-Plus
QWQ Plus
Claude-3.5 3laude-3.5-halku 5-halku
Claude 7 --sonnet
DeepSeek-R) DeepSeek-V3
GPT-40 mini
GPT-40
GPT-4.1
ol-mini
o3-mini
04-mini
03
Token Count for Solver Refinement & Debugging for Burgers Equation
Average Tokens
10000
2600.25
8940.17
5000
353875
4244.7
3753.00
4694.32
4650.67
2406.08
2616.67
2111.75
2349.33
2036.92
2754.83
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2 Owen-2.5.Max 5.x
QwQ Plus
Claude daude-3.5-haiku 55-halku
Claude-3. --nnnett
DeepSeek-Ri DeepSeek-V3
-R1
GPT-4o mini
GPT-40
GPT-4.1
ol-mini
03-mini
04-mini
03
Token Count for Solver Refinement Debugging for CNS Equation
Average Tokens
20000
10804-83
10000
9330.42
6050.33
7951.75
$420.67
8411.67
6767.83
5885.42
5487.75
$398.83
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5-Mai
QwQ-Plus
Claude-3 3audee3.5-halku 5-haiku
Claude-3.7 7.7 -onnnet
DeepSeek-R1 DeepSeek-V3
GPT-40 mini
GPT-40
GPT-4.1
0l-mini
03-mini
04 04-mini
o3
Token Count for Solver Refinement & Debugging for Darcy Flow Equation
Average Tokens
14693
10000
6393.58
5719.08
7127.87
$450.00
$392.33
2925.75
3337.00
2839.92
3437.00
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5.M Maxx
Claude-3.5 5-haiku
Plus
3.5haku
Claude-3.1 7 sonnett
DeepSeek-V3
DeepSeek-R1
GPT-40 mini
GPT-4o
GPT-4.1
01-mini
03-mini
04-mini
03
Token Count for Solver Refinement & Debugging for Reaction Diffusion Equation
Average Tokens
5585.83
5000
4970.25
4866.17
4206.42
3029.55
2762.58
2791.58
2500
1556.58
1846.50
2159.92
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.5.MM MaxaMax
QWQ QwQ-Plus Plus
Claude-3 Saudd-3.5-halkir h-hakku
3.a 7 sonnet
DeepSeek.V3
DeepSeek-Ri
GPT-40 mini
GPT-40
GPT-4.1
01-mini
03-mini
04-mini
03
Figure 16:
Average output token count in the refinement stage for each LLM each PDE family.
27

---

Published in Transactions on Machine Learning Research (02/2026)
B Implementation Details
B.1 PDE Descriptions
We share below the descriptions of PDEs considered in our work which are used to prompt LLMs.
Advection Equation
Your task is to solve a partial differential equation (PDE) using Python in batch mode.
The PDE the 1D advection equation, given by
(0,,(())))0 E (0,1), L (0,2)
u(0, =,)(r),
I (0,
where 8 is constant representing the advection speed. In our task, we assume the periodic boundary condition.
Given the discretization of uo (z) of shape [batch_size] NJ. where is the number of spatial points, you need to implement solver predict u(t,-) for the specified subsequent. time steps (t 4,,-..: (r). The solution of shape [batch_ size, T+1, N] (with the initial time frame and the subsequent steps). Note that although the required time steps are specified, you should consider using smaller time steps internally to obtain more stable simulation.
In particular, your code should be tailored to the case where 8 = =0.1, i.c., optimizing it particularly for this use case. You will be completing the following code skeleton:
inport numpy
def solver (uo, batch, t_coordinate, beta): Bolves the Advection equation for all tines in t_coordinate.
Args:
ub_batch ndarray:: Initial condition [batch_size, X0, where batch_size the number of different initial conditions, and N is the number af spatial grid points. coordinate idprraarr))):
coordinates shape [T+1].
It begins with t_0-0 and follows the time steps t_i, ... t.T. speed.
beta (flost):
Returne:
solutions (ap ndarray): Shape (batch_aize, T+1, N].
solutions [t, 0, !] containa the initial conditions (u0 batch),
solutions 1, :] containa solutions
t_coordinate [1].
Implement
Advection equation
Binto:
1. Consider using PyTorch or JAX with CPU acceleration
Remenber to handle data types and device placanent appropriately
raturn solutions
28

---

Published in Transactions on Machine Learning Research (02/2026)
Burgers Equation
Your task is to solve a partial differential equation (PDE) using Python in batch mode. The PDE the burgers equation, given by
Ocu(s, t) 02 ((())))) vd f), x x (0, t (0, 0) =00)))
x x (0,1)
where is a constant representing the viscosity. In our task, we assume the periodic boundary condition.
Given the discretization of (z) of shape [batch size, NJ where N is the number of spatial points, you need to implement solver to predict u(-(() for the specified subseqent time steps (1 = 41, ... IT). The solution is of shape [batch size, T+1, N] (with the initial time frame and the subsequent steps). Note that althongh the required time steps are specified, you should consider using smaller time steps internally to oblain more stable simulation.
You will be completing the following code skeleton:
inport numpy
def solver (uo, batch, t_coordinate, mu): "n"Solves the Burgers' equation for all tines in t_coordinate.
Arge:
bathatch (ap adarray): where batch_eize is the number of different initial conditions, and is the number spatial grid points.
condition [batch_aize, X),
"_coordinate Inp.ndarray): Time coordinates of ahape [T+1].
It begins with t_0-0 and follovs the tine steps t.l, a. t_T. nu (float): Viecoeity coefficient.
solutions (ap ndarray): Shape [batch_size, T+1, N].
solutions [:, 0, 3] contains the initial conditions (a0, (u0, batch), solutions[;, contains
at tine t_coordinate(1) [1].
Implement
Burgers' equation
Bints:
1. Conaider uaing PyTorch
acceleration
2. Remember to bandle data typea and device placement appropriately
return solutions
29

---

Published in Transactions on Machine Learning Research (02/2026)
Reaction Diffusion Equation
The PDE diffusion-reaction equation, given by
{ @,u((, x)) v0,n -(0x pm(11 u) 0, E (0, te (0.7) (0, - uo(),
I E (0,
where u and p are coefficients representing diffusion and reaction terms, respectively. In our task, we assume the periodic boundary condition.
Given the discretization of (z)) of shape [batch size, NJ where N the number spatial points, you need to implement solver to predict u(.,t) t) for the specified subsequent time steps (t........ tr). The solution is of shape [batch_ size, T+1, N] (with the initial time frame and the subsequent steps). Note that although the required time steps specified, you should consider using smaller time steps internally to obtaín more stable simulation.
In particular, your code should be tailored to the case where y = 0.5, = 1.0, i.e., optimizing it particularly for this use casc.
You will be completing the following code skeleton provided below:
inport numpy as
def solver (uo_ batch, E_coordiaate, rho):
"""Solves the reaction diffusion equation for all tines in t_coordinate.
Arga:
u0_batch (up ndarray):
condition [batch_size, X),
where batch_size is the number of different initial conditions, and is the number spatial grid points. t_coordinate Inpradyrryy): Time coordinates of shape [T41].
It begins with t_0=0 and follovs the time steps t,l, ... t.T. (float): paraneter rho (floRt): The paraneter rhe in the equation.
the equation.
Returns:
solutione hap.aderray) Shape baze T+1, 8].
solutions(:: 0, contains the initial conditions (u0_batch),
solutions[;, 1, contains the solutions at time t_coordinate [i].
**.
TODO: Implement the solver for 1D reaction diffusion equation
Hints:
1. Consider using PyTorch JAX with GPU acceleration
2. Remember to handle data types and device placement appropriatoly
solutions
30

---

Published in Transactions on Machine Learning Research (02/2026)
Compressible Navier Stokes Equation
Your task is to solve a partial differential equation (PDE) using Python in batch mode. The PDE the 1D compressible Navier-Stokes equations, given by
0,p + 0,(pv)
p(0,v + vdzu) = -0,p + +1jmnu + (C 1/3)0,(Om)) 0 c+ mo² + +02 [( < <e+p+ pa² 2 vo' =0
where is the mass density, v is the velocity, is the gas pressure, p/ - (I' with = 5/3, o' = (S án)ā_v v is the viscous stress tensor, and 17, are the shear and bulk viscosity coefficients, respectively. In our lask, we assume periodic boundary conditions. The spatial domain is Ω [-,
is the internal energy
Given the discretization of the initial velocity, density, pressure, each of shape [batch size, NJ where N is the number of spatial points, you need to implement solver to predict the state variables for the specified subsequent time steps (4 l1,).., The solver outputs velocity, density, pressure, each of shape [batch size, T+1, NJ (with the initial time frame and the subsequent steps). Note that although the required time steps are specified, you should consider using smaller time steps internally to obtain more stable simulation.
In particular, your code should be tailored to the case where fj = <= 0.1, ie., optimizing it particularly for this use case. You will be completing the following code skeleton provided below:
inport numpy aa up
def volver (Vxo, densityo, pressured, t_coardinate, ata, zata): """Solves
compressible Navier-Stoker all tines t_coordinate.
Arge:
Vx0 (np. adarray:: Initial velocity of shape bbatch_size X]
where is the nunber evenly spaced spatial points. density0 ndprrdar)))): Initial density [batch,siza 8].
preasured adarraa))))): Initial pressure [batch_size, X]. t_coordinste (ap. adarray):
coordinates It begins with t_0=0 and follovs the time steps t.l, .... t.T.
shape [T+1].
eta (float): rho (float): The bulk viacosity contricients in the equation.
viscosity coafficients
equation.
Returna:
Vx_pred (np. ndarray); Shape [batch_eize Iti, X).
tinefrane identical The aubsequent timeframes are the solutions at the corresponding tine steps.
deneity_pred (mp. ndarray): Shape [batch_ats size, T+1, N]. pressure_pred (dar Shape [batch_size, Tel, E]. .*.
TODO: Implement the solver for 1D compressible Navier-Stoker equations
Bints:
1. Consider using PyTorch JAX with CPU acceleration
Remenber handle data typea device placenent appropriately return Vx_pred, density_pred, pressure_pred
2.
31

---

Published in Transactions on Machine Learning Research (02/2026)
Darcy Flow
Your task to solve a partial differential equation (PDE) using Python in batch mode. The PDE is the 2D Darcy flow equation, given by:
(a(x)Vu(x)) 1, 102
with the boundary condition:
u(x)= = 0, x E (E, 1)2
where u(x) is the solution function, and a(x)) is a batch of coefficient function. Given the distionization of the coefficient. function a(x) of shape [batch size, N, NJ, where the number of spatial grid points in each direction, you need implement solver predict u(z) for the
specified subsequent time steps. The solution should be of shape [batch_ size, N, N].
You will be completing the following code skeleton provided below:
inport numpy as ap
def solver (a):
"S"Sa the Darcy equation.
Args:
(ap. adarray): Shape obatch_oize, и, N), confficient aquation, where ll denotes the nunber of apatial grid points in each direction.
Returna:
udpruaar)))) Shape [batch_size,
Inploment
Darcy equation
Binta:
1. Consider using PyTorch or JAX with GPU acceleration
2. Alternativaly,
NumPy SciPy for CPU-based computation
3. Remenber to bandle data typea and device placenent appropriately
return
32

---

Published in Transactions on Machine Learning Research (02/2026)
B.2 Instructions to the language model
In this subsection, present our instructions to the language model, including system prompl and instructions for debugging and refinement.
System Prompt
You are an intelligent Al researcher for coding, numerical algorithms, and scientific computing. Your goal is to conduct cutting-odge research in the field of PDE solving by leveraging and creatively improving existing algorithms to maximize performances based feedbacks. Follow the user's requirements carefully and make sure you understand them. Always document your code as comments to explain the reason behind them. Always use Markdown cells to present your code.
Example Debugging Prompt
Thank you for your implementation! When running the code, 1 got the following output and error message:
Code output:
Error message: {error_message)
Can you think step-by-step to identify the root cause of the error and provide solution to fix it? Please provide detailed explanation of the error and the solution you propose. You can refer to the code implementation you provided carlier and analyze the error message to identify the issue.
Your response should be in the following format: [Your rationale for debugging (think step-by-ste by-tep]
*x python
[Your bug-Tre bug-free implenentation]
33

---

Published in Transactions on Machine Learning Research (02/2026)
Example Refinement Prompt.
Your task is to solve a partial differential equation (PDE) using Python in batch mode. cade_deecription
You will be improving the following existing code samples. The code samples are provided below: {code_sanee sandless
Your tasks are:
1. Understand the above code samples. Compare their techniques and performances.
2. Identify the parts that could potentially be improved.
3. Plan on how you can improve the function.
4. Improve the function.
The goal to get very low nRMSE (normalized RMSE) and make the code as fast as possible. You must analyze the implementation and test results of the examples provided in the code template, and think about how you can improve them reduce the nRMSE.
If the RMSE is much higher than 1e-2 or becomes Nan, it is likely that there is a bug in the implementation and If the rumning time is much longer than 600s, you must prioritize making it more efficient.
must debug it or think about completely different approaches.
The convergence rate is the empirical order of convergence with respect to spatial resolution. also a good indicator of the performance of the algorithm which you may consider.
You can implement auxiliary functions add additional arguments the function if needed. You use PyTorch JAX with GPU acceleration. You can consider known techniques and analyze their effectiveness based existing results. You should also consider combining existing techniques even developing new techniques since you are doing cutting-cdge rescarch.
Your generated code will be executed to evaluated. Make sure your 'solver' function runs correctly and efficiently. You must use print statements for to keep track of intermediate results, but do not print too much information. Those output will be useful for future code improvement and/or debugging. Your output must follow the following structure:
1. Summary of the existing implementation along with their performances. Identify what techniques work well, what could be buggy (due to high error Nan), and what could be further improved.
2. Rationale for the new implementation (think step-by-step) based on the summary of the existing
implementation.
3. Your python implementation (modularized with appropriate auxiliary functions):
"e "python
[Your improved implenentation]
In the refinement prompl, Scoddesamplle contains multiple samples, each of which is structured as follows:
Code Sample
Code Sample {id}:
"1 ''python
(code)
:er
Running the above code leads to the following output: (code_stdous
The {code_ soodeestdout> is the running log of the code, which include the nRMSE, rumning time, and empirical convergence rate (e.g., nRMSE: 1.355e-02 | Time: 4.56s represented as texts.
Empirical convergence rate: 0.913),
34

---

Published in Transactions on Machine Learning Research (02/2026)
B.3 Detailed Experiment Setup
Datasets. We focus on 5 represenlative PDE families that are commonly employed testbeds for ML methods, including Advection, Burgers, Reaction- Disfiuiion, Compressible Navier-Stokes (CNS), and Darey Flow. These span a range of spatial dimensions, boundary conditions, and numerical stiffness. The datascts are drawn from PDEBench (Takamoto et al., 2022) and FNO paper (Li et al., 2022b) (MIT license). We randomly sample 100 instances for each family for testing. We also sample another 50 instances as the development set, which is used to provide the execution feedback to LLMs without leaking the true test data.
Models. We benchmark 16 LLMs, including proprietary models such as GPT-4o (OpenAI, 2024a), Claude- 3.5/3.7 (Anthropic, 2024), Gemini 2.0/2.5 (Team et al., 2023) Qwen-2.5-Max (Team, 2024), and open- weights models such as DeepSeek-Ve V3 (DeepSeek-Al et al., 2024). We also include reasoning models such as o3 (OpenAT, 2024b), DeepSeek-R1 (DeepSeek-AT et al., 2025), QwQ-Plus (Team, 2025), and Gemini 2.0 Flash Thinking (Google, 2025). Additionally, we evaluate agentic workflows: FunSearch (Romera-Pattee et al., 2024) and AIDE (Jiang al., 2025).
Baselines. We compare LLM-generated solvers against naive numerical solvers based finite difference methods with the central difference scheme (and forward Euler time integration for time-dependent PDEs), well as high-quality solvers based on established softwares or hand-crafted implementations. The high-quality reference solvers are included in the supplementary materials. A brief summary of the reference solver implementation are presented in Tables 7 & 8.
Table 7: Summary of reference solver types and spatial discretization.
PDE
Reference Solver Type
Spatial Discretization
Advection
Analytical
Lincar Grid Interpolation
Burgers
Pseudo-Spectral Method
Fourier Series (N modes, 2/3 dealiasing) Analytical + Central Finite Difference (2nd Order)
React-Diff
Strang Splitting
CNS
Explicit Finite Difference
Central Finite Difference (2nd Order) + Artificial Viscosity
Darcy
Iterative Linear Solver
Central Finite Difference (5-point stencil)
Table 8: Summary of reference solver time integration.
PDE
Reference Solver Time Integration
Advection
Exact Shift (x Bt)
Burgers
Semi-Implicit (Integrating Factor)
React-Diff
Forward Euler
CNS
Stabilized Forward Euler
Darcy
N/A (PDE is time-independent; Conjugate Gradient for the linear solver)
We also include standard neural network solvers and foundation models as baselines, including FNO (Li et al., 2022b), a spectral neural operator; U-Net (Ronneberger et al., 2015), a convolutional architecture for solution rogression; PINNs (Raissi et al., 2019), models that encode PDE constraints into the loss; ORCA (Shen et al., 2023), a cross-modal fine-tuning approach which exlends general lext/ /vision foundation models to specialized domains; PDEformer (Ye et al., 2024), a Traph-Tranffommc et al., 2024b), an UNN--Transfomee--bae PDE foundation model.
PDE foundation model; and UPS (Shen
Evaluation Metrics. We use normalized root mean squared error (nRMSE) as the primary measure of solution quality. We also report (1) dcbug success rate: the proportion of generations successfully fixed after error-driven feedback; (2) convergence rate: empirical convergence rale with increased spatial resolution; and (3) execution time: runtime of generated solvers. All solvers are evaluated on a NVIDIA GeForce RTX 2080 Ti GPU with 11GB memory. We enforce runtime limits to manage computational resources: a timcout of 1200 seconds is applied to the Advection, Burgers, and Darcy Flow Equation solvers, while limit of 2400 Equation solvers.
seconds is set for the more complex CNS and
35

---

Published in Transactions on Machine Learning Research (02/2026)
B.4 Method configurations
Repeated sampling and debugging. For OpenAT's series models, set the decoding temperature to 1.0 due to the API constraints that do not allow using other temperature settings. For the other models, We set the temperature to 0.7. We sample 32 solutions from each model. If the solver runs into an exccution error or Inf /NaN values, we let the LLM debug using the prompt in Appendix B.2.
Refinement.
the refinement stage, we select the best programs generated the sampling & debugging stage for each PDE and use them "seed" programs for refinement. Note that we use the same set of "seed" directly compare the LLMs' ability refining given set
programs for each LLM in this stage that of strong programs. We consider refinement based on 3/4/5 seed implementations, and generate 4 samples in each setting. We then report the best solver performance among the 12 samples.
FunSearch. FunSearch features algorithm discovery based on a program database. The program database consists of a few "islands" of programs. In our experiment, we set the number of islands to and the island reset period to 3600s. We use 2 existing programs the prompt and instruct the LLM to generate an improved version. We run the FunSearch process for 32 iterations. Specifically, in the first 12 iterations, LLMs generate solvers from scratch to warm start the program database; in the subsequent 20 iterations, LLMs improve the existing solvers following an evolutionary algorithm. In each iteration, the language model decoding temperature is set to 0.7.
AIDE. AIDE uses tree search to explore in the space of code. We run AIDE for 96 steps. In the search process, the max debug depth, debug probability, and number drafts are set to 5, 0.9, and 24, respectively. Following Jiang et al. (2025), the language model decoding temperature set to 0.5 for code generation.
36

---

Published in Transactions on Machine Learning Research (02/2026)
C Case Study on Burgers' Equation
We present detailed analysis on the numerical schemes and time stepping strategies used in the generated solvers for Burgers' Equation. We investigate the solver choices from a few axes: spatial discrelization methods, numerical schemes for the convective flux term, the time integration methods, and the time stepping strategies. We prompt OpenAl's o3 and o4-mini models to classify solver choices based on the code. We visualize the distribution of different methods in Figures 6, 17, & 18, respectively.
Figure 5 shows the distribution of spatial discretization methods used by different. LLMs for PDE solvers. Finite Difference dominates across most models, with some variation. QWQ Plus and DeepSeek-Vi V3 incorporate significant Spectral method usage, while GPT-4.1 and o3 employ Finite Volume methods more frequently than other models.
Figure displays time integration method preferences across LLMs. Explicit Euler is the dominant approach for most models. DeepSeek-RI uniquely prefers IMEX methods, while Gemini 2.5 Pro, Qwen-2.5-Max, and o3 utilize higher-order Kutga methods more frequently.
Figure 17 illustrates how different LLMs implement convective flux schemes. Central differencing is the predominant method for most models, while Slaude-3.--oommee uniquely favors Upwind schemes. GPT-4.1 shows the most diverse approach, and o3 heavily utilizes Lax-Friedrichs local method.
Figure 18 presents the time stepping strategics used by LLMs. Typically, strong LLMs like DeepSeck- R1 and tend to choose adaptive time steps, balancing numerical stability and efficiency.
Convective Flux Scheme Distribution
100%
5%
6%
22N
28%
0%
10%
80%
Proportion
6.2
60%
90%
40%
69%
90%
50%
81%
90%
41%
54%
20%
3.1%
23%
17%
16%
10%
Gemini 2.0 Flash
Gemini 2.0 Thinking
0%
Qwen-2.5 Max
QwQ OwQ-Plus Plus
Claude Claude-3.5-haikt h-haiiuu
Claude-3.7 soude-3.7-somes 7-sonnet
DeepSeek-V3
DeepSeek-R1
GPT-4o mini
GPT-40
GPT-4.1
03
ol-mini
03-mini
04-mini
Gemini 2.5
Upwind
Central
Lax-Friedrichs (global)
Lax-Friedrichs (local)
Godunov
Others
Figure 17: Distributions of the numerical schemes for the convective flux term employed by each LLM. "LF" in the legend refers to Lax-Friedrichs schemes, with global and local variants.
Time Step Strategy Distribution
100%
2.8
27%
80%
Proportion
3.2%%
69%
6.%
60%
84%
91.%
22%
8.4%
9-4%
100%
94%
40%
84%
73%
20%
9%
3.4%
0%
6%
Gemini 2.0 Flash
Gemini 2.0 Thinking
Gemini 2.5 Pro
Qwen-2.544 M-x
QwQ-Plus ha5-u
Claude-3.5 Claude-3.5-halku
Claude-3.1 3.7 7-sonnet
DeepSeek-V3
DeepSeek-R1
GPT-4o mini
GPT-4o
GPT-4.1 1
o3
el-mini
03-mini
04-mini
Adaptive (CFL Condition based)
Adaptive (error control based)
Fixed
Others
Figure Distributions of the time stepping strategies employed by each LLM.
37

---

Published in Transactions on Machine Learning Research (02/2026)
D Case Study on Reaction-Difusion Equation
As discussed in Section 5.7, LLMs consistently miss the key trick that involves directly using the analytical solution of the reaction term in the solver. In this section, compare the best two LLM generated solver with the hand-crafted reference solver. We first present concise comparison in Table 9, and then present. the full implementation.
Table 9: Comparisons of Reaction-Diffusion Equation Solvers. Explicit Exponential Time Differencing.
"IMEX-ETD" is short for Implicit-
Reference Solver
LLM Solver 1
LLM Solver
Library Time Integration
JAX Strang splitting Analytical
PyTorch
PyTorch
Sequential splitting
IMEX-ETD
Reaction Term
Discretized
Discretized
Diffusion Term
Central difference
Crank-Nicolson
Spectral
Time Stepping
Adaptive (CFL)
Fixed
Fixed
Reference solver.
inport nampy as up import jax inport jaz. numpy as jup
def solver batch t_coordinate, nu, rho):
an Solooees the reaction diffusion equation all times in c.cordiinate.
Arge:
_0_ batch (ap ndarray): Initial condition [batch_size, }), batch_size
number different initial conditions, and R is the nunber of apatial grid pointa. coordiante (np.ndarray, Time coordinates of shape [7+1].
t_0=0 and fallous the time ataps t_i, ... nu (float): The paraneter nu in the equation. (float): paraneter
bagins
equation.
solutions (ap adarray): Shape [batch_size T+1, X].
solutions [:, 0, :] contains the initial conditions (u0, (u0_batch) batch),
solutions [;, contains solutions time t_coordinate [1]
available
if jax. devicess ''gpa)):
print accaleration enabled")
elae:
print GPU not available, using GPU*)
Convert inpute to JAX arrays
u_batch jap array (a0_batch, pappe=jjpp.fllattt t_ caordinate arpay (traaytt__coreee (tooo(dinatt))))
batch_size, F u_batch.shap can(t_ccordinata) -
Spatial discretization 1.0 Donain [0,
Calculate internal tine step for stability (CFL condition for diffusion) dt_internal 0.25 dz++2
print(f (f"Sinnatinn
paraneters: nu-(bu), rho-(rbo)")
print Gfid size: (-(X), Time points: T=(T+1), Batch size: (batch_size)") print (f"xx((xx:.66 dt, rs_internais(de int.daaliitrral.e...........
Initialize solutions array solutions up. ((bacch_ba T W), fleat5)). solutions o, :] (r (u_batch)
Pieceuise Exact Solution for reaction ters vectorized across batch jjax.
38

---

Published in Transactions on Machine Learning Research (02/2026)
def reaction_step (u_batch,
8olvea the reaction part u.t rho (1 analytically for the entire batch in parallel
Avoid division apailon = 1e-10 return (1.0 oxp dt) u_batch) (u_batch opeilon))
with small epsilon
Diffusion step with periodic boundaries vectorized across batch Qjax jit def dtffurior_ dt):
...
Solvea the diffusion part u_t nu u,xx using central differencing periodic boundary
in parallel
Implement periodic boundary conditions using roll slong the spatial dimension u_left jup. rep. (u_batch, axis=1)
u_right (op.rooll -1, axis=1)
Central difference for Laplacian (vectorized across batch)
laplacian (u_right 2eu_batch u_left) / (dx*+2)
return batch dt nu laplacian
cot-conplled full time step with Strang splitting vectorized acroas batch ejax.jit
def time_ (tt__sttppu_ dt): sun
Perforne tine using Strang splitting for the entire batch: 1. Half step reaction
Full step diffusion
a. Half atep of reaction
...
First half reaction
b_batch reaction_ stec:itnnttp dt/2)
Full diffusion
b_batch diffusion_stop (t_ babchhh dt)
Second saattion.toouuu_bttec dt/2)
b_batch
u_batch
Procass
entire batch simultaneously
for 1 in range (1, T 1):
t_start t_coordinste [i-1] t_end = c_cooddinate()) current_t t_atart
Evolve to next sare point using snall internal steps atep_
while current_t t_end:
sure dt (in (dd__internal, t_end current_t)
step beyond tLend
batch time atun_attp(u_batt)) dt) current_t += dt stap_count +
a Print progrese periodically max(1, T//10) step (13/(7} completed (internal steps: count]_cont]))
Save solution for all batches (u_batch)
aolutions [i, 1, :] battc))
print p ("iimllation completed")
return
LLM-generated solver 1. Generated by GPT-4o.
import numpy as up inport torch
def crank_nicolson_ntte naw.iixxatttxxx (3I, dx, dt. nu):
"""Construcss the and B natrices for the Crank-Nicolses nethod, dx++2)
Identity natrix
torch. (N, davicos/cuda")
39

---

Published in Transactions on Machine Learning Research (02/2026)
Off oft-dnngonal shifts (periodic boundary) ofr_diag -I forch.roll(1. shifts=1, dins=1) off_diag +* torch.roll(I, sbifts--1, dims-1)
and B natrices == (1 2+r)
ofr_diag
off_diag
return A, B
def www.rreettonteennceee. rho, dt): """Applies the reaction tern using explicit Buler. return dt rho
solver (uO batch, c,coodiiatee nu, rho):
"""Solves the 10 reaction diffusion equation for all tines in t_coordinate, Extract the dinensione batch_eize, F c0_baach.ehaa. Lea(t_coordsmatt)
Convert to torch tensors for GPU operations
10_batch tensor
doype=torch.float222 device='cude
Spatial step dx 1.00
step stability dt dt_internal - 0.1 dx++2 / nu auz_internal_ster int (et.ccc(e...ooriiaatteeeee [1] +__coordinate (0)) st_interral))
dt_internal (t_coordinate(11 [1] t_coordinate [o]) aun_internel_stepe
Precompute Crank-Nicolson natrices A, a crank_nicaaiiccccc srank_nicolaon, natrix (M, dr, dt at_intarnal, nu)
i_i torch.linslg inv (A). to('eua)))
array Initialize weroe ((batch_size, T+1, H), device='cuda") solutions solutions[:, 0, :] u0_batch
Time-stepping loap
range T+1): - solutions [:. t-1, c).clone()
range (nang (nt__nintrraalatteeeeooooo ssteps):
apply_reactioon
rho, dt_internsl) Apply reaction
= norchl (B, u.T).T # Now a has shape (batch_ B)
torch. (A_iav, u.T).T
* Apply diffusion tern, msintaining correct shape
solutione
return solutions. () numpy
LLM-generated solver 2. Generated by Clando-3.7-Someet.
Inport numpy
inport torch
typing inport Tuple
solver (u0_ ""oooves the 1D reaction diffusion equation for all times in t_coordinate.
t_coordinate, , au=0.6, r)o=1.0):
Arge:
uo batch (np- Initial condition baatch_iie,, R], batch_eize
initial conditions, and F is the sumber of spatial grid points. coordinate ndarray):
coordinates of shape [7+1].
It begins with t_a=0 and follous the time steps t.l, t.T. nu (float): The paraneter nu the squation. Default 0.5. rho (float): The paraneter
the equation. Default is 1.0.
Returns:
solutions (ap adarray): Shape (batch_size, T+1, E]. solutions [;, 0, contains initial
(u0
solutional:, 1, contains the solutions at time t_coordinate [1].
Check if CPU available davice torch.devica ('cuda' if torch. print (f
in a
else 'cpu?)
(device)")
input PyTorch
40

---

Published in Transactions on Machine Learning Research (02/2026)
batch_ ub shape
u0 tensor torch.tansor
tansor (u0_ batch, dtypemtorch. devicendevice)
Set the spatiel grid and wavenunbers for epectral nethod dx 1.0 /
2.0 mp-pi torch.
fft. ftt. fXtfrag d=dx)
(.to@daviiceddvvice, dtypestorch. faaattt
Linear operator -au k"2 (diffusion tern in Fourier apace) 22
Uae nuch smaller tine step stability dt it_ientrnal a 1e-6 # Very conservative tine step
print(f In time step: (dt_intornal)")
Initialize solution array solutions torch. torch.zoros(
(batch_aize, cotft_oorddaate)) N), dtypa=torch. device-device
solutionsl:, :] ub_tensor Set initial condition
Current solution u. a_currntt = 10_ tensor
current_tine
stepping loop
for in range (1, len(t_coordinate)):
target_tine s_coordinate
(f tine starget_tina:..e
Integrate
target
while current_tine time target_tine:
Adjust
step exactly target needed
cin(dt_internal, target_tine current_time)
Take tine step a_current Siaggrat__topuu_uurree L, dt, rho, device)
simpler, stable nethod
# Check for Ball Inf values
torch. isnan (u, current)
any anarre
(u__enreett....e .atrrent).anyy. ():
(f" Varning: Nall Inf detected tine (current_tine))) Try to recover by clanping values
torch.nan_to_nnu
u_current,
mam=0.5,
po,
neginf=0.0
Clamp to phynically neaningful range u_current torch. clamp ca_crree,,,,, 0.0, 1.0)
current_ting
Store solution solutions [;,, :]
targat
numpy array
retura solutions, cpu () nunpy ()
def intagrate_ategy
u: Forch.Tenaor L: torch. Tensor, dt: float,
fleat, device: torch, device torch.Tenae Tensor:
Take one tine step using stabilized semi iipliepp nethod.
Args:
(torch. Teneor): Current solution.
L. (torch. Tensor): Linear operator in Fourier space. (float): step.
rho (rloat): Reaction confficient. device (torch.device) Device to perfora conputations on.
Returns:
step.
41

---

Published in Transactions on Machine Learning Research (02/2026)
Transforn to Fourier apace
u_hat torch. ffth. dim=1)
Apply 2/3 rule for dealiaring (cet highest 1/3 frequencios to zero) u.shape cutoff = ///
torch. ones (N, dtypertorch. bool, devicandevica) [aat[aatoff:K-cctoo = False
Evaluste nonlinear tern
tine
N_a nonlinear teranrtterenn rho)
N_u_hat torch. dim=i)
Apply dealiasing
N_u_ N_u_bat[:, mask]
Compute the exponentisl tern for the diffueion part Use safer expenentisl calculation with clipping exp_tern torch. torch.exp(L dt)
Prevent underflov by setting very small values exp_ torch. chere(
torch. (exp_etoorral
1e-15,
torch. zeros_ zeros_like (exp_teer),, exp_tern
diffusion implicitly, reaction explicitly u_nev_bat exp_term u_hat + texp_tern 1.0) / L + M_u_hat
Senii Seli-ieplicit update:
node separately avoid division by zero zero_mode (torch.abe(L) < 1e-14)
Handle
if torch. any (zero_node): u_nov_bat [:,, zero_node] u_hat(;, zaro_made] át N_u_hat zero_noda.]
space u_aex torch. torch. ffffffffffff..... (u_nev_bat, dim=1).real
Transforn
Ensure solution stays within physically zeaningful bounds u_nex torch. (u.ve 0.0, 1.0)
return u_nes
def nonlinear_ tern (u: torch. Tensor, rho: float) -> torch, Tensor:
ane
Compute the nonlinear reaction tern rho (1 u))
Arga:
(torch, Tensor): Current solution.
(flost):
coefficient.
Returns:
torch. Tensor: Nonlinear
evaluation.
vithin [o, 13 prevent nuserical issues u_eafe torch, torch. 0.0, 1.0) u_safe)
retura u_safe
To investigate potential mitigation strategies for this failure mode, conduct preliminary experiment utilizing human-in-the-loog guidance. Specifically, we simulate the provision of domain-expert feedback by incorporating the following instruction into the prompt:
Solver Hint
You should consider operator splitting techniques to solve this PDE efficiently. Also note that when you split the operators, the reaction term has analytical solution and you can use 1.0 (1.0 e-rho * dt) (1.0 - u_batch) / (u_batch + epsilon)) for the reaction step.
We evaluate this approach using 04-mini and present the comparative performance metrics in Table 10. With numan-in-the-loop guidance, the LLM-ggnerated solver achieves an nRMSE of 1.69 x 10-3, representing significant improvement accuracy over the vanilla sctting (nRMSE of 1.44 x 10-1). Notably, the guided model outperforms even the reference solver (nRMSE 2.29 x 10-3) while being slightly faster. These promising results demonstrate the effectiveness and promise of human-in-the-loop corrections. More
42

---

Published in Transactions on Machine Learning Research (02/2026)
systematic study the human-in-the-loop PDE solver code generation, well as other potential approaches like targeted model refinement, can be interesting directions of future research.
Table 10: Enhancing Sraction-Diffusion Equation Solvers via Human Hints. For CodePDE with and without hints, the reported uRMSE and running time correspond to the solver selected via best-of-32 sampling.
Method
nRMSE
Running time
Reference solver
2.29 x 10
1329.899
CodePDE w/o hint (o4-mini)
1.44 x 10
191.61s
CodePDE w/ hint (o4-mini)
1.69 x 10-3
1221.578
43

---

Published in Transactions on Machine Learning Research (02/2026)
E Additional Discussions
E.1
When should we trust LLM generated solvers?
While presenting promising results, we note that LLM-generated solvers are not inherently trustworthy in their raw stale. They can be prone to numerical instability, logical errors in complex physics, and suboptimal scheme selection when faced with nuanced mathematical structures. We identify three primary failure modes:
Inabiiiyy Exploit Mathematical Nuances. One of the most distinct failure modes is the inability recognize when parts of PDE possess analytical solutions that can stabilize the numerical method. This phenomenon is detailed in Appendix D)
Ligh-Diffeeulty Physics. As the physical complexity of the PDE increases, LLM performance notably degrades. For example, in single-shot generation (without self self-debugging), only 16.6% of generated solvers for the Compressible Navier-Stokes (CNS) equation are bug-free free.
Complenity Overload
The "Reliability us. Sophistication" Trap. "working" solver is not always a "trustworthy" high- fidelity solver. In Scction 5.5, we highlight a trado-off where some models appear reliable only because they conservatively default simplistic numerical schemes.
However, reliability can be substantially improved through structured inference framework like CodePDE. By incorporating automated debugging and iterative refinement, numerical instability can be significantly mitigated; if a solver encounters NaN or Inf values, CodePDE is designed to autonomously repair the code. Furthermore, the convergence test serves as critical filter: if a solver fails to converge with grid refinement, is not reliable in practice, regardless of its nRMSE performance on isolated test cases. We also observe that deploying frontier models, rather than cost-efficient alternatives, essential for minimizing these failures.
Finally, for high-stakes applications, human verification serves the ultimate safeguard. We note that this stands in sharp contrast to neural network solvers in prior works, which are typically opaque black boxes. Because CodePDE generates human-readable code, it uniquely enables expert review to function as the final layer of quality assurance.
44

---

Published in Transactions on Machine Learning Research (02/2026)
E.2 Fundamental Differences Between Neural PDE Solvers and LLM-Generated Solvers
In our evaluation, we use neural PDE solvers (e.g,, FNO, PINNs, foundation models) baseline for comparison. We note that while both aim to solve differential equations, they differ significantly in their mechanism, data requirements, and interpretability. In this subsection, we discuss the main distinctions and clarify the fairness and scope of our comparisons:
Mcchanism and Output. Neural Solvers typically act as direct function approximators (Raissi al., 2019; Lu al., 2021; Li et al., 2022b). Physics-Informee Neural Networks (Raissi et al., 2019) are trained to approximate solution functions by minimizing physics-ifformed loss. Operator learning such FNO and DeepONet (Lu et al., 2021; Li et al., 2022b) are trained to map imput conditions (initial conditions, parameters) directly to a solution tensor. The "solver" is the neural network itself, and the process is often a "black box" operation where the internal reasoning is opaque. CodePDE does not predict the solution directly. Instead, it frames PDE solving as a code generation task. The output is an executable Python script implementing established numerical methods. The actual "solving" is performed by the generated code, leveraging the rigorous mathematical guarantees of classical nunerical analysis.
Data Requirements and Generalization. Neural Solvers are typically dala-intensive. Task-specific models like FNO require large training datasets of pre-computed solutions to learn the operator. Foundation models (e.g., UPS, PDEformer) rely on expensive offline pre-training and may require additional finetuning when applied in the downstream applications (Shen et al., 2024b; Ye et al., 2024). Their generalization is generally limited to the distribution of the training data (in-distribution generalization). CodePDE operates in zero-shot manner regarding the specific PDE instance. It relies on the LLMs' pre-trained knowledge of code and mathematics, requiring no specific training data for the PDE at hand. It generalizes via reasoning, allowing to construct solvers for novel equations or boundary conditions as long as they can be described in natural language and solved via algorithms.
Interpretability and Debuggability. Neural Solvers often lack transparency; understanding why a particular solution was produced is difficult. CodePDE generates human-readable code. This provides full transparency into the method used. Furthermore, this interprotability enables automated debugging.
Given these fundamental differences, we view neural network methods and CodePDE two distinct paradigms in machine learning for scientific computing. Our evaluation demonstrates that, is capable of achieving competitive or even superior accuracy compared to existing neural network methods. We further emphasize that this comparison is fair if not advantageous to the neural baselines since neural networks typically benefit. from large-scale data requirements and lask-specific training, whereas CodePDE operates without such prerequisites. That being said, we acknowledge that CodePDE requires access to frontier large language models. Ultimately, practitioners should evaluate these tradc-offs and select the paradigm that best aligus with the specific characteristics of the target applications.
paradigm, CodePDE
45

---

Published in Transactions on Machine Learning Research (02/2026)
F
Broader Impact
Our work highlights the potential for LLMs to democratize access to scientific computing tools. As these models improve, they could reduce the expertise barrier for implementing numerical PDE solvers, allowing domain scientists without extensive programming or numerical methods backgrounds to leverage computational approaches their research.
However, this democratization necessitates responsible deployment. Incorrect or numerically unstable implementations might lead to erroneous scientific conclusions if users lack the expertise to critically evaluate the generated solutions. Furthermore, the automated execution of model-generatedd scripts requires robust sandboxing to ensure code-execution safety. Future work should focus not only on improving solution quality but also on enhancing explainability and providing appropriate disclaimers about solution reliability.
As the capabilities of LLMs scientific computing continue to advance, we anticipate a gradual shift in how compulational science is conducted from writing code from scratch lo specifying problems in natural language and refining LLM-generated implementations. This transition also invites consideration of computational sustainability, balancing the cost of largo-scale model inference against the efficiency gains of the generated solvers. CodePDE provides a loundation for measuring progress along this trajeclory and identifying the areas where further research is most needed to realize the full potential of LLMs as scientific computing agents.
46


Reviews:
[{"review_id": "tahlx1OFOk", "reviewer": "Reviewer_vcD2", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper frames PDE solving as a code-generation problem and introduces CodePDE, an inference-time framework that combines (1) chain-of-thought–guided solver generation, (2) iterative automated debugging, (3) feedback-driven solver refinement, and (4) test-time scaling (best-of-n sampling) to produce executable numerical solvers from large language models. The authors evaluate 16 LLMs on five PDE families (Advection, Burgers, Reaction–Diffusion, Compressible Navier–Stokes, Darcy Flow) using metrics including normalized RMSE (nRMSE), debug success rate, empirical convergence behavior, and runtime. Key empirical findings: LLM-generated solvers with debugging + refinement outperform handcrafted reference solvers on 4 of 5 tasks after refinement; automated debugging raises average bug-free rates from 41% to 84%; test-time scaling yields steady gains up to moderate budgets (notably between n=4 and n=16); Reaction–Diffusion remains a consistent failure mode, and CNS shows high convergence-test failure rates.\n\n**Key strengths**\n\n- Clear framing of PDE solving as code generation and a practical, modular pipeline that unites generation, execution feedback, and refinement.\n\n- Thorough empirical study across multiple LLMs and representative PDE families with several useful metrics (nRMSE, convergence tests, debug-success).\n\n- Ablation evidence that debugging, refinement, and test-time scaling are individually important \n\n**Key Weaknesses**\n\n- The Reaction–Diffusion failure mode is identified, but the diagnosis and attempted mitigations are relatively narrow.\n\n- Limited discussion of failure modes and when LLM-generated solvers should not be trusted\n\n- Comparison with neural solvers and foundation models is somewhat unfair, as these methods are designed for different use cases (e.g., learning from data, zero-shot generalization)\n\n- No analysis of the exact costs(LLM+GPU) involved while comparing different methods", "claims_and_evidence_yn": "Yes", "claims_explanation": "The paper provides substantial empirical evidence to support its main claims", "audience_yn": "Yes", "audience_explanation": "This paper addresses an important intersection of machine learning and scientific computing that is highly relevant to the TMLR audience", "requested_changes": "The additions below would further strengthen the work: \n\n- Discuss when LLM-generated solvers should not be trusted and their failure cases\n\n- Clarify the fairness and scope of comparisons with neural PDE solvers and foundation models (Add an explicit section in the Appendix/Discussions to discuss the fundamental differences between PDE solvers and LLMs)\n\n- Provide a transparent accounting of compute and monetary costs for each method (Exact LLM inference costs, GPU hours etc)", "broader_impact_concerns": "None", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "DGlAcfTkFe", "reviewer": "Reviewer_nLdX", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper introduces CodePDE, a new inference-time framework that uses large language models (LLMs) to automatically generate, debug, and refine numerical solvers for partial differential equations (PDEs).\n\nInstead of training specialized neural PDE solvers, CodePDE treats PDE solving as a code-generation problem: given a natural-language description of a PDE and its boundary/initial conditions, an LLM writes executable solver code (e.g., finite-difference, spectral, or IMEX methods).\n\nThe framework adds three key inference-time components that make this process practical and reliable:\n\nAutomated Debugging — runs generated code, captures runtime errors, and prompts the LLM to fix bugs iteratively.\n\nSelf-Refinement — uses numerical feedback (nRMSE) to improve solver accuracy through feedback-driven rewriting.\n\nTest-Time Scaling — employs best-of-n sampling to enhance solver quality by leveraging additional inference compute.\n\nAcross five canonical PDE families (Advection, Burgers, Reaction–Diffusion, Compressible Navier–Stokes, and Darcy Flow), and 16 different LLMs, CodePDE produces high-quality, convergent solvers — often matching or exceeding classical hand-crafted methods and outperforming existing neural solvers that require expensive offline training.\n\nThe work also presents a systematic analysis of LLM behaviors in scientific programming, uncovering:\n\na trade-off between solver reliability (safe, low-order schemes) and sophistication (high-order or diverse methods);\n\ndistinct model capabilities for generation versus refinement; and\n\nclear evidence that debugging + refinement + scaling are the essential ingredients for any scientific “agentic” workflow.\n\nFinally, because the generated solvers are human-readable, CodePDE enables transparent error diagnosis and interpretability — offering a practical bridge between symbolic scientific reasoning and automated code synthesis.", "claims_and_evidence_yn": "Yes", "claims_explanation": "🧩 1. Claim:\n\nLLMs can generate high-quality PDE solvers that rival or outperform traditional and neural solvers.\n\nEvidence:\n\nQuantitative comparisons in Table 1 show lower normalized RMSE (nRMSE) for LLM-generated solvers (with refinement) on 4 of 5 PDE families.\n\nExample: Burgers Equation — Gemini 2.0 Flash achieves \n1.06\n×\n10\n−\n4\n1.06×10\n−4\n vs \n3.55\n×\n10\n−\n4\n3.55×10\n−4\n for the reference solver.\n\nNeural baselines (PINN, FNO, UPS foundation model) perform worse despite training.\n\nAssessment: ✅ Convincing\n\nThe use of multiple PDEs and baselines supports the claim well.\n\nThe accuracy metric (nRMSE) and convergence tests are appropriate.\n\nMinor caveat: reference solver tuning is not deeply discussed, but the performance margin is large enough to be credible.\n\n🧩 2. Claim:\n\nStructured inference — combining debugging, refinement, and test-time scaling — is essential for reliable solver generation.\n\nEvidence:\n\nSection 5.2 → 5.4 + 5.6 (Ablation Study) directly quantify this.\n\nDebugging boosts executable-code success from 41 % → 84 %.\n\nRefinement converts underperforming solvers into surpassing references (e.g., CNS 2.41 × 10⁻² → < 1.89 × 10⁻²).\n\nBest-of-n scaling improves accuracy up to n ≈ 16.\n\nRemoving these components (naïve prompting) yields unusable nRMSE ≈ 1.5 × 10⁻¹.\n\nAssessment: ✅✅ Very convincing and quantitatively backed.\n\nClear causal link between framework components and outcome.\n\nMultiple metrics (nRMSE, convergence, success rate) triangulate the effect.\n\n🧩 3. Claim:\n\nThere exists a trade-off between solver reliability and sophistication (simple vs. high-order schemes).\n\nEvidence:\n\nSection 5.5 Convergence Analysis:\n\nGPT-4o mini → low failure (8.6 %) but first-order methods.\n\nGemini 2.5 Pro and o3 → higher-order methods but more failures.\n\nThe authors visualize convergence-order distributions across LLMs.\n\nAssessment: ✅ Supported and well interpreted.\n\nQuantitative data plus qualitative explanation (safety vs. exploration) match the claim.\n\n🧩 4. Claim:\n\nFrontier models exhibit solver diversity and exploratory capability, important for scientific discovery.\n\nEvidence:\n\nSection 5.7: model-wise analysis of generated numerical schemes.\n\nMost models choose finite-difference + Euler;\n\nDeepSeek-R1 → spectral + IMEX;\n\no3 → broadest scheme diversity.\n\nVisualizations (Figures 5–6) substantiate this behavioral difference.\n\nAssessment: ✅ Supported, though qualitative.\n\nDiversity is clearly demonstrated; “scientific discovery” phrasing is aspirational but grounded in observed exploration behavior.\n\n🧩 5. Claim:\n\nLLM-generated solvers are interpretable and allow insight into errors and reasoning.\n\nEvidence:\n\nReaction–Diffusion case study (§ 5.7):\n\nThe paper inspects code to show that LLMs incorrectly discretize the reaction term (missed analytical sub-solution).\n\nThis post-hoc inspection is only possible because outputs are human-readable code.\n\nAssessment: ✅ Credible and clearly demonstrated.\n\nConcrete example substantiates interpretability advantage over black-box neural solvers.", "audience_yn": "Yes", "audience_explanation": "Relevance to TMLR Audience\n1. Strong alignment with machine learning research interests\n\nTMLR (Transactions on Machine Learning Research) targets readers who care about:\n\nAdvances in learning paradigms (not only model training but also inference, reasoning, and agentic workflows),\n\nIntersections of ML with scientific domains, and\n\nEmpirical rigor and reproducibility.\n\nCodePDE fits squarely within these aims because it:\n\nExtends ML into scientific reasoning and simulation, a rapidly expanding area (ML for science).\n\nProposes a novel inference framework rather than yet another model architecture — an idea-driven contribution consistent with TMLR’s scope.\n\nProvides extensive, reproducible experiments across 16 LLMs and 5 PDE families.\n\nSo it would interest readers who follow both scientific ML and LLM reasoning trends.\n\n2. Bridges two active research communities\n\nScientific Machine Learning / Physics-Informed ML: Researchers developing PINNs, FNOs, or operator-learning methods would be intrigued that inference-only LLMs can sometimes outperform trained neural PDE solvers.\n\nLLM Reasoning & Agentic Systems: Those studying chain-of-thought, self-refinement, and test-time scaling will see this as a concrete, quantitative domain where those mechanisms demonstrably matter.\n\nThis bridge between symbolic scientific computing and LLM reasoning is rare — and therefore valuable to a broad TMLR audience.\n\n3. Provides empirically grounded insights, not hype\n\nThe study contributes:\n\nRigorous benchmarks (PDEBench, FNO datasets),\n\nQuantitative analyses (nRMSE, convergence orders),\n\nAblation results demonstrating which inference-time components truly matter.\n\nTMLR reviewers and readers appreciate clear evidence over novelty for its own sake, and this work delivers that.", "requested_changes": "1. Clarify Reference Solvers and Ground Truth (Critical)\n\nAlthough the paper lists which solvers were used for each PDE (e.g., analytical, spectral, Strang splitting), it doesn’t specify:\n\nThe numerical order, grid size, or time-stepping details of each reference.\n\nHow “ground truth” solutions were generated for error computation (analytical vs. high-resolution numerical vs. dataset-provided).\n\n➡ Proposed adjustment: Add a concise table summarizing for each PDE family:\n\nThe reference solver type, discretization order, and grid/time step used;\n\nHow the ground truth reference is defined and validated.\n\nThis is critical for fair comparison and reproducibility.\n\n2. Quantify Inference and Compute Costs (Critical)\n\nThe paper reports code execution runtime (solver speed), but not the total inference cost:\n\nNo data on token usage, number of LLM calls per debugging/refinement loop, or wall-clock inference time.\n\nTest-time scaling (best-of-n) increases cost, but the trade-off is not quantified.\n\n➡ Proposed adjustment: Include a short table summarizing inference cost per PDE family:\n\nAverage tokens or latency per solver (single-shot vs. 4-round debug vs. refined vs. best-of-n=16).\n\nCompare total cost vs. classical or neural baselines.\n\nThis is essential for evaluating real-world feasibility.\n\n3. Address Oracle Selection in Best-of-n Scaling (Critical)\n\nThe best-of-n scaling experiments select the lowest nRMSE among candidates, assuming access to ground truth at inference time — effectively an oracle.\n\n➡ Proposed adjustment:\nAdd a variant where solver selection uses unsupervised proxies (e.g., residual norm, stability, runtime errors) instead of true nRMSE, or at least acknowledge this limitation explicitly.\n\nThis change is important for methodological transparency and external validity.\n\n4. Include Details on Execution Environment and Safety (Important but Non-Critical)\n\nThe framework executes LLM-generated Python code, yet the paper omits details on:\n\nWhether the runs were sandboxed,\n\nTimeout/memory protections,\n\nHow unsafe imports or file operations were handled.\n\n➡ Proposed adjustment:\nAdd a short paragraph (even in supplement) explaining code execution safety and reproducibility—e.g., “All generated solvers were executed in Docker containers with restricted I/O and 2 GB memory limits.”\n\nThis isn’t a blocker for acceptance but strengthens trust and replicability.", "broader_impact_concerns": "The ethical risks associated with this work are relatively low, as it focuses on automated code generation for scientific computing rather than on human data or decision-making. However, a few broader-impact aspects could be better addressed. First, since the framework executes LLM-generated code, there are potential safety and security concerns (e.g., malicious or unstable code execution, unbounded resource use) that merit explicit discussion of sandboxing, auditing, or containment strategies. Second, the environmental cost of extensive inference-time sampling (best-of-n and multi-round refinement) could be acknowledged, as repeated LLM calls substantially increase compute usage. Finally, the authors might comment briefly on the reliability and verification of automatically generated scientific code if such methods are adopted in critical applications (e.g., engineering or climate modeling).\n\nThese are not severe ethical shortcomings, but adding a short Broader Impact paragraph addressing code-execution safety, computational sustainability, and responsible deployment of automated scientific solvers would strengthen the submission.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "cvzkMplvLq", "reviewer": "Reviewer_vX6Y", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "## Summary\nThis paper proposes CodePDE, an inference framework for generating numerical PDE solvers using LLMs. The framework consists of five stages—task specification, code generation, automated debugging, evaluation, and solver refinement—and is applied to 16 LLMs across five representative PDE families.\n\nThe authors further analyze model behaviors including solver correctness, numerical scheme diversity, convergence properties, and the effects of test-time scaling. The work offers several interesting empirical insights, such as the trade-off between solver reliability and solver sophistication, and differing strengths across reasoning-oriented vs refinement-oriented LLMs.\n\n## Strengths\n1. Framing PDE solving as LLM-driven code generation is a fresh and potentially impactful perspective. This paradigm opens useful avenues for automation, interpretability, and hybrid symbolic–numerical reasoning in scientific computing.\n2. The experiments cover 16 prominent LLMs, 5 PDE families with varying difficulty and characteristics, Multiple metrics, and Test-time scaling curves. This is a very comprehensive evaluation of LLM-based PDE solver.\n3. The CodePDE pipeline is clear, modular, and pragmatic. The inclusion of iterative debugging and refinement significantly increases solver reliability and performance, demonstrating realistic agentic capabilities.\n\n## Weaknesses\n1. While useful, the framework largely combines generic inference techniques—prompting, best-of-n sampling, iterative debugging, and LLM refinement—without introducing PDE-specific algorithms or new inference methods.\n2. Missing ablations on PDE-specific prompting or model knowledge. It would be highly informative to explore: (1) Providing PDE textbooks or numerical methods summaries to the model (2) Enforcing numerical constraints (CFL, positivity, mass conservation) within prompts (3) The impact of adding symbolic tools\n3. Some key components lack precise algorithmic specification. For example, The refinement stage does not clearly describe how the nRMSE feedback is structured for the model.", "claims_and_evidence_yn": "Yes", "claims_explanation": "N/A", "audience_yn": "Yes", "audience_explanation": "N/A", "requested_changes": "See weaknesses.", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Decision: Accept (poster)