Monte Carlo Methods for Calculating Some Physical Properties of Large Molecules

Abstract

Abstract. In this paper we describe Monte Carlo methods for solving some boundary-value problems for elliptic partial differential equations arising in the computation of physical properties of large molecules. The constructed algorithms are based on walk on spheres, Green's function first passage, walk in subdomains techniques, and finite-difference approximations of the boundary condition. The methods are applied to calculating the diffusion-limited reaction rate, the electrostatic energy of a molecule, and point values of an electrostatic field.Key words. Monte Carlo method, random walk, diffusion, reaction rate, electrostatic energy, molecule AMS subject classifications. 65C05, 65N99, 78M25, 92C451. Introduction. Elliptic partial differential equations such as the Laplace, Poisson, Poisson-Boltzmann, etc. are effectively used as mathematical models in different branches of computational biophysics and chemistry. Calculation of the diffusionlimited reaction rate, the electrostatic potential and field, the internal energy -are all problems that can be reduced to the solution of a diffusion equation with some conditions on the boundary and at the infinity. The intrinsic analogy between diffusion and electrostatics makes it possible to apply the same computational techniques to solve problems coming from these different fields. It is worth noting that the solution domain in this class of problems is usually infinite. So, it is natural that since the days of Maxwell, the boundary-element method has been used as an effective tool for solving such problems. In particular, one can calculate the capacitance and the diffusionlimited reaction rate as a surface integral [18]. The boundary-element method, as well as finite-difference and finite-element methods, are still commonly used for solving electrostatics and diffusion problems arising in biophysics. Review of these and other techniques used in the computation of molecular electrostatic properties is given in [2].Another possible way of computationally treating these problems comes from the probabilistic representation of solutions to elliptic and parabolic partial differential equations as functionals of diffusion process trajectories [15,7,8]. Direct computational simulation of physical diffusion in this case coincides with the approximation of Brownian motion as the solution to a stochastic differential equation via a firstorder Euler scheme [19,16]. This approach was applied [4] to simulation studies of diffusion-limited reactions. Though computationally far from being optimal, it allows one to include different physical phenomena (hydrodynamic, electrostatic, etc.) into the computational scheme, and, what is essential, it is efficient enough to be competitive with deterministic methods. Later, this algorithm was modified [30] by incorporating different boundary conditions [28,26]. In addition, other variants of the algorithm were suggested. In [31], in particular, the Brownian dynamics simulation method was compared to the algorithms based on survival probab...

Simulating diffusions with piecewise constant coefficients using a kinetic approximation

“…The values of u ε and v ε at the reference points are accurate from 2 to 3 digits. We now study the case a 1 = 80 and a 2 = 4 where the physical constants are taken from a problem on the computation of molecule physical properties [36] and where now ε = 0.0005 using the same scaling than in dimension one. The values of u ε and v ε at the reference points are as accurate as previously and the CPU times similar which indicates that our scaling is efficient.…”

“…This kind of problem arises in a lot of modelling problems, such as molecular electrostatics [36], geophysics [40], ecology [7,39], astrophysics [43], magneto/electro-encephalography (MEG, EEG) [19],...…”

“…An alternative method derived from finite differences. In [36], M. Mascagni and N. Simonov proposed another way to solve a PDE of type (1). When the particle reaches a point x of S that separates two zones with diffusion coefficients a 1 and a 2 , it is reinjected at a point x + ε n(x) with probability p + = a 1 /(a 1 + a 2 ) and x − ε n(x) with probability p − = a 2 /(a 1 + a 2 ), where n(x) is the unit vector normal to S and point toward the direction where the diffusion coefficient is a 1 .…”

An introduction to probabilistic methods with applications

“…Besides the fact that these mathematical results provide a significant contribution to the stochastic analysis of linearized Poisson-Boltzmann equations, the authors also provide a rigorous analysis, with detailed numerical comparisons and refined stochastic algorithms, of the Monte Carlo approximation method proposed by Mascagni and Simonov [29].…”

Random Walks for Solving Boundary-Value Problems with Flux Conditions

“…This algorithm works well for the Dirichlet boundary conditions. Common way of treating flux conditions is to simulate reflection from the boundary in accordance with the finite-difference approximation of normal derivative [5][6][7][8]. Such approach introduces bias into the estimate.…”