1 The method of random walk on spheres for solving boundary-value problems from molecular electrostatics Nikolai A. Simonov Institute of Computational Mathematics and Mathematical Geophysics SB RAS, Novosibirsk, Russia and School of Computational Science, Florida State University, Tallahassee, FL, 32306, USA, E-mail: simonov@csit.fsu.edu Michael Mascagni Department of Computer Science and School of Computational Science, Florida State University, Tallahassee, FL, 32306, USA, E-mail: mascagni@fsu.edu Abstract - We consider the boundary-value problem for el- liptic equation in the whole space with continuity boundary conditions on the surface of a compact body and describe the advanced walk-on-spheres algorithm for solving this problem. Keywords— Monte Carlo, walk on spheres, random walk, Poisson-Boltzmann, molecule, electrostatics, point value, free energy I NTRODUCTION The method of random-walk-on-spheres is an efficient Monte Carlo algorithm that makes it possible to calculate values of the solution of an elliptic partial differential equa- tion equation and its derivatives at given points without the need to solve the boundary-value problem in the whole domain [1], [2]. This method is well established for the Dirichlet problem, whereas problems involving other types of boundary conditions that utilize the normal derivative of the solution are still the subject of current research. Here, we consider one of the topical problems of electrostatics that involves such boundary conditions, and propose an advanced approach to constructing an appropriate Monte Carlo algorithm. The article is organized as follows. The first section is devoted to the description of the non-classical boundary- value problem we consider. In the second part we give a brief review of random walk algorithms that we use here and ways of their efficient computational implementation. In the third section we present the Monte Carlo estimate for point value of the solution. This estimate is based on the approach described in section 3. The fourth part of the article is devoted to the description of the new way of treat- ing flux-type boundary conditions. Finally, we present the results of some model computational experiments and dis- cuss further directions of our future research. I. STATEMENT OF THE PROBLEM Let u(x) satisfy the following elliptic partial differential equations inside bounded domain, G, and in its exterior: ǫ i Δu(x)+4πρ(x)=0 ,x ∈ G, (1) ǫ e Δu(x) - ǫ e κ 2 u(x)=0 ,x ∈ G 1 ≡ R 3 \ G, (2) where ǫ e ≥ ǫ i > 0, ρ(x) is a fixed function (the static charge density), and κ 2 (the inverse Debye length) is a pos- itive constant. These equations are coupled by the continu- ity conditions on the boundary of the domain: u i (y)= u e (y) ,ǫ i ∂u i ∂n(y) = ǫ e ∂u e ∂n(y) ,y ∈ ∂G. (3) Here, for convenience, u i denotes the solution to (1) in the interior of G, and u e is the solution of (2) in G 1 . We also assume that u e (x) → 0 as |x| goes to infinity; a common assumption in electrostatics. Boundary-value problems of this kind arise in electrostat- ics problems for dielectric media. In such settings, the so- lution, u(x), gives the electrostatic potential for the system of a dielectric body, G, immersed into another dielectric. In molecular biophysics applications, G can be thought of as a molecule in aqueous solution. In the framework of the widely used implicit solvent model, water and ions dis- solved within are treated as a continuous medium, whose properties are characterized by the dielectric permittivity, ǫ e , whereas the molecule under investigation is described explicitly. The solute (large molecule) is thought of as a cavity with dielectric constant, ǫ i , which is usually much less than that of the exterior environment. To describe the geometry of this cavity different models can be applied [3]. One of the simplest and most com- monly used in electrostatics computations is the model in which the molecule in question is thought of as a union of a large number of intersecting spheres (atoms): G =  M m=1 B(x (m) ,r (m) ). In other models, the molecular sur- face is defined as the contact surface formed between the boundary of G (so-called van der Walls envelope) and a solvent molecular probe (sphere) of appropriate radius. While in Chemistry, one learns that electrical charge is quantized in integer multiples of the charge of the elec- tron, in molecular electrostatics it is common to assign to every spherical atom, with radius r (m) , a partial electrical charge, q m , which is positioned at its center, x (m) . This takes into account the partial effective charges of the large molecule that occur due to electrostatic screening, bond polarization, and other common chemical phenomena re- lated to electrostatics. Depending on the model, some ad-