Solving Partial Differential Equations on Closed Surfaces with Planar Cartesian Grids

J. Thomas Beale

Solving Partial Differential Equations on Closed Surfaces with Planar Cartesian Grids

T0 review · 0 major / 4 minor · reviewed 2026-08-14 · deepseek-v4-flash

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗ Review JSON

Pith's one-line read Representing a closed surface by cut points on a planar grid, with values at secondary points fixed by quadratic interpolation, lets ordinary finite differences solve Laplace-Beltrami diffusion, advection, and shallow-water equations on…

desk verdict A solid, honest methods paper: the cut-point equilibration idea is adapted to surface PDEs with convincing second-order numerics, though the load-bearing interpolation step remains unproved. read the letter →

arxiv 1908.01796 v2 pith:VSGA4EHX submitted 2019-08-05 math.NA cs.NA

J. Thomas Beale This is my paper · ORCID

classification math.NAcs.NA MSC 65M0665M5058J3558J4535Q86

keywords CartesiangridssurfacePDEsLaplace-Beltramioperatorequilibrationlevelsetmethodshallowwaterequationsfinitedifferencemethods

verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

Most PDE solvers on curved surfaces need a mesh or an extension of the equation into a neighborhood of the surface. This paper proposes a simpler route: if the surface is given as a level set, it is represented by the points where the surface crosses the intervals of a regular three-dimensional Cartesian grid. Some of these cut points carry independent unknowns (primary points); the rest (secondary points) get their values from the primary ones through a quadratic-interpolation 'equilibration' step. The paper shows that with this representation, the Laplace-Beltrami operator, surface gradients, divergence, advection, and the shallow-water equations can be discretized by ordinary finite differences in the three coordinate planes, and the observed accuracy is second order on spheres, ellipsoids, and a nonconvex Cassini oval.

What carries the argument

The machinery is the equilibration operator E_h. Given values at primary cut points, E_h fills in values at secondary cut points by solving the system u_s = Π_sp u_p + Π_ss u_s, whose rows are the quadratic interpolation condition (2.1) along the surface; since I − Π_ss is strictly diagonally dominant, the system is invertible. Composing E_h with the finite-difference discretization of the surface differential operator — e.g., the nine-point Laplace-Beltrami stencil Δ_h in (3.11) or (3.15), or the coordinate-plane formulas for surface gradient and divergence — produces the reduced operator Δ_h E_h on primary points. The equilibration is what lets each primary point use the coordinate plane in which the surface normal is nearly aligned with a coordinate axis, and it is what transfers updated values back to secondary points between time steps.

What would settle it

Run the diffusion test of Section 3 on a Cassini oval with a near-self-intersecting value such as b/a = 1.02, where the surface has large curvature; if halving the grid size h no longer reduces the relative L2 error by roughly a factor of four, the claimed second-order accuracy fails.

Watch

Extended reading notes

Core claim

The central claim is that a single, fixed discretization of any smooth closed level-set surface reduces surface PDEs to standard planar finite-difference calculations. The surface is sampled at cut points where it meets grid intervals; primary points are those closest to grid points, and secondary points are tied to primary points by solving a small interpolation system. The discrete Laplace-Beltrami operator is built as a nine-point stencil in the appropriate coordinate plane, and the equilibration operator E_h passes values from primary to secondary points so that stencils can be applied seamlessly as the coordinate system changes across the surface. The paper demonstrates second-order convergence for diffusion and Poisson problems, for a linear advection equation with a known solution, and for the shallow-water equations on a sphere using a MacCormack Lax-Wendroff scheme with artificial viscosity; it also computes low eigenvalues of the discrete spherical Laplacian and shows they converge to the exact values at second order. The author is careful to state that convergence is observed but not proved, and a proof is left for the future.

Load-bearing premise

The scheme's accuracy depends on the quadratic interpolation that fills in secondary-point values being accurate enough not to spoil the second-order finite-difference stencils it is composed with, and this is verified numerically but not proven for surfaces in three dimensions.

Editorial extensions

If this is right

Reaction-diffusion systems on surfaces are a direct extension, since the paper's Laplace-Beltrami discretization can be combined with pointwise reaction terms (Section 5).
Poisson's equation on a closed surface can be solved by augmenting the reduced Laplace-Beltrami operator with one scalar constraint, and the paper shows second-order convergence for this formulation on the sphere.
The method's coordinate-plane formulas for surface gradient, divergence, and material derivative mean that planar finite-volume or finite-difference schemes can be ported to smooth closed surfaces with only the equilibration step added.
Surface integrals can be computed directly from the cut points with a high-order quadrature rule, allowing conserved quantities such as mass and energy to be monitored in long-time simulations without building a separate surface mesh.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

A natural next step, not taken in the paper, is to use higher-order interpolation in the equilibration step; if the accuracy barrier is really the quadratic interpolation, this would open the way to fourth-order surface PDE solvers with the same cut-point structure.
The paper's proof of resolvent positivity covers only closed curves in the plane; if a similar argument can be devised for surfaces in R^3, it would supply the missing maximum-principle stability theory for implicit time stepping.
The surface-discretization machinery could plausibly be combined with a semi-Lagrangian time step by interpolating from neighboring primary points, a possibility the paper mentions but does not develop.
The restriction to smooth level-set surfaces is explicit; applying the method to point clouds or to surfaces with rapidly varying curvature would require additional machinery, and those settings are outside the paper's claim.

Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Desk Editor's Note

read the letter

Here's the short version: the paper does something genuinely new — it takes Ying-Wang's cut-point surface discretization, built for boundary integral equations, and turns it into a general finite-difference toolbox for surface PDEs. The Laplace–Beltrami, gradient, divergence, and shallow-water formulations are derived cleanly in three coordinate patches, and the numerical evidence is consistently second order across diffusion, Poisson, advection, and the shallow-water equations. The computations look honest: no fitted parameters in the central claim, Taylor expansions are given for the discrete operators, and the convergence tables are detailed across multiple surfaces including a nonconvex Cassini oval. The author explicitly states that a convergence proof 'seems difficult' and gives a resolvent-positivity proof only for curves in R2. That is the soft spot, and it is real: the equilibration step uses quadratic interpolation at secondary cut points, and a naive error count suggests the O(h^3) interpolation error could enter the discrete Laplacian as O(h) after division by h^2. The paper does not isolate and measure this contribution directly, so the observed second-order accuracy relies on cancellations that are not analyzed. This is a genuine gap, but not a fatal one: the empirical convergence is consistent across many runs, and the method is simple enough that others can and probably will test it. The shallow-water max-norm data are less clean — for ν=.5, the max error in Φv stops decreasing between N=320 and N=640 — but the L2 errors continue to improve, and the authors themselves present this as a feasibility study, not a mature scheme. Who is this for? Anyone working on surface PDEs who wants an alternative to closest-point, level-set extension, or point-cloud methods, especially if they already use Cartesian grids and want something simpler. It is not a finished theory paper — the missing convergence proof matters — but it is a solid, clearly written methods paper that deserves a serious referee and, I think, publication after revision. My recommendation: send it to review. It will be read, used, and probably extended.

Referee Report

0 major / 4 minor

Summary. The paper presents a finite-difference method for solving PDEs on smooth closed surfaces, built on the Ying--Wang cut-point representation of a level-set surface. Unknowns are carried at primary cut points, values at secondary cut points are obtained by solving a quadratic-interpolation equilibration system, and differential operators are written in each of the three coordinate planes after projection. The method is applied to surface diffusion and Poisson problems with the Laplace--Beltrami operator, to a linear advection equation with an exact solution on the sphere, and to the shallow water equations on a sphere using a Lax--Wendroff-type scheme with artificial viscosity. The paper reports second-order convergence in extensive numerical tables, computes eigenvalues of the discrete spherical Laplacian, proves a resolvent-positivity result for the analogous discretization on curves in two dimensions, and explicitly states that a full convergence proof for surfaces is not yet available.

Significance. If the observed second-order accuracy holds in the advertised generality, the method is a simple and attractive alternative to triangulated, closest-point, or RBF-based surface PDE solvers for smooth closed surfaces given by level sets: it reduces surface computations to standard finite differences on planar Cartesian grids without boundary conditions. The central discretization is parameter-free in the sense that the finite-difference coefficients are obtained from Taylor expansions rather than from data fitting, and the paper provides unusually complete numerical evidence, including convergence tables for diffusion, advection, and shallow water, eigenvalue tests for the discrete Laplace--Beltrami operator on the sphere, and a high-order quadrature rule. The author's explicit acknowledgment of the missing convergence proof and the candid reporting of less regular shallow-water maximum-norm results are also strengths of the presentation. The main reservation is that the abstract's phrasing of the accuracy claim is broader than what the shallow-water tables support in the maximum norm.

minor comments (4)

[Sections 3 and 5, Eqs. (2.1)-(2.4), (3.23)-(3.25)] The paper states in Section 5 that proving convergence 'seems difficult,' and the consistency analysis in Section 3 is performed for the primary stencil assuming exact values at all stencil points. Since the O(h^3) quadratic-interpolation error at secondary points could in principle enter a second-difference stencil as O(h), it would be informative to include a direct numerical consistency test of the composed operator Delta_h E_h (for example, apply Delta_h E_h to a known smooth function on a nontrivial surface and tabulate the discrete Laplace--Beltrami error as a function of h). This would isolate the equilibration contribution from the primary stencil truncation and would further substantiate the observed second-order convergence.
[Abstract, Section 4, Tables 4.2 and 4.3] The abstract's claim of 'second order accuracy in examples with a variety of equations, including ... shallow water equations' is broader than the reported data. In Table 4.2 with nu=1, the maximum-norm errors for Phi v do not decrease between N=320 and N=640 (e.g., day 1 values 2.26e-3 and 2.17e-3), and the nu=0.5 results in Table 4.3 are non-monotone, with the day-1 maximum error for Phi v increasing from 4.48e-3 at N=320 to 1.37e-2 at N=640. The L2 errors do show roughly second-order behavior, so I suggest qualifying the accuracy statement for the shallow water example as L2 accuracy and restating the feasibility-only goal for nonlinear hyperbolic problems.
[Table 3.2, Cassini oval row] In the BDF2 L2 column for N=320, the entry printed as '2.33-5' should presumably be '2.33e-5'; please correct the typographical formatting.
[Section 4, Eqs. (4.33)-(4.34)] The artificial-viscosity parameter nu takes the values 0.5 and 1 in the shallow water experiments, and the convergence behavior depends noticeable on this choice. A brief sentence on how these values were selected, and on the sensitivity of the reported rates to nu, would improve reproducibility.

Circularity Check

0 steps flagged · score 0.0 of 10

No significant circularity: the discrete operators are derived from Taylor expansions against external benchmarks, and the only self-citations support auxiliary quadrature, not the central claim.

full rationale

The paper's central derivation is self-contained. The discrete Laplace-Beltrami operator is constructed by Taylor-expanding the metric coefficients and centered differences, with the consistency calculation shown explicitly in Eqs. (3.2)-(3.12): the paper derives Δ_h u = Δu + O(h^2) directly from Taylor expansions, not from the numerical results. The equilibration map E_h of Eq. (2.4) is an interpolation-defined linear operator, not a fit to target data. The claimed second-order accuracy is presented as observed convergence in Tables 3.1-3.2, 4.1-4.3, against exact solutions (spherical harmonic diffusion, the analytic advection solution (4.5), and the Williamson test [29]); no fitted parameter is renamed as a prediction. The only self-citations are to [2] and [31] for the surface quadrature rule used in integral diagnostics, which is a supporting tool and not the load-bearing PDE discretization; the surface representation itself is credited to the external reference [32]. The paper explicitly states that a full convergence proof is missing (Section 5) and that the resolvent-positivity proof is only for curves in R^2 (Appendix A); these are admitted limitations and correctness risks, not imported conclusions. I find no equation, parameter, or theorem whose derivation reduces to its own input, so the circularity score is 0.

Assumptions & free parameters 3 free parameters · 5 assumptions · 0 invented entities

The central numerical claim rests on standard smoothness assumptions and on an empirical stability assumption. The only hand-set parameters are discretization thresholds and viscosity, none of which are fitted to the reported errors.

free parameters (3)

η (admissible cut-point threshold) = 0.45
Cut points with |nν| ≥ η are retained; η must be <1/√3. Chosen for the examples, not fitted to improve accuracy.
θ (quadrature partition-of-unity angle) = 62.5 degrees
Used only in the surface-integral quadrature (Appendix B); paired with η=.45.
ν (artificial viscosity) = 0.5 or 1
Added to the shallow-water Lax-Wendroff scheme for stability (Section 4); not optimized.

assumptions (5)

domain assumption The surface is C^2 and is given as a level set φ=0 with |∇φ|≥c0 on the surface.
Stated at the start of Section 2; required for cut-point construction and Lemmas 2.1-2.2.
domain assumption Grid spacing h is small enough that each grid interval contains at most one admissible cut point and primary points are quasi-uniform.
Invoked in Lemma 2.2 and used throughout the derivations.
standard math Standard Taylor expansions for smooth functions on the surface patches give O(h^2) truncation of the discrete operator.
The derivation in Section 3 assumes u and metric coefficients are smooth enough.
standard math The equilibration matrix I - Π_ss is nonsingular and the solve (2.3)-(2.4) is well-defined.
Proven via strict diagonal dominance from |θ| ≤ 1/2; used to define E_h.
domain assumption The time-stepping methods (forward Euler, BDF2, MacCormack) are stable for the chosen k and h.
Not proved; observed numerically in Sections 3 and 4.

how reviews work

0 comments

Cite this review

Pith. "Pith review of Solving Partial Differential Equations on Closed Surfaces with Planar Cartesian Grids." pith.science (2026). https://pith.science/paper/VSGA4EHX

@misc{pith2026190801796,
  author       = {Pith},
  title        = {Pith review of: Solving Partial Differential Equations on Closed Surfaces with Planar Cartesian Grids},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/VSGA4EHX}},
  note         = {Machine review of arXiv:1908.01796}
}

read the original abstract

We present a general purpose method for solving partial differential equations on a closed surface, based on a technique for discretizing the surface introduced by Wenjun Ying and Wei-Cheng Wang [J. Comput. Phys. 252 (2013), pp. 606-624] which uses projections on coordinate planes. Assuming it is given as a level set, the surface is represented by a set of points at which it intersects the intervals between grid points in a three-dimensional grid. They are designated as primary or secondary. Discrete functions on the surface have independent values at primary points, with values at secondary points determined by an equilibration process. Each primary point and its neighbors have projections to regular grid points in a coordinate plane where the equilibration is done and finite differences are computed. The solution of a p.d.e. can be reduced to standard methods on Cartesian grids in the coordinate planes, with the equilibration allowing seamless transition from one system to another. We observe second order accuracy in examples with a variety of equations, including surface diffusion determined by the Laplace-Beltrami operator and the shallow water equations on a sphere.

Figures

Figures reproduced from arXiv: 1908.01796 by the authors.

Discussion (0). Continue with ORCID to comment.

Reference graph

Works this paper leans on

32 extracted references · 31 canonical work pages

[1]

Aris, Vectors, Tensors, and the Basic Equations of Fluid Mechanic s, Dover, New York, 1962

R. Aris, Vectors, Tensors, and the Basic Equations of Fluid Mechanic s, Dover, New York, 1962

work page 1962
[2]

J. T. Beale, W.-J. Ying, and J. R. Wilson, A simple method for computing singular or nearly singular integrals on closed surfaces , Commun. Comput. Phys. 20 (2016), 733-53

work page 2016
[3]

M. J. Berger, D. A. Calhoun, C. Helzel and R. J. LeVeque, Logically rectangular ﬁnite volume methods with adaptive reﬁnement on the sphere , Phil. Trans. R. Soc. A 367 (2009), 4483-4496

work page 2009
[4]

Bertalmio, L.-T

M. Bertalmio, L.-T. Cheng, S. J. Osher and G. Sapiro, Variational problems and partial diﬀerential equations on implicit surfaces , J. Comput. Phys. 174 (2001), 759-780

work page 2001
[5]

Bouchon, Monotonicity of some perturbations of irreducibly diagona lly dominant M- matrices, Numer

F. Bouchon, Monotonicity of some perturbations of irreducibly diagona lly dominant M- matrices, Numer. Math. 105 (2007), 591-601

work page 2007
[6]

J. H. Bramble and B. E. Hubbard, On a ﬁnite diﬀerence analogue of an elliptic boundary problem which is neither diagonally dominant nor of non-neg ative type, J. Math. Phys. 43 (1964), 117-132

work page 1964
[7]

D. A. Calhoun, C. Helzel and R. J. LeVeque, Logically rectangular grids and ﬁnite volume methods for PDEs in circular and spherical domains , SIAM Rev. 50 (2008), 723–52. 19

work page 2008
[8]

Chen and C

Y. Chen and C. B. Macdonald, The closest point method and multigrid solvers for elliptic equations on surfaces , SIAM J. Sci. Comput. 37 (2015), A134–A155

work page 2015

Show all 32 references

[9]

Chu and R

J. Chu and R. Tsai, Volumetric variational principles for a class of partial di ﬀerential equations deﬁned on surfaces and curves , Res. Math. Sci. 5 (2018), Paper No. 19

2018
[10]

M. P. DoCarmo, Diﬀerential Geometry of Curves and Surfaces , second ed., Dover, New York, 2016

2016
[11]

, D. R. Durran, Numerical Methods for Fluid Dynamics With Applications to G eo- physics, Second Ed., Springer, New York, 2010

2010
[12]

Dziuk and C

G. Dziuk and C. M. Elliott, Finite element methods for surface PDEs , Acta Numerica 22 (2013), 289-396

2013
[13]

Fornberg and N

, B. Fornberg and N. Flyer, A Primer on Radial Basis Functions with Applications to the Geosciences, SIAM, Philadelphia, 2015

2015
[14]

F. X. Giraldo, J. S. Hesthaven, and T. Warburton, Nodal high-order discontinuous Galerkin methods for the spherical shallow water equations , J. Comput. Phys. 181 (2002), 499-525

2002
[15]

Hackbusch, Elliptic Diﬀerential Equations: Theory and Numerical Trea tment, sec- ond ed., Springer, Berlin, 2017

W. Hackbusch, Elliptic Diﬀerential Equations: Theory and Numerical Trea tment, sec- ond ed., Springer, Berlin, 2017

2017
[16]

Leung, J

S.-Y. Leung, J. Lowengrub, and H.-K. Zhao, A grid based particle method for solving partial diﬀerential equations on evolving surfaces and mod eling high order geometrical motion, J. Comput. Phys. 230 (2011), 2540-61

2011
[17]

Liang and H.-K

J. Liang and H.-K. Zhao, Solving partial diﬀerential equations on point clouds , SIAM J. Sci. Comput. 35 (2013), A1461-86

2013
[18]

C. B. Macdonald and S. J. Ruuth, The implicit closest point method for the numerical solution of partial diﬀerential equations on surfaces , SIAM J. Sci. Comput. 31 (2009), 4330–4350

2009
[19]

C. B. Macdonald, B. Merriman and S. J. Ruuth, Simple computation of reaction- diﬀusion processes on point clouds , Proc. Natl. Acad. Sci. USA 110 (2013), 9209-14

2013
[20]

Petras and S

A. Petras and S. J. Ruuth, PDEs on moving surfaces via the closest point method and a modiﬁed grid based particle method , J. Comput. Phys. 312 (2016), 139–56

2016
[21]

Peyret and T

R. Peyret and T. D. Taylor, Computational Methods for Fluid Flow , Springer, New York, 1983

1983
[22]

S. J. Ruuth and B. Merriman, A simple embedding method for solving partial diﬀerential equations on surfaces , J. Comput. Phys. 227 (2008), 1943-61

2008
[23]

Schwartz, D

P. Schwartz, D. Adalsteinsson, P. Colella, A. P. Arkin, and M. On sum, Numerical computation of diﬀusion on a surface , Proc. Natl. Acad. Sci. USA 102 (2005), 11151-56. 20

2005
[24]

Shankar, G

V. Shankar, G. B. Wright, R. M. Kirby and A. L. Fogelson, A radial basis function (RBF)-ﬁnite diﬀerence (FD) method for diﬀusion and reactio n-diﬀusion equations on surfaces, J. Sci. Comput. 63 (2016), 745–768

2016
[25]

Staniforth and J

A. Staniforth and J. Thuburn, Horizontal grids for global weather and climate prediction models: a review , Q. J. R. Meteorol. Soc. 138 (2012), 1–26

2012
[26]

P. N. Swarztrauber, D. L. Williamson, and J. B. Drake, The Cartesian method for solv- ing partial diﬀerential equations in spherical geometry , Dyn. Atmos. Oceans 27 (1997), 679-706

1997
[27]

P. A. Ullrich, C. Jablonowski and B. van Leer, High-order ﬁnite-volume methods for the shallow-water equations on the sphere , J. Comput. Phys. 229 (2010), 6104-34

2010
[28]

Wang, S.-Y

M. Wang, S.-Y. Leung and H.-K. Zhao, Modiﬁed virtual grid diﬀerence for discretizing the Laplace-Beltrami operator on point clouds , SIAM J. Sci. Comput. 40 (2018), A1– A21

2018
[29]

D. L. Williamson, J. B. Drake, J. J. Hack, R. Jakob and P. N. Swar ztrauber, A stan- dard test set for numerical approximations to the shallow wa ter equations in spherical geometry, J. Comput. Phys. 102 (1992), 211–24

1992
[30]

D. L. Williamson, The evolution of dynamical cores for global atmospheric mod els, J. Meteorol. Soc. Jpn., 85B (2007), 241–269

2007
[31]

J. R. Wilson, On computing smooth, singular and nearly singular integral s on implicitly deﬁned surfaces, Ph.D. thesis, Duke University (2010), http://search.proquest.com/docview/744476497

2010
[32]

Ying and W.-C

W.-J. Ying and W.-C. Wang, A kernel-free boundary integral method for implicitly deﬁned surfaces, J. Comput. Phys. 252 (2013), 606-624. 21

2013

Pith tools

Reviewed August 14, 2026 · model on record in the stance chip above.

[1] [1]

Aris, Vectors, Tensors, and the Basic Equations of Fluid Mechanic s, Dover, New York, 1962

R. Aris, Vectors, Tensors, and the Basic Equations of Fluid Mechanic s, Dover, New York, 1962

work page 1962

[2] [2]

J. T. Beale, W.-J. Ying, and J. R. Wilson, A simple method for computing singular or nearly singular integrals on closed surfaces , Commun. Comput. Phys. 20 (2016), 733-53

work page 2016

[3] [3]

M. J. Berger, D. A. Calhoun, C. Helzel and R. J. LeVeque, Logically rectangular ﬁnite volume methods with adaptive reﬁnement on the sphere , Phil. Trans. R. Soc. A 367 (2009), 4483-4496

work page 2009

[4] [4]

Bertalmio, L.-T

M. Bertalmio, L.-T. Cheng, S. J. Osher and G. Sapiro, Variational problems and partial diﬀerential equations on implicit surfaces , J. Comput. Phys. 174 (2001), 759-780

work page 2001

[5] [5]

Bouchon, Monotonicity of some perturbations of irreducibly diagona lly dominant M- matrices, Numer

F. Bouchon, Monotonicity of some perturbations of irreducibly diagona lly dominant M- matrices, Numer. Math. 105 (2007), 591-601

work page 2007

[6] [6]

J. H. Bramble and B. E. Hubbard, On a ﬁnite diﬀerence analogue of an elliptic boundary problem which is neither diagonally dominant nor of non-neg ative type, J. Math. Phys. 43 (1964), 117-132

work page 1964

[7] [7]

D. A. Calhoun, C. Helzel and R. J. LeVeque, Logically rectangular grids and ﬁnite volume methods for PDEs in circular and spherical domains , SIAM Rev. 50 (2008), 723–52. 19

work page 2008

[8] [8]

Chen and C

Y. Chen and C. B. Macdonald, The closest point method and multigrid solvers for elliptic equations on surfaces , SIAM J. Sci. Comput. 37 (2015), A134–A155

work page 2015