PaperlayerView exact revision on arXiv

Source-linked AI summary

Grid-Free Monte Carlo for Time-Dependent Diffusion

Zihong Zhou, Rohan Sawhney, Eugene d'Eon, Wojciech Jarosz

arXiv:2609.12306v1cs.GRcs.NAmath.NA

TL;DR

Transient diffusion requires resolving time-dependent behavior, but conventional approaches face meshing and time-stepping burdens, and grid-free WoS and WoSt have largely been limited to steady-state problems. This paper extends those walks to heat equations with initial conditions and time-dependent inputs using finite time budgets and improved kernel samplers. The method directly estimates requested times; short-time walks often take fewer steps, and shared walks improve efficiency across multiple target times.

  • Problem

    Transient heat-equation solvers must handle spatial meshing and sequential time stepping, while grid-free WoS and WoSt remain largely limited to steady-state problems.

    Source: §1 · p. 2

  • Method

    The method extends WoS and WoSt with finite-time random walks that sample exit times and use low-bias, tabulation-free exit-time and rejection samplers for time-dependent inputs.

    Source: §1 · p. 2

  • Results

    At short times, walks require fewer steps than steady-state counterparts because they often exhaust their time budgets before reaching the boundary; shared walks improve efficiency when evaluating multiple target times.

    Source: §1 · p. 3

  • Takeaways & Limitations

    The solver estimates arbitrary spacetime queries without volumetric meshing or time stepping, while retaining parallel, progressive, and output-sensitive evaluation.

    Source: §7 · p. 15

    Monte Carlo variance, particularly at low walk counts, remains the primary limitation.

    Source: §7 · p. 15

Abstract

from arXiv · show

Many scientific applications require modeling how diffusive systems evolve over time, not merely their eventual steady states. While conventional steady-state analysis of partial differential equations (PDEs) on complex geometries is already hindered by costly volumetric meshing, transient analysis further requires sequential time stepping and careful step size selection. Grid-free Monte Carlo solvers such as walk on spheres (WoS) and walk on stars (WoSt) avoid this meshing bottleneck but remain largely limited to steady-state problems. We generalize WoS, for pure Dirichlet problems, and WoSt, for mixed Dirichlet--Neumann problems, to heat equations with initial conditions and time-dependent source and boundary data. We equip each random walk with a finite time budget and sample an exit time at every spatial step. If the exit time exceeds the remaining budget, the walk samples an interior point and evaluates the initial condition; otherwise, it continues with a reduced budget, accumulating source and boundary contributions. Our main technical contribution is a suite of kernel sampling and variance reduction techniques, including a low-bias, tabulation-free exit time sampler and efficient rejection samplers. Unlike grid-based transient solvers, our method directly estimates the solution at any requested time without volumetric meshing or sequential time marching. It also retains the parallel, progressive, and output-sensitive evaluation of WoS and WoSt while eliminating time step selection and temporal discretization bias entirely. Finally, we show how sharing walks enables efficient estimates at multiple target times.

1 Introduction

Transient diffusion matters because its evolution toward equilibrium can reveal when and where heat accumulates, but conventional solvers face costly spatial resolution and sequential time stepping. The paper extends WoS and WoSt to transient heat equations using finite-time random walks, with kernel sampling techniques and shared walks supporting direct spacetime estimates.

1 Introduction

  • Finite element methods for transient problems must resolve localized features while selecting sequential time steps, with under-resolution and step-size choices introducing spatial and temporal error.The initial condition and time-dependent sources can create localized interior features that the mesh must resolve.Source: §1 · p. 2
  • WoS and WoSt avoid volumetric meshing but have remained largely limited to steady-state PDEs.WoS addresses Dirichlet problems, while WoSt extends the framework to mixed Dirichlet–Neumann problems.Source: §1 · p. 2
  • Finite-time WoS and WoSt walks sample an exit time at each spatial step, continuing with a reduced time budget or terminating inside the domain to evaluate the initial condition.Walks reaching the Dirichlet boundary evaluate boundary data; source and Neumann contributions are sampled along the walk.Source: §1 · p. 2, Figure 2 · p. 2
  • A tabulation-free, low-bias exit-time sampler and efficient rejection samplers support time-dependent kernels, alongside variance reduction and shared walks for multiple target times.The samplers exploit distinct short-, intermediate-, and long-time behavior of the heat kernel and its derivatives.Source: §1 · p. 2
  • At short times, walks often exhaust their budgets before reaching the boundary and require fewer steps than steady-state walks; at long times, they approach steady-state behavior with modest additional sampling cost.Efficiency depends on both average walk length and exit-time sampling overhead.Source: §1 · p. 3, Figure 3 · p. 3
  • Monte Carlo noise decreases with more walks, but noisy estimates make the method less competitive for dense global solutions, where conventional grid-based solvers remain effective.The method is particularly attractive for sparse or localized queries, short horizons, and progressive evaluation.Source: §1 · p. 3

2 Related Work

Related work contrasts discretized deterministic solvers with grid-free Monte Carlo methods, which have largely addressed steady-state problems. The paper positions its extensions and kernel samplers against the limitations of both families.

2.1 Grid-Based Deterministic Solvers

  • Finite difference and finite element methods discretize space, while boundary element methods reduce spatial discretization to the boundary but can incur dense-system costs and domain integrals for volumetric inputs.Transient FEM also requires time stepping, with accuracy, stability, and cost depending on integrator and step size.Source: §2.1 · pp. 3–4

2.2 Grid-Free Monte Carlo Solvers

  • The paper extends grid-free evaluation to continuous spacetime and develops robust heat-kernel and exit-time sampling, including short-time method-of-images approximations and long-time spectral expansions.Its exit-time sampler uses sampled exponential variables and a moment-matched Gamma approximation of the tail; targeted strategies also address localized initial conditions and time-independent source and Neumann data.Source: §2.1 · p. 4, §2.2 · pp. 4–5
  • Grid-free WoS avoids volumetric meshes and Brownian-path time discretization, and WoSt extends the approach to mixed boundary conditions, but existing methods remain focused primarily on steady-state problems.Earlier transient approaches have distinct constraints: WoHB uses a temporal epsilon-shell, while other methods restrict inputs or can exhibit high variance.Source: §2.2 · p. 4

3 Background

The background reviews grid-free random walks for steady-state and heat-equation problems, then introduces the notation and kernel-based steps underlying WoS, WoSt, and floating random walks.

3 Background

  • The background surveys WoS, WoSt, and floating random walks for steady-state Dirichlet, mixed-boundary, and heat-equation problems.Source: §3 · p. 5

3.1 Notation

  • The notation distinguishes query and interior points, remaining time budgets, elapsed time, normalized time, and the assumed unit diffusivity.Other positive constant diffusivities are addressed later in the paper.Source: §3.1 · p. 5

3.2 Walk on Spheres

  • WoS recursively samples exits from the largest ball centered at the current point, accumulating source contributions and evaluating Dirichlet data near the boundary.Uniform sphere sampling cancels the boundary weight.Source: §3.2 · p. 5

3.3 Walk on Stars

  • WoSt uses a star-shaped region bounded by a spherical portion and Neumann boundary, sampling rays to advance the walk or continue from Neumann data.Its radius is chosen using distances to the Dirichlet boundary and a Neumann-boundary silhouette point.Source: §3.3 · pp. 5–6

3.4 Floating Random Walks

  • Floating random walks extend WoS to the homogeneous heat equation with Dirichlet data, sampling an exit time alongside each spatial step.If the sampled exit time fits the remaining budget, the walk advances; otherwise it returns the constant initial value.Source: §3.4 · p. 6
  • The heat kernel describes heat spread at a given time, whereas the harmonic Green’s function integrates that kernel over all elapsed times.For a ball, exit location and exit time distributions factor into a uniform sphere location and an exit-time density.Source: Figure 5 · p. 6, §3.4 · p. 6

4 Method

The method develops a grid-free solver for the heat-equation initial-boundary value problem, extending WoS and WoSt to time-dependent inputs.

4 Method

  • The solver progresses from time-dependent WoS with general initial conditions to kernel samplers, mixed-boundary WoSt, time-dependent source and Neumann data, and positive constant diffusivities.Source: §4 · p. 6

4.1 Walk on Spheres with General Initial Conditions

For general initial conditions, the method extends WoS by pairing spatial steps with time-budget and heat-kernel sampling. It offers both time-exhaustion evaluation and per-step initial-condition connections.

4.1 Walk on Spheres with General Initial Conditions

  • For nonconstant initial conditions, the heat-kernel integral cannot factor out the initial value, so the walk uses Monte Carlo integration and a recursive estimator.Source: §4.1 · p. 6

4.1.1 Initial Condition at Time Exhaustion.

  • Sampling the full exit-time density selects whether the walk exits before its budget or runs out of time inside the ball and evaluates the initial condition.On exit, it continues with reduced time; on exhaustion, it samples an interior point from the time-conditioned heat kernel.Source: §4.1.1 · p. 7

4.1.2 Initial Condition Connections.

  • The alternative initial-condition connection strategy evaluates the initial condition at every walk step, with a contribution weighted by the probability of remaining inside the ball.A single walk therefore accumulates initial-condition contributions from every visited ball.Source: §4.1.2 · p. 7, Figure 6 · p. 7
  • Initial-condition connections are typically most effective for localized initial conditions at intermediate times.Source: §4.1.2 · p. 7

4.2 Exit Time and Heat Kernel Sampling

The section develops exit-time and heat-kernel samplers needed to evaluate transient walks across time scales. Its exit-time sampler avoids tabulation, while rejection sampling supports efficient interior-point draws.

4.2.1 Exit Time Sampling.

  • The exit-time sampler draws the first K exponential variables and approximates the remaining sum with a Gamma variable matched to its mean and variance.Dimensionless unit-ball samples are scaled by R^2 for each walk step.Source: §4.2.1 · pp. 7–8
  • 5.8 × 10^-6 maximum relative error over moments 3–12, versus 5.6 × 10^-3 for a 1000-entry CDF table; the sampler is also an order of magnitude faster than WoMS.The sampler exactly preserves the exit time’s mean and variance.Source: §4.2.1 · p. 8, Figure 7 · p. 9

4.2.2 Heat Kernel Sampling.

  • For stable heat-kernel evaluation, the method uses images for τ ≤ 0.01 and a 16-term spectral expansion for longer times.The spectral representation converges rapidly at long times but can require many terms and suffer floating-point cancellation at short times.Source: §4.2.2 · pp. 8–9
  • Rejection samplers for the time-conditioned heat-kernel radial density maintain low rejection rates across short-, intermediate-, and long-time regimes, adding no bias beyond numerical kernel approximations.The long-time proposal uses a Beta(3, 2) density that closely matches the limiting radial shape.Source: §4.2.2 · p. 9

4.3 Walk on Stars with General Initial Conditions

The section extends time-dependent walks on stars to mixed Dirichlet–Neumann problems with general initial conditions. A finite time budget supplies a stopping rule even for pure Neumann problems, and the initial-condition estimate remains unbiased.

4.3 Walk on Stars with General Initial Conditions

  • For the initial-condition formulation, the section assumes zero source and Neumann data and applies the heat-equation Green representation to a star-shaped region.The time-dependent Poisson kernel on the Neumann portion measures flux through an interior surface, not a first-exit distribution.Source: §4.3 · pp. 9–10

4.3.1 Recursive Walk.

  • In pure Neumann problems, time-budget exhaustion terminates the walk inside the domain and the initial condition anchors the solution, unlike the steady-state solution defined only up to an additive constant.Shorter budgets retain high-frequency initial-condition detail, while longer budgets capture a more averaged profile.Source: §4.3.1 · p. 10, Figure 8 · p. 10

4.3.2 Initial Condition Contribution.

  • To avoid a geometry-dependent heat-kernel normalization over the star, the method samples from the enclosing ball’s normalized kernel and extends the initial condition by zero outside the star.Assigning zero to samples outside the star preserves an unbiased estimate of the initial-condition integral.Source: §4.3.2 · p. 10

4.4 Time-Dependent Source and Neumann Data

The section adds time-dependent source and Neumann contributions to the walk-on-stars representation. It samples source terms from an enclosing-ball space-time kernel and estimates Neumann terms by sampling boundary locations and elapsed times.

4.4 Time-Dependent Source and Neumann Data

  • For source contributions, the method samples elapsed time from a normalized enclosing-ball density, then samples an interior location from the time-conditioned heat kernel.In WoSt, the source is extended by zero outside the star region; in WoS, all sampled source locations lie inside the ball.Source: §4.4 · p. 11
  • 15% efficiency improvement is reported for the time-independent source optimization in a source-dominated manufactured-solution test.The estimator exhibits predictable Monte Carlo convergence as walks per point increase.Source: Figure 10 · p. 11, §4.4 · p. 11

4.4.2 Neumann Boundary Condition.

  • For Neumann contributions, the method samples a boundary location and then an elapsed time proportional to the heat kernel, accumulating the contribution at each walk step.The resulting estimator supports spatially and temporally varying Neumann data and is demonstrated on a pure Neumann problem.Source: §4.4.2 · pp. 11–12

4.5 Constant Non-Unit Diffusivity

The section handles any positive constant diffusivity by rescaling time and running the unit-diffusivity solver.

4.5 Constant Non-Unit Diffusivity

  • For positive diffusivity κ, the solver uses initial time budget κt and converts sampled elapsed times back by dividing by κ, without changing spatial sampling or geometric queries.The rescaling multiplies the source by 1/κ, reparameterizes boundary data in time, and leaves the initial condition unchanged.Source: §4.5 · p. 12

5 Improving Efficiency

The efficiency improvements target localized initial conditions, time-independent data, and repeated target-time queries. Together, they reduce unnecessary sampling and amortize work without adding bias to shared-walk estimates.

5 Improving Efficiency

  • The efficiency suite also skips zero initial-condition contributions and combines sampling strategies for localized data.It further omits temporal samples for time-independent source and Neumann data and shares walks across target times.Source: §5 · p. 12

5.1 Localized Initial Conditions

  • Multiple importance sampling combines heat-kernel and region proposals, working across time budgets and initial-condition support sizes.A support-bounding-volume test skips initial-condition evaluations when the support is disjoint from the current ball.Source: §5.1 · p. 12, Figure 11 · p. 12

5.2 Time-Independent Source and Neumann Data

  • For time-independent source and Neumann data, integrating over elapsed time removes temporal samples while preserving each estimator’s expectation and variance.The source is sampled from a density based on the time-integrated heat kernel; the Neumann contribution reuses the hierarchical boundary sampler.Source: §5.2 · p. 12

5.3 Shared Walks Across Multiple Times

  • Sharing spatial steps and exit-time samples across target budgets amortizes geometric queries without additional bias, yielding a 3× runtime-efficiency improvement in the illustrated five-budget example.Each unfinished estimate updates its own remaining budget and evaluates time-dependent contributions separately.Source: §5.3 · pp. 12–13, Figure 12 · p. 13

6 Evaluation

The evaluation compares the method with other grid-free Monte Carlo solvers and FEM, then demonstrates it on a transient ITER heat-conduction problem. Results cover sampling behavior, discretization tradeoffs, and runtime for a sparse temperature-history query.

6.1 Implementation

  • Early time-budget exhaustion often makes transient solves faster at short times, while initial-condition connections favor localized data when their variance reduction offsets their cost.At longer times, time-dependent kernel sampling adds modest overhead relative to steady-state WoS and WoSt.Source: §6.1 · p. 13, Figure 13 · p. 13

6.2.1 Walk on Boundary.

  • WoB’s sample variance was approximately 107 for the larger time budget, versus 10^-3 for our method, despite WoB tracing roughly 5× more walks per point.The equal-runtime comparison used WoB on an NVIDIA RTX 3060 Laptop GPU and the method on the CPU configuration described in the implementation section.Source: §6.2.1 · p. 14, Figure 14 · p. 14

6.2.2 Walk on Heat Balls.

  • For the smaller time budget, our method ran approximately 3× faster than WoHB, whose smaller spacetime heat balls require more walk steps.WoHB’s temporal epsilon-shell termination retains bias, whereas the method can sample a terminal point from the time-conditioned heat kernel.Source: §6.2.2 · p. 14, Figure 14 · p. 14

6.3 Comparison with the Finite Element Method

  • Coarse FEM meshes under-resolve localized initial conditions, while aggressive refinement increases runtime and memory yet leaves visible aliasing artifacts.The method operates on the boundary mesh and avoids resolving localized interior data on a volume mesh.Source: §6.3 · p. 15
  • FEM time-step choice trades accuracy, stability, and cost; even an unconditionally stable method needs small steps to resolve rapidly changing data.The method evaluates arbitrary requested times without time stepping or temporal discretization bias.Source: §6.3 · p. 15

6.4 Thermal Analysis

  • The ITER example substitutes fixed-temperature Dirichlet pipe walls for unsupported Robin conditions and uses a triangulated boundary mesh rather than FEM’s tetrahedral volume mesh.The authors characterize it as a representative transient thermal analysis, not a complete thermal-hydraulic model.Source: §6.4 · p. 15, Figure 15 · p. 15
  • One million shared walks reconstructed 5,000 probe-time values over [0, 250] s in about three minutes; coarse-step FEM took about one hour and overestimated temperature past the safety limit.The FEM prediction improved with a smaller time step, at substantially greater cost.Source: §6.4 · p. 15

7 Conclusion and Future Work

The paper extends WoS and WoSt to heat-equation initial-boundary value problems, estimating arbitrary spacetime queries without volumetric meshing or time stepping. Its scope remains limited by Monte Carlo variance and constant-diffusivity linear heat equations, with broader boundary conditions and variable coefficients left for future work.

7 Conclusion and Future Work

  • The time-dependent WoS and WoSt solver handles general initial conditions and time-dependent source and boundary data, with shared walks supporting estimates at multiple target times.It evaluates arbitrary spacetime queries without volumetric meshing or time stepping while retaining parallel, progressive, and output-sensitive evaluation.Source: §7 · p. 15
  • Monte Carlo variance is the primary limitation, and Neumann-dominated problems with large time budgets can be costly because walks repeatedly interact with the reflecting boundary.Additional walks reduce noise predictably and run efficiently in parallel; adapting variance-reduction methods and distributed or GPU implementations are future directions.Source: §7 · pp. 15–16
  • On the gear example, FEM exhibits spatial aliasing and mesh artifacts even after refinement to 4 million tetrahedra, while Monte Carlo avoids volume meshing and time stepping.FEM runtimes exclude mesh generation, and further refinement fails due to integer overflow; Monte Carlo estimates converge as walks per point increase.Source: Figure 16 · p. 16
  • On the highly refined mesh, RK4 requires a 10^-7 time step for stability; Monte Carlo avoids time-discretization bias, with sampling variance decreasing as walks increase.Backward Euler is unconditionally stable, but large steps introduce numerical diffusion and small steps require tuning.Source: Figure 17 · p. 16
  • The formulation assumes a linear heat equation with constant diffusivity and mixed Dirichlet–Neumann conditions; variable coefficients and transient Robin or radiative conditions remain extensions.Efficient estimators for spatial, temporal, and parameter derivatives could support sensitivity analysis, optimization, and inverse problems.Source: §7 · p. 16
  • For a moving Gaussian source, Monte Carlo captures the overall heat profile without time stepping, while reducing FEM step size recovers the profile at greater cost.Increasing walks per point predictably reduces Monte Carlo sampling noise.Source: Figure 18 · p. 17
Loading 2609.12306v1…