Source-linked AI summary

LLM-Assisted Generation of Transparent, Open-Source Multiphysics Models of Electrochemical Devices

Sebastian Castro, Maya F. Schuchert, Spencer A. McCluskey, Eric W. Lees, Justin C. Bui

arXiv:2610.10320v1physics.chem-phcs.AIphysics.comp-ph

TL;DR

Thorough, precise implementation guides are needed to constrain agent-built multiphysics models and prevent ambiguity. This paper specifies such a guide for electrochemical modeling and uses it to build open-source Julia models verified against an independently implemented COMSOL model. The work also shows that silent physics errors can escape staged checks, underscoring the importance of human-specified physics and verification.

  • Problem

    Implementation guides must specify equations, parameters, boundary conditions, and build plans precisely to prevent ambiguity that can lead to model errors and inconsistencies.

    Source: The implementation guide · p. 17
  • Method

    The workflow uses a machine-readable guide specifying model physics and numerical methods, together with staged verification criteria, to direct agent construction of multiphysics models.

    Source: A specification-driven, staged build · p. 5, The implementation guide · p. 17
  • Results

    The one-dimensional, 1,480-degree-of-freedom electrochemical model reproduced expected physics, was calibrated to experimental data, and was verified against an independently implemented COMSOL model.

    Source: Discussion · p. 15
  • Takeaways & Limitations

    Open, inspectable models let researchers who understand the device but are not numerical-method specialists build, troubleshoot, extend, and share multiphysics models while retaining control of physics and implementation.

    Source: Discussion · p. 16

    Staged checks can miss numerically silent physics errors, as a planted diffusivity defect converged and passed global charge conservation.

    Source: Error handling with multiphysics agents · p. 13, Discussion · p. 15

Abstract

from arXiv · show

Multiphysics continuum models are powerful tools for studying electrochemical devices, enabling in silico reactor design and resolution of local pH, potential, and concentration fields that govern device performance but are difficult to measure experimentally. However, constructing such models requires substantial numerical expertise or reliance on proprietary software. Here, we show that frontier large language model agents can remove this implementation burden while keeping the underlying physics under researcher control. Using one-dimensional electrochemical CO2 reduction to CO in a porous gas diffusion electrode as a test case, we develop a machine-readable, human-specified modeling harness containing governing equations, parameters, numerical methods, logical build stages, and human-verifiable checkpoints. From this specification, the agent reproducibly constructs complete multiphysics models in open-source Julia. Independently built models, including fully autonomous agent-built models, agree with an equivalent COMSOL implementation to within 0.7% of the peak CO partial current density, and with one another to within 0.04%. Systematically planted errors demonstrate the importance of explicit specifications for reproducibility and reveal the agent's capabilities and limitations in debugging model physics. This framework establishes a more transparent approach to multiphysics modeling in which physical descriptions and governing equations, rather than specialized code, become the primary inputs for computational model development.

Results

The staged guide adds coupled physics through five build stages, with human checkpoints to assess and correct implementation. The resulting polarization behavior and spatial fields provide checks against expected electrochemical behavior.

A specification-driven, staged build

  • Five stages increase physical complexity from transport, buffer reactions, and kinetics to a fully electrostatically coupled system, with human-reviewed stop gates between stages.The model grows from 400 to 1,480 degrees of freedom on the production mesh; gates prescribe quantitative physicality checks and diagnostic plots.Source: A specification-driven, staged build · pp. 5–6
  • The polarization curves show FECO peaking at intermediate potentials and declining at high potentials, where hydrogen evolution outcompetes CO2 reduction.The behavior reflects coupled CO2 reduction kinetics, CO2 transport, and acid–base consumption in the porous electrode.Source: A specification-driven, staged build · pp. 6–7
  • Stage 2 mass-transfer boundaries reduce CO partial current through CO2 depletion and pH buildup, while Stage 3 gas-to-liquid transfer raises it again.Stage 4 gas momentum adds a small correction; fitting at Stage 4 brings iCO toward experimental data, and Stage 5 migration raises it by mitigating CO2 salting-out.Source: A specification-driven, staged build · p. 7
  • Stage-5 profiles show pH highest near the gas-diffusion-layer side, where CO2 consumption generates OH−, and dissolved CO2 lowest there before increasing toward the electrolyte.Gas-phase CO2 enters from the flow channel, diffuses through the GDL, and is consumed in the catalyst layer, while generated CO diffuses out; these microscale fields are difficult to measure experimentally.Source: A specification-driven, staged build · p. 7

Agreement with COMSOL and across builds

Three independent Julia builds closely reproduce a physically equivalent COMSOL model and agree closely with one another. The reported residual differences are small relative to peak CO current, though small solver or mesh differences may contribute to the offset.

Agreement with COMSOL and across builds

  • Across three independent Julia builds, the average of per-build iCO RMSDs against COMSOL was 2.27 mA cm−2 at 44 sweep potentials, 0.7% of the approximately 330 mA cm−2 peak CO current density.Source: Agreement with COMSOL and across builds · p. 7
  • Julia iCO values were up to 2.2% above COMSOL above 5 mA cm−2 and about 0.9% above at peak.The authors suggest slight mesh-generation or numerical-implementation differences as possible sources of the discrepancy.Source: Agreement with COMSOL and across builds · p. 8
  • Across three builds on a common 89-point voltage sweep, pooled RMSDs were 0.12 mA cm−2 for CO current and 0.11 mA cm−2 for H2 current, with CO faradaic efficiency agreeing within 0.01 percentage points.The fitted exchange current density was 106.47 ± 0.11 A m−2 and reorganization energy was 1.6111 ± 0.0004 eV.Source: Agreement with COMSOL and across builds · p. 8

Solver implementation and agent cost–benefit deviation

The Julia implementation uses a customizable solver stack and performs competitively in the reported sweep, although differences in discretization and computing environment preclude a direct speed comparison. Agents also substituted simpler solver components when they judged them equivalent, so human specification and checking remain important.

Solver implementation and agent cost–benefit deviation

  • The 89-potential Julia sweep took about 13 s, versus 21 s for COMSOL's 61-potential sweep at similar problem size and compute budget, but the comparison is not direct.The implementations differ in discretization, linear-solve library, and solve environment.Source: Solver implementation and agent… · p. 9, Figure 4 · p. 9
  • Exploiting model sparsity provided nearly all solver speed gains, with linear solving and Jacobian assembly each about two orders of magnitude faster than dense alternatives.The solver ablation configurations returned the same polarization curve to the precision reported, so the bars compare cost.Source: Solver implementation and agent… · p. 9, Figure 4 · p. 9
  • Agents replaced the specified analytical block Jacobian and MUMPS solver with forward-mode automatic differentiation and UMFPACK, judging the simpler alternatives at least as fast.Those substitutions were non-detrimental in solve time and were adopted, but the authors caution that agents may make such changes unless humans explicitly require and check the specified numerics.Source: Solver implementation and agent… · pp. 9–10

Diagnostic transparency during model development

Model development exposed both convergence problems and silent physics errors, with human review helping identify issues that agents did not consistently diagnose. In mature builds, transparent diagnostics let agents connect numerical behavior to local physics and solver residuals.

Diagnostic transparency during model development

  • As the guide was refined, early critical failures and substantial human input gave way to final replicates requiring minimal intervention and having no critical failures.The development challenges and agent responses are summarized in Fig. 5.Source: Diagnostic transparency during model… · p. 10, Figure 5 · p. 11
  • A sign error in the Stage-2 CO2 mass-transfer boundary condition recurred in five of six early builds because the guide omitted its sign convention.The error reversed the flux direction; some agents found it after nonconvergence, while others needed human prompting.Source: Diagnostic transparency during model… · p. 10
  • Two late-stage builds admitted spurious, partially converged polarization points because convergence tolerances were loose; humans spotted them visually and agents tightened the absolute residual tolerance.The agents independently converged on the same fix.Source: Diagnostic transparency during model… · p. 12
  • Agents struggled to identify hidden physics bugs that did not block convergence, including an ambiguous diffusivity definition that differed by a factor of 3.54.Cross-comparison of divergent replicates helped expose the conflicting Bruggeman correction definitions, after which the guide was revised.Source: Diagnostic transparency during model… · p. 12
  • An agent could inspect the sparse Stage-5 Jacobian, report its condition number, locate reaction-rate and HER-dominated cells, and check current continuity against boundary flux.It could also explain these diagnostics through local CO2 availability and pH, and trace variables and residuals to the governing physics.Source: Diagnostic transparency during model… · p. 12

Error handling with multiphysics agents

In a controlled study, agents resolved three of four planted errors, but a numerically silent constitutive error escaped detection and persisted through fitting. Agent review caught common implementation mistakes, while the silent failure shows why independent references remain necessary.

Study design

  • The controlled study tested four naturally encountered defects—a boundary-condition sign reversal, omitted conservation term, omitted constitutive correction, and unit inconsistency—against an unmodified-guide control.Builds used the same prompt and Claude Opus 4.8 at maximum effort.Source: Error handling with multiphysics agents · pp. 12–13

Results

  • Three of four planted defects were resolved before final development; two without human intervention and one after a human hint, while Case 3 went undetected.The sign reversal made the boundary-value problem unsolvable, and the unit inconsistency produced transport coefficients thirty-fold too high; the agent rejected the latter on magnitude.Source: Error handling with multiphysics agents · p. 13
  • Case 3 reproduced experimental polarization despite converged residuals and charge conservation of 2.4 × 10^-11, while its fitted exchange current density was 41% below every other build.The fitting step absorbed most of the error; an independent COMSOL reference helped resolve silent guide errors.Source: Error handling with multiphysics agents · p. 13
  • Agent review caught all syntax errors and most code-implementation convergence errors without human input, alongside common sign, unit, and math mistakes.These checks provide a safety net, but do not catch every physical-model error.Source: Error handling with multiphysics agents · pp. 13–14

Fully autonomous builds

Six fully autonomous builds produced virtually the same model output after substantial development of the implementation guide. The workflow shifts numerical implementation toward explicit physical specification and verification, while leaving researchers responsible for the model’s assumptions.

Fully autonomous builds

  • Six autonomous builds produced virtually the same output, with fitted i0 = 106.56 ± 0.07 A m−2 and λ = 1.6113 ± 0.0002 eV for Opus and i0 = 106.51 ± 0.01 A m−2 and λ = 1.6108 ± 0.0003 eV for Sonnet.Builds were instructed to bypass stop gates, use only the guide, and receive no human input; the figure compares polarization curves and iCO deviations with COMSOL and experimental data.Source: Fully autonomous builds · p. 14, Figure 6 · p. 14

Discussion

  • Autonomous builds became possible only after substantial iteration on the implementation guide, and a consistently stated but physically incorrect assumption can still be reproduced and obscured by parameter fitting.Independent builds, external references, and human physical intuition remain safeguards.Source: Fully autonomous builds · p. 14, Discussion · p. 15
  • The model couples transport, carbonate chemistry, kinetics, flow, and electrostatics in a stiff 1,480-degree-of-freedom system calibrated to experimental data and verified against an independently implemented COMSOL model.The demonstrated model is one-dimensional and uses a macrohomogeneous porous-electrode description selected for studying polarization, selectivity, and concentration profiles.Source: Discussion · p. 15
  • The implementation guide records equations, boundary conditions, conventions, numerical methods, verification criteria, and expected physical magnitudes as an inspectable artifact independent of a particular implementation.Staged checks catch many implementation and convergence failures but cannot ensure that the specified physics is correct.Source: Discussion · p. 15
  • The workflow makes researchers’ understanding of the physics the primary interface to model building while keeping the resulting models open-source and fully inspectable.Researchers remain responsible for stating and evaluating the physics precisely.Source: Discussion · p. 16

Methods

The study models a one-dimensional silver-catalyzed CO2 reduction cathode and uses a rebuilt COMSOL model as its benchmark. The model couples electrochemical reactions and buffer chemistry with liquid and gas transport under specified operating and structural assumptions.

Model system and governing equations

  • The one-dimensional CO2 reduction model uses a porous-electrode domain and GDE formulation from Weng et al. [7], MHC kinetics from Lees et al. [36], and a rebuilt COMSOL reference based on Muhieddine Orfali et al. [37].Only the cathode catalyst layer and gas diffusion layer are modeled; the anode and membrane are excluded, with full-cell and higher-dimensional models left for future work.Source: Model system and governing equations · p. 16
  • Two catalyst-surface reactions generate hydroxide, which carbonate buffer equilibria consume, linking current, local pH, dissolved CO2 depletion, and reduced CO2-reduction kinetics.The stated feedback includes mass-transport limitation as local pH rises and CO2 is depleted.Source: Model system and governing equations · p. 16
  • The steady-state, isothermal model at 50 °C couples liquid Nernst–Planck transport, gas Stefan–Maxwell transport with Henry’s law phase transfer, Darcy flow, MHC CO2-reduction kinetics, and Butler–Volmer hydrogen-evolution kinetics.It assumes fixed porosity and liquid saturation, represents the catalyst layer as solid, ionomer/liquid, and gas-filled pore phases, and treats the GDL as dry without anolyte crossover.Source: Model system and governing equations · p. 17

The implementation guide

The implementation guide is a human-specified, machine-readable harness that translates physical understanding into a constrained model-building process. Its precision, staged refinement, and collaborative versioning support auditable implementation.

  • The guide specifies the problem’s parameters, equations, boundary conditions, and scaffolded build plan, making the agent’s role implementation rather than physical inference.The case-study guide spans 18 sections covering equations, parameters, and solver information.Source: The implementation guide · p. 17
  • Explicit equations, complete parameter tables, and specified numerical methods constrain the agent’s algorithmic choices and reduce ambiguity.The methods include finite-volume discretization, a damped Newton solver, and solver initialization and preconditioning protocols for stiff problems.Source: The implementation guide · p. 17
  • Current frontier models require substantial scaffolding because they may take unphysical shortcuts or hallucinate parameters, making human specification and verification essential.The authors say purpose-built agents or future frontier models may handle greater specification flexibility, but no specification is ambiguity-free.Source: The implementation guide · p. 18
  • Iteratively refining a shared, version-controlled guide lets independent builds expose ambiguities or errors, while agreement across builds supports confidence in the model and specification.The authors contrast this collaborative artifact with an opaque, version-locked solver project.Source: The implementation guide · p. 18

Staged construction and stop gates

The model is built through stages of increasing physical complexity, with human-reviewed stop gates before progression. This staged process helps isolate errors and refine the guide before complexity compounds.

  • Stage 5 completes the 1,480-DOF system at the production mesh, with migration enabled and Stage 4-fitted kinetic parameters held fixed.The production mesh has 80 CL and 120 GDL cells; migration raises the simulated CO partial current.Source: Staged construction and stop gates · p. 18
  • The five-stage build progresses from a 400-DOF liquid-phase problem to gas transport, phase transfer, and pressure balance before full electrochemical migration is enabled.A fitting step estimates CO2R and HER kinetic parameters against experimental data at Stage 4, before charged-species coupling increases stiffness.Source: Staged construction and stop gates · p. 18
  • Human reviewers use prescribed quantitative pass/fail criteria and diagnostic plots at stop gates to check physicality and authorize each next stage.Checks include electroneutrality, mass conservation, current continuity, and Newton residual convergence; plots show trends such as rising pH with cathodic current.Source: Staged construction and stop gates · p. 19
  • Testing a coupled-model error in a 400-DOF subproblem is simpler than finding it after it contaminates the 1,480-DOF system, and stop gates also reveal guide ambiguities.The authors say building a specification of this scale from scratch in one prompt is not practical at this stage.Source: Staged construction and stop gates · p. 19

Agents and models

Multiple researchers used the shared guide to build models, while the reported builds also used persistent project memory. The authors compare builds to examine reproducibility despite the probabilistic nature of LLMs.

  • The authors used frontier Claude models for different build stages and chose them over a custom agent to demonstrate feasibility with an established, capable agent.They report that weaker Sonnet-tier models made more errors than Opus in earlier builds from less mature guides.Source: Agents and models · p. 19
  • Researchers across laboratories and with varying experience built models from the same guide, instructing agents to consult only it apart from project memory.The authors found comparing replicates helpful for compensating for LLMs’ probabilistic nature.Source: Agents and models · pp. 19–20
  • All nine reported builds used persistent project memory containing some earlier results, but it lacked λ, HER parameters, polarization curves, and COMSOL output.A tenth build using only the guide and no project memory fitted i0 = 106.59 A m−2 and λ = 1.6109 eV.Source: Agents and models · p. 20
  • The untargeted λ, HER parameters, and full 89-point polarization curves agreed across builds, and the authors report that agreement did not rest on memory’s expected values.Each fit started from the guide’s prescribed initial point and was deterministic.Source: Agents and models · p. 20

Planted-error study

The section describes how planted defects in the specification were assessed and sets out the matched-compute conditions for Julia–COMSOL solve-time benchmarking and mesh-convergence checks.

Planted-error study

  • A planted defect counted as detected only when the build log identified it, named the affected relation, and proposed a correction.Defects were placed in the guide rather than code, and reviewers approved or rejected stage progression by inspecting plots without auditing source code.Source: Planted-error study · p. 20
  • Each build was instructed to follow the guide exactly, stop for human review at each gate, and log errors and discrepancies without being told a defect existed.Explanatory text naming the affected quantity was edited along with a defect to prevent detection from an isolated self-contradiction.Source: Planted-error study · p. 20

Solve-time benchmarking and mesh convergence

  • Julia and COMSOL were timed on the same workstation at matched core budgets and comparable degrees of freedom, with both converging at every swept potential.COMSOL swept 61 potentials at 50 mV spacing and Julia 89 at 25 mV; startup and compilation time were excluded, and Julia ran under WSL.Source: Solve-time benchmarking and mesh… · pp. 20–21
  • The authors document Stage 5 mesh convergence separately and report the fastest of three Julia timing repetitions, with variability given in supplementary material.The two recorded COMSOL runs had the same time.Source: Solve-time benchmarking and mesh… · p. 21

Reporting summary

The reporting summary is linked to the article, while the repository provides the implementation materials and modified guides; funding sources and a funding disclaimer are also stated.

Reporting summary

  • Research-design information is available in the Nature Portfolio Reporting Summary linked to the article.Source: Reporting summary · p. 21

Code availability

  • The repository provides the implementation guide, five-stage Julia source, build logs, planted-error guides, and COMSOL reference materials; running COMSOL requires a v6.4 license and specified modules.Source: Code availability · p. 21

Funding

  • The authors disclose support from the NSF Graduate Research Fellowship Program, Schmidt Science Fellows, the Rhodes Trust, the Peter Wall Legacy Award, and NSERC Discovery Grants; the views need not reflect NSF views.Source: Funding · pp. 21–22
Loading 2610.10320v1…