pith. sign in

arxiv: 2511.21713 · v1 · submitted 2025-11-18 · ⚛️ physics.comp-ph

A Self-Adjusting FEM-BEM Coupling Scheme for the Nonlinear Poisson-Boltzmann Equation

Mauricio Guerrero-Montero , Michal Bosy , Christopher D. Cooper This is my paper

Pith reviewed 2026-05-17 20:48 UTC · model grok-4.3

classification ⚛️ physics.comp-ph
keywords nonlinear Poisson-Boltzmann equationFEM-BEM couplingrelaxation parameter optimizationNewton-Raphson methodmolecular electrostaticsRNA structures
0
0 comments X

The pith

A coupled FEM-BEM scheme for the nonlinear Poisson-Boltzmann equation automatically selects its own optimal relaxation parameter to converge without manual tuning.

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

The nonlinear Poisson-Boltzmann equation models electrostatics around molecules but the sinh term makes convergence difficult in highly charged cases such as nucleic acids. This work replaces the usual trial-and-error search for a relaxation parameter with an automatic optimization step inside a finite-element / boundary-element coupling. The method is tested on a spherical cavity against APBS and then on several RNA structures, where Newton-Raphson combined with a cubic start for the nonlinearity performs best. Automatic choice of the relaxation factor cuts iterations by roughly forty percent and yields a measured 1.37-fold speed-up over the best hand-tuned value on the most charged test molecule.

Core claim

The paper shows that embedding an optimization routine for the relaxation parameter inside a Newton-Raphson iteration on a FEM-BEM discretization lets the solver reach the nonlinear Poisson-Boltzmann solution reliably and faster than any fixed manual choice, without requiring the user to guess the parameter in advance.

What carries the argument

The automatic relaxation-parameter optimizer embedded in the Newton-Raphson loop of the FEM-BEM coupling.

If this is right

  • The nonlinear form of the equation can be used routinely for nucleic acids and other highly charged biomolecules where the linear approximation is inaccurate.
  • Newton-Raphson with a cubic approximation on the first iteration becomes the preferred nonlinear solver inside the coupling.
  • Eliminating manual tuning removes a common barrier to adoption of the nonlinear Poisson-Boltzmann model in structural biology workflows.

Where Pith is reading between the lines

These are editorial extensions of the paper, not claims the author makes directly.

  • The same self-adjusting idea could be transferred to other nonlinear elliptic problems that currently rely on user-tuned relaxation or continuation parameters.
  • If the optimizer proves robust on larger proteins or membrane systems, the method would lower the computational cost of including explicit nonlinear solvent effects in drug-design calculations.

Load-bearing premise

The automatic procedure for choosing the relaxation parameter will remain effective and stable for any molecular geometry, charge pattern, or solvent condition that might be encountered.

What would settle it

A new molecular structure with an extreme charge distribution or irregular shape on which the automatic method either diverges or requires more iterations than the best hand-selected relaxation value.

Figures

Figures reproduced from arXiv: 2511.21713 by Christopher D. Cooper, Mauricio Guerrero-Montero, Michal Bosy.

Figure 1
Figure 1. Figure 1: Two-region domain split according to Equation (2). [PITH_FULL_IMAGE:figures/full_fig_p005_1.png] view at source ↗
Figure 2
Figure 2. Figure 2: Domain split into three regions according to Equation (3). [PITH_FULL_IMAGE:figures/full_fig_p006_2.png] view at source ↗
Figure 3
Figure 3. Figure 3: Different surfaces around the molecule. The van der Waals (VdW) and solvent [PITH_FULL_IMAGE:figures/full_fig_p007_3.png] view at source ↗
Figure 4
Figure 4. Figure 4: Cross section of an example FEM mesh (in grey), with the SES identified within [PITH_FULL_IMAGE:figures/full_fig_p016_4.png] view at source ↗
Figure 5
Figure 5. Figure 5: Solvation energy of a sphere for different values of [PITH_FULL_IMAGE:figures/full_fig_p018_5.png] view at source ↗
Figure 6
Figure 6. Figure 6: Function |gd| for the first 4 iterations of Picard (solid lines) and Newton-Raphson (dashed lines) methods on the sphere. Note this is a semilog plot of the absolute value, hence, zeros are located at the cusps, where gd changes sign [PITH_FULL_IMAGE:figures/full_fig_p020_6.png] view at source ↗
Figure 7
Figure 7. Figure 7: Function fd for the first 4 iterations of Picard (solid lines) and Newton-Raphson (dashed lines) methods on the sphere. 20 [PITH_FULL_IMAGE:figures/full_fig_p020_7.png] view at source ↗
Figure 8
Figure 8. Figure 8: Evolution of the norm [PITH_FULL_IMAGE:figures/full_fig_p023_8.png] view at source ↗
Figure 9
Figure 9. Figure 9: Number GMRES iterations in each nonlinear iteration ( [PITH_FULL_IMAGE:figures/full_fig_p023_9.png] view at source ↗
read the original abstract

The Poisson-Boltzmann equation is widely used to model molecular electrostatics; however, it is usually solved in linearised form because the sinh nonlinearity is challenging, limiting its applicability in highly charged systems such as nucleic acids. This work presents a solution method for the nonlinear Poisson-Boltzmann equation based on a coupled finite/boundary element scheme that automatically finds an optimal relaxation parameter, ensuring fast and reliable convergence of the nonlinear solver without user intervention. We validated our solver against APBS for a spherical cavity, and used RNA-based structures to perform a thorough study of the different algorithmic choices, and to test our implementation. We found that the best alternative to solve the Poisson-Boltzmann equation was using a Newton-Raphson method where the nonlinearity was gradually introduced with a cubic approximation in the first iteration. Newton-Raphson was also the best method to find the optimal relaxation factor, reducing the number of iterations by 40%. Including other optimisation techniques, we were able to obtain a 1.37x speed-up with respect to the best hand-picked relaxation factor for 1HC8 (molecule with highest charge in our tests), avoiding any trial-and-error process to find the relaxation factor.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Referee Report

2 major / 2 minor

Summary. The paper presents a coupled finite-element/boundary-element (FEM-BEM) scheme for the nonlinear Poisson-Boltzmann equation that automatically determines an optimal relaxation parameter. The central algorithmic contribution is the use of a Newton-Raphson solver combined with a cubic ramp for gradual introduction of the sinh nonlinearity in the first iteration, together with additional optimization techniques, to achieve reliable convergence without manual tuning of the relaxation factor. Validation is performed against APBS on a spherical cavity, and performance is assessed on RNA-based structures, reporting a 40% reduction in iterations and a 1.37× speedup relative to the best hand-picked relaxation factor for the highest-charge test case (1HC8).

Significance. If the automatic relaxation procedure proves robust, the method removes a practical barrier to using the full nonlinear Poisson-Boltzmann equation for highly charged systems such as nucleic acids, where linearized approximations are currently preferred due to solver instability. The reported empirical gains on RNA structures and the parameter-free character of the relaxation search constitute a concrete advance in computational molecular electrostatics.

major comments (2)
  1. [Abstract and Results (RNA-based structures)] The central claim that the scheme 'automatically finds an optimal relaxation parameter, ensuring fast and reliable convergence ... without user intervention' rests on the assumption that the Newton-Raphson-plus-cubic-ramp procedure (and unspecified additional techniques) remains effective for arbitrary geometries, charge distributions, and solvent conditions. The only supporting evidence is agreement with APBS on a single spherical cavity and iteration counts on RNA structures; no analytic convergence proof, no counter-example search, and no results on qualitatively different systems (e.g., membrane proteins or extreme charge densities) are provided. This leaves the robustness claim load-bearing yet only partially substantiated.
  2. [Abstract] The abstract states that 'including other optimisation techniques, we were able to obtain a 1.37x speed-up', yet the manuscript supplies only summary performance numbers without detailed error metrics, convergence plots, or full algorithmic pseudocode for the self-adjusting procedure. This absence makes it difficult to verify reproducibility and to assess whether the reported speedup is attributable to the automatic relaxation search or to other implementation details.
minor comments (2)
  1. [Methods] Notation for the relaxation parameter and the cubic ramp function should be defined explicitly in the methods section rather than introduced only in the abstract.
  2. [Results] The manuscript would benefit from a table summarizing iteration counts, wall-clock times, and final residuals for all tested RNA structures under both automatic and hand-tuned relaxation, to allow direct comparison.

Simulated Author's Rebuttal

2 responses · 2 unresolved

We thank the referee for the constructive feedback and the opportunity to address the concerns raised. We respond point by point to the major comments and indicate the revisions we will incorporate.

read point-by-point responses

  1. Referee: [Abstract and Results (RNA-based structures)] The central claim that the scheme 'automatically finds an optimal relaxation parameter, ensuring fast and reliable convergence ... without user intervention' rests on the assumption that the Newton-Raphson-plus-cubic-ramp procedure (and unspecified additional techniques) remains effective for arbitrary geometries, charge distributions, and solvent conditions. The only supporting evidence is agreement with APBS on a single spherical cavity and iteration counts on RNA structures; no analytic convergence proof, no counter-example search, and no results on qualitatively different systems (e.g., membrane proteins or extreme charge densities) are provided. This leaves the robustness claim load-bearing yet only partially substantiated.

    Authors: We agree that broader validation would strengthen the robustness claim. Our tests deliberately target the spherical cavity for direct comparison with APBS and RNA structures because these represent the highly charged biomolecular cases where the nonlinear Poisson-Boltzmann equation is most relevant and where solver instability currently favors linearized approximations. The Newton-Raphson solver with cubic ramp for gradual nonlinearity introduction, together with the self-adjusting relaxation search, produced reliable convergence without manual tuning across the tested RNA cases (including the highest-charge structure 1HC8). An analytic convergence proof lies outside the scope of this numerical-methods contribution. We will revise the manuscript to add an explicit discussion of the method's intended scope and acknowledged limitations for other geometries such as membrane proteins. revision: partial

  2. Referee: [Abstract] The abstract states that 'including other optimisation techniques, we were able to obtain a 1.37x speed-up', yet the manuscript supplies only summary performance numbers without detailed error metrics, convergence plots, or full algorithmic pseudocode for the self-adjusting procedure. This absence makes it difficult to verify reproducibility and to assess whether the reported speedup is attributable to the automatic relaxation search or to other implementation details.

    Authors: We accept that the current presentation of performance results is insufficient for full reproducibility. The reported 1.37× speedup for 1HC8 combines the automatic relaxation search (via Newton-Raphson) with the cubic-ramp initialization and other solver optimizations described in the methods. In the revised manuscript we will insert the complete pseudocode for the self-adjusting relaxation procedure, add convergence plots for the RNA test cases, and include tabulated error metrics (maximum potential difference and iteration counts) relative to both APBS and the best hand-tuned relaxation factor. These additions will make clear the contribution of the automatic parameter selection. revision: yes

standing simulated objections not resolved
  • An analytic convergence proof for the Newton-Raphson-plus-cubic-ramp scheme across arbitrary geometries and charge distributions
  • Empirical results on qualitatively different systems such as membrane proteins or extreme charge densities

Circularity Check

0 steps flagged

Numerical solver with external validation; no circular reduction in algorithmic claims

full rationale

The paper presents a computational algorithm for solving the nonlinear Poisson-Boltzmann equation via FEM-BEM coupling, with an automatic relaxation parameter found via Newton-Raphson optimization plus cubic ramping during the solve itself. Performance is demonstrated by direct comparison to APBS on a spherical test case and by iteration counts/speedups versus hand-picked factors on RNA structures. No step in the method description or results reduces a claimed prediction or uniqueness result to a fitted input or self-citation by construction; the optimization procedure is independent of the final electrostatic solution and the validation uses external benchmarks.

Axiom & Free-Parameter Ledger

1 free parameters · 1 axioms · 0 invented entities

The paper introduces no new physical axioms or entities; the contribution is algorithmic. The only notable free parameter is the relaxation factor, which the method is designed to determine automatically rather than fix by hand.

free parameters (1)
  • relaxation parameter
    Central to the nonlinear solver convergence; the paper claims an automatic procedure removes the need to choose it manually.
axioms (1)
  • domain assumption Iterative relaxation methods can be applied to the nonlinear Poisson-Boltzmann equation with guaranteed convergence under suitable parameter choice.
    Standard assumption underlying all relaxation-based solvers for the nonlinear PBE.

pith-pipeline@v0.9.0 · 5519 in / 1288 out tokens · 35876 ms · 2026-05-17T20:48:47.368965+00:00 · methodology

discussion (0)

Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.

Reference graph

Works this paper leans on

52 extracted references · 52 canonical work pages

  1. [1]

    Shen and F

    J. Shen and F. A. Quiocho, Journal of Computational Chemistry16, 445 (1995)

    work page 1995

  2. [2]

    N. A. Baker, Methods in Enzymology383, 94 (2004)

    work page 2004

  3. [3]

    H. Fu, H. Chen, M. Blazhynska, E. Goulard Coderc de Lacam, F. Szczepaniak, A. Pavlova, X. Shao, J. C. Gumbart, F. Dehez, B. Roux, et al., Nature Protocols17, 1114 (2022)

    work page 2022

  4. [4]

    M. K. Gilson, A. Rashin, R. Fine, and B. Honig, Journal of Molecular Biology184, 503 (1985)

    work page 1985

  5. [5]

    N. A. Baker, D. Sept, M. J. Holst, and J. A. McCammon, Proceedings of the National Academy of Sciences of the USA98, 10037 (2001)

    work page 2001

  6. [6]

    Rocchia, E

    W. Rocchia, E. Alexov, and B. Honig, The Journal of Physical Chem- istry B105, 6507 (2001)

    work page 2001

  7. [7]

    A. H. Boschitsch and M. O. Fenley, Journal of Chemical Theory and Computation7, 1524 (2011)

    work page 2011

  8. [8]

    Jurrus, D

    E. Jurrus, D. Engel, K. Star, K. Monson, J. Brandi, L. E. Felberg, D. H. Brookes, L. Wilson, J. Chen, K. Liles, et al., Protein Science27, 112 (2018)

    work page 2018

  9. [9]

    C. M. Cortis and R. A. Friesner, Journal of Computational Chemistry 18, 1591 (1997)

    work page 1997

  10. [10]

    L. Chen, M. J. Holst, and J. Xu, SIAM journal on numerical analysis 45, 2298 (2007)

    work page 2007

  11. [11]

    Xie and S

    D. Xie and S. Zhou, BIT Numerical Mathematics47, 853 (2007). 39

    work page 2007

  12. [12]

    S. D. Bond, J. H. Chaudhry, E. C. Cyr, and L. N. Olson, Journal of Computational Chemistry31, 1625 (2010)

    work page 2010

  13. [13]

    Shaw, Physical Review A32, 2476 (1985)

    P. Shaw, Physical Review A32, 2476 (1985)

    work page 1985

  14. [14]

    B. J. Yoon and A. Lenhoff, Journal of Computational Chemistry11, 1080 (1990)

    work page 1990

  15. [15]

    Juffer, E

    A. Juffer, E. F. Botta, B. A. van Keulen, A. van der Ploeg, and H. J. Berendsen, Journal of Computational Physics97, 144 (1991)

    work page 1991

  16. [16]

    A. H. Boschitsch, M. O. Fenley, and H.-X. Zhou, The Journal of Physical Chemistry B106, 2741 (2002)

    work page 2002

  17. [17]

    B. Lu, X. Cheng, J. Huang, and J. A. McCammon, Proceedings of the National Academy of Sciences of the USA103, 19314 (2006)

    work page 2006

  18. [18]

    Geng and R

    W. Geng and R. Krasny, Journal of Computational Physics247, 62 (2013)

    work page 2013

  19. [19]

    C. D. Cooper, J. P. Bardhan, and L. A. Barba, Computer physics com- munications185, 720 (2014)

    work page 2014

  20. [20]

    S. D. Search, C. D. Cooper, and E. Van’t Wout, Journal of Computa- tional Chemistry43, 674 (2022)

    work page 2022

  21. [21]

    L. E. Felberg, D. H. Brookes, E.-H. Yap, E. Jurrus, N. A. Baker, and T. Head-Gordon, Journal of computational chemistry38, 1275 (2017)

    work page 2017

  22. [22]

    S. V. Siryk and W. Rocchia, The Journal of Physical Chemistry B126, 10400 (2022)

    work page 2022

  23. [23]

    A. Jha, M. Nottoli, A. Mikhalev, C. Quan, and B. Stamm, The Journal of Chemical Physics158(2023)

    work page 2023

  24. [24]

    Jha and B

    A. Jha and B. Stamm, arXiv preprint arXiv:2309.06862 (2023)

    work page arXiv 2023

  25. [25]

    A. H. Boschitsch and M. O. Fenley, Journal of Computational Chemistry 25, 935 (2004)

    work page 2004

  26. [26]

    Ying and D

    J. Ying and D. Xie, Applied Mathematical Modelling58, 166 (2018). 40

    work page 2018

  27. [27]

    M. Bosy, M. W. Scroggs, T. Betcke, E. Burman, and C. D. Cooper, Journal of Computational Chemistry45, 787 (2024)

    work page 2024

  28. [28]

    M. D. Altman, J. P. Bardhan, J. K. White, and B. Tidor, Journal of Computational Chemistry30, 132 (2009)

    work page 2009

  29. [29]

    Nicholls and B

    A. Nicholls and B. Honig, Journal of computational chemistry12, 435 (1991)

    work page 1991

  30. [30]

    C. Li, M. McGowan, E. Alexov, and S. Zhao, Mathematical biosciences and engineering: MBE17, 6259 (2020)

    work page 2020

  31. [31]

    M. J. Holst and F. Saied, Journal of computational chemistry16, 337 (1995)

    work page 1995

  32. [32]

    Cai, M.-J

    Q. Cai, M.-J. Hsieh, J. Wang, and R. Luo, Journal of Chemical Theory and Computation6, 203 (2010)

    work page 2010

  33. [33]

    M. L. Connolly, Journal of applied crystallography16, 548 (1983)

    work page 1983

  34. [34]

    Johnson, J

    C. Johnson, J. Nédélec, et al., Mathematics of Computation35, 1063 (1980)

    work page 1980

  35. [35]

    Guerrero-Montero, M

    M. Guerrero-Montero, M. Bosy, and C. D. Cooper, Journal of Compu- tational Chemistry46, e70036 (2025)

    work page 2025

  36. [36]

    O. Steinbach,Numerical approximation methods for elliptic boundary value problems(Springer, New York, 2008), ISBN 978-0-387-31312-2, fi- nite and boundary elements, Translated from the 2003 German original

    work page 2008

  37. [37]

    B. J. Yoon and A. M. Lenhoff, Journal of Computational Chemistry11, 1080 (1990)

    work page 1990

  38. [38]

    M. K. Gilson and B. Honig, Proteins: Structure, Function, and Bioin- formatics4, 7 (1988)

    work page 1988

  39. [39]

    Betcke and M

    T. Betcke and M. Scroggs, Journal of Open Source Software6, 2879 (2021)

    work page 2021

  40. [40]

    Logg and G

    A. Logg and G. N. Wells, ACM Transactions on Mathematical Software (TOMS)37, 1 (2010). 41

    work page 2010

  41. [41]

    et al.,trimesh, URLhttps://trimsh.org/

    D.-H. et al.,trimesh, URLhttps://trimsh.org/

    work page

  42. [42]

    M. F. Sanner, A. J. Olson, and J.-C. Spehner, Biopolymers38, 305 (1996)

    work page 1996

  43. [43]

    Decherchi and W

    S. Decherchi and W. Rocchia, PloS one8, e59744 (2013)

    work page 2013

  44. [44]

    C. T. Lee, J. G. Laughlin, J. B. Moody, R. E. Amaro, J. A. McCammon, M. Holst, and P. Rangamani, Biophysical Journal118, 1003 (2020)

    work page 2020

  45. [45]

    T. Wang, C. D. Cooper, T. Betcke, and L. A. Barba, arXiv preprint arXiv:2103.01048 (2021)

    work page arXiv 2021

  46. [46]

    J. S. Kieft and I. Tinoco, Structure5, 713 (1997)

    work page 1997

  47. [47]

    J. W. Ponder and D. A. Case, Advances in Protein Chemistry66, 27 (2003)

    work page 2003

  48. [48]

    T.J.Dolinsky, J.E.Nielsen, J.A.McCammon, andN.A.Baker, Nucleic acids research32, W665 (2004)

    work page 2004

  49. [49]

    Dock-Bregeon, B

    A. Dock-Bregeon, B. Chevrier, A. Podjarny, J. Johnson, J. De Bear, G. Gough, P. Gilham, and D. Moras, Journal of molecular biology209, 459 (1989)

    work page 1989

  50. [50]

    Westhof and M

    E. Westhof and M. Sundaralingam, Biochemistry25, 4868 (1986)

    work page 1986

  51. [51]

    L. Guo, J. Ding, R. Guo, Y. Hou, D.-C. Wang, and L. Huang, Journal of Biological Chemistry289, 1478 (2014)

    work page 2014

  52. [52]

    G. L. Conn, A. G. Gittis, E. E. Lattman, V. K. Misra, and D. E. Draper, Journal of molecular biology318, 963 (2002). 42

    work page 2002