Abstract
For solving some elliptic boundary value problems, Monte Carlo diffusion algorithms are often the most efficient ones. Among them, last-passage algorithms are good for obtaining charge density at a specific point on a conducting surface. In our previous research, we developed the last-passage Monte Carlo algorithm for charge density on a flat conducting surface. In the research, we used the Laplace Green’s function only on a flat surface. In this paper, we further develop the last-passage algorithm on a spherical surface. We demonstrate the last-passage algorithm by obtaining charge density on a sphere held at unit potential. In addition, using the last-passage algorithm we compute the mutual capacitance and charge distribution of two conducting spheres. We compare them with the analytic results of J. Lekner to find an excellent agreement.
Similar content being viewed by others
1 Introduction
According to the Feynman-Kac formula of probabilistic potential theory [1, 2], supposing that the function coefficients of the differential operator L and the function, \(c(\vec {x} )\le 0\), are bounded and Lipschitz continuous, the boundary function is bounded and continuous, and the characteristic form is non-negative definite, for the following Dirichlet’s boundary value problem,
the solution in the domain D is
where \(X_s\) and \(X_t\) are Itô diffusion processes corresponding to the differential operator L, \(\tau = \inf \{t:X_t \in \partial D \}\) is the first-passage time and \(X(\tau \)) the first-passage point on the boundary.
Based on the above theorem, for the simple homogeneous and isotropic diffusion process (one of the Itô diffusion processes, that is, the Brownian diffusion process) corresponding to the Laplacian differential operator, the following boundary value problem
with the Dirichlet’s boundary condition,
from Eq. (3), has the solution in the domain D given by
Here, \(\tau = \inf \{t:X_t \in \partial D \}\) is the first-passage time, and \(X(\tau \)) the first-passage point of the homogeneous and isotropic diffusion on the boundary. If we introduce the first-passage probability distribution, \(p(\vec {x},\vec {y})\), where \(\vec {x} \in D\) and \(\vec {y} \in \partial D\), the solution Eq. (6) can be rewritten as an integral equation,
And when we reinterpret the Laplace diffusion problem with the boundary condition, \(\psi (\vec {y})=1\) as an electrostatic problem, that is, the boundary as a conducting surface, we obtain the correspondence, Fig. 1. Besides, in electrostatics using Green’s function \(G(\vec {x},\vec {y})\) the solution is
where \(\vec {n}\) is the normal vector to the boundary and \({\partial G(\vec {x},\vec {y})}/{\partial \vec {n}}\) the charge density. From Eqs. (7) and (8), we can conclude the correspondence between the first-passage probability on the absorbing boundary and the surface charge density on the conducting surface,
This correspondence is the basis for the Monte Carlo diffusion capacitance algorithm [3].
Along with the similar correspondences, the point charges in the electrostatic Poisson problem correspond to the sources of diffusion, the conducting surface to the absorbing boundary, the electrostatic induced charge density on the conducting surface to the first-passage probability density on the boundary and so on [4, 5]. These correspondences are of great value to diffusion Monte Carlo studies of electrostatic potential problems involving complex geometries and/or singular boundary conditions because it allows massively parallel computation and local implementation of boundary conditions. Using the correspondences, diffusion Monte Carlo algorithms have been developed [6].
In diffusion Monte Carlo algorithms, first- and last-passage algorithms have been used for computing charge distribution on a conducting surface [4, 5, 7,8,9]. The first-passage algorithms allow one to decompose a complex diffusion path into geometrically simple diffusion jumps, using the Markovian nature of the diffusion path. To expedite diffusion simulations in the first- and last-passage algorithms, rather than “Walks-on-Spheres” (WOS) [7], other diffusion simulation techniques such as Green’s function based first-passage (GFFP) methods [10, 11], “Walks-on-Planes” (WOP) [12, 13], off-centered WOS [9], and parallel-plates method [14] have been developed.
Here, we briefly summarize a first-passage algorithm, Kai-Lai Chung’s last-passage algorithm and their relation [15] (see Fig. 2). The first-passage algorithm [5, 12] can give us the overall surface charge distribution on a conductor held at unit potential with respect to infinity,
Here, \(F(\infty ;d\vec {y})\) is the first-passage probability distribution of the diffusion starting from infinity to \(d\vec {y}\). In contrast, the Kai-Lai Chung’s last-passage algorithm [16] gives the charge density as
Here, \(L(\vec {x}; d\vec {y})\) is the last-passage density of the diffusion starting from \(\vec {x}\), to \(d\vec {y}\) and going to infinity and \(u(\vec {x}, \vec {y})\) the electrostatic potential difference between \(\vec {x}\) and \(\vec {y}\).
There is a time-reversal symmetry of the Brownian diffusion processes of the first-passage algorithm starting from infinity to \(\vec {x}\) and the Chung’s last-passage algorithm starting from \(\vec {x}\) to infinity (See Fig. 2). The first-passage point of the first-passage algorithm becomes the last-passage point of the Chung’s last-passage algorithm. Equating the two algorithm equations, the last-passage probability distribution of the Chung’s algorithm can be related to the first-passage probability distribution starting from infinity to \(\vec {x}\). That is,
The first-passage algorithm is good for overall charge distribution on the conducting surface while the Chung’s last-passage algorithm can compute charge density around a point on the conducting surface as an importance sampling around the point, and Given-Hwang’s last-passage algorithm [5] can compute charge density at a specific point on the conducting surface. The Given-Hwang’s last-passage algorithm can be considered as a special last-passage algorithm of the Chung’s last-passage algorithm [4]. In the Given-Hwang’s last-passage algorithm, the point \(\vec {x}\) inside the domain in Fig. 2 is on the boundary, that is, the limit case of the Chung’s last-passage algorithm to the boundary. Basically, both last-passage algorithms are the time-reversals of the first-passage algorithm and they are the same kind [4]. The implementation of the Chung’s last-passage algorithm is under investigation now.
At first [5], Given-Hwang’s last-passage algorithm was developed with a half-sphere on a flat surface (see Fig. 3). According to the isomorphism of Fig. 1, in Fig. 3 the electrostatic potential \(V(\vec {x}+\vec {\varepsilon })\) by a point dipole at \(\vec {x}\), very near the point \(\vec {x}\) on the conducting surface held at zero potential by the point dipole, is
Here, \(g(\vec {x} +\vec {\varepsilon },\vec {y})\) is the probability density associated with a Brownian particle (random walker, diffusing particle or diffusion) initiating at the point \(\vec {x}+\vec {\varepsilon }\) and making first-passage on the surface \(\partial \Omega _{\vec {y}}\) at the point \(\vec {y}\) [5], and \(p(\vec {y},\infty )\) the probability associated with the diffusing particle initiating at the point \(\vec {y}\) on the upper first-passage surface and diffusing to infinity. Applying Gauss’ law gives the charge density
where the term
is the Laplacian Green’s function for a point dipole centered on the flat conducting surface. Later [17], the last-passage algorithm was developed and tested for the mutual capacitance of a system of flat conductors.
In this paper, we further develop the Given-Hwang’s last-passage algorithm not on a flat surface but on a spherical surface. We verify the last-passage algorithm by obtaining charge density on a single sphere held at unit potential (see Fig. 4). Then, as a useful demonstration, using the last-passage algorithm on the spherical surface we compute the mutual capacitance and charge distribution of two spheres and compare with the previous research by Lekner [18].
2 Last-Passage Algorithm on Sphere
In this paper, following the previous research [4], heuristically we can obtain the Laplacian Green’s function \(G(\vec {x}, \vec {y})\) for a point dipole centered on the spherical surface at point \(\vec {x}\) and normal to the surface (see Fig. 4). The Green’s function is based on the isomorphism, provided by probabilistic potential theory [1, 2, 5], between electrostatic Dirichlet Poisson equation problem on a conducting surface and the corresponding Brownian diffusion process expectation.
At first, we introduce a spherical coordinate centered at \(\vec {x}\) (see Fig. 4), where the z-axis is normal to the sphere with radius R and \(\theta \) the angle between the line perpendicular to the surface and \(\vec {y}\). Then, the potential \(V({\vec {y}})=V(r,\theta , \phi \)) on the spherical surface \(\partial \Omega _{\vec {y}}\) due to the dipole at \(\vec {x}+\vec {\varepsilon }\) and \(\vec {x} - \vec {\varepsilon } R/(R+\varepsilon )\) is given by
Here, the vector \(\vec {\varepsilon }\) is \(\varepsilon \hat{n}\), where \(\hat{n}\) is the unit vector normal to the conducting surface. Differentiating partially with respect to r, we get
And differentiating again partially with respect to \(\varepsilon \) and setting \(\varepsilon =0\), we obtain
After normalizing and putting \(\varepsilon _0=1\) and \(q=1\), it becomes the Green’s function,
To verify the Green’s function on the spherical surface, we compute the charge density on a sphere held at unit potential (see Fig. 4). Here, we use the same coordinate system as before and the correspondence between the electrostatic problem and the diffusion one. By the correspondence, the probability of going to infinity from \(\vec {y}\) is \(1-R/r'\) because the electrostatic potential at \(\vec {y}\) is \(R/r'\) (see Fig. 1). Here, \(r'\) is the distance between \(\vec {y}\) and the center of the sphere with radius R:
Using (14), the charge density on the spherical conducting surface is
where \(\cos \alpha =-r/(2R)\). Changing the variable a little, using \(\theta =\alpha p\), the density becomes
The analytical charge density of the sphere is \(q/(4\pi R^2)\). Here, the sphere should be at unit potential on its surface, \(V(r')=q/r'\), so \(q=R\). Then, the analytical charge density is \(1/(4\pi R)\). Actually, Eq. (22) is integrable and the integration gives the exact value, \(1/(4 \pi R)\).
In addition, using the Monte Carlo integration we compute the charge density. Comparing with the analytical value \(1/(4\pi R)\), we verify the last-passage algorithm. In Table 1, we find that the evaluated \(\sigma \) coincides with the analytical charge density very well.
3 Mutual Capacitance and Charge Distribution of Two Spheres
In the previous research [17], we computed the mutual capacitance of the circular-disk parallel-plate capacitor. In that case, the geometry of the capacitor is flat.
Here, to demonstrate the last-passage algorithm on a spherical surface, we compute the mutual capacitance of two spheres with radii \(R_1\) and \(R_2\) and separation d between the two spheres (see Fig. 5). The mutual capacitance matrix \(C_{ij}\) of a set of N conductors is defined by the relation [17]
where \(V_i\) and \(Q_i\) are the potential and total charge on the i-th conductor, respectively. \(C_{ij}\) is the total charge on conductor i, when unit potential is applied to conductor j while grounding all the other conductors. Here, \(C_{ii}\) and \(C_{ij}\) are given by
Correspondingly,
where \(P(\vec {y} \rightarrow C^i)\) [\(P(\vec {y} \rightarrow C^j)]\) is the probability that a diffusing particle starts at the point \(\vec {y}\) will be absorbed on the surface \(C^i\) [\(C^j\)] of the i-th [j-th] conductor [17] (for the definition of \(\vec {y}\), you can refer to Fig. 4). In this way, we can compute the total charges on the spherical surfaces using the probabilities of going to infinity or going back to the diffusion-starting sphere or going to the other sphere. This is another correspondence between an electrostatic problem and the diffusion one.
There are analytic expressions for capacitance and charge distribution on the two conducting spheres using bispherical coordinates [18]. The capacitance coefficients are
with
The charge distribution on the two spheres is
with
Here, \(P_n(x)\) is the Legendre polynomial. \(\theta _1\) and \(\theta _2\) are polar angles from the central axis that passes the centers of the two spheres.
In Fig. 6, we demonstrate the Monte Carlo error convergence of \(C_{11}\), as an example, with the epsilon layer, \(\varepsilon = 10^{-6}\), with \(R_1 = 1.0\), \(R_2=1.0\), and \(d=3.0\). The graph shows the typical Monte Carlo error convergence. In Table 2, mutual capacitance matrices computed by the Monte Carlo method in various situations are compared with analytic results. Here, the number of random walks per run is \(10^9\).
In Fig. 7, for the sample case of \(R_1 = 1.0\), \(R_2=1.0\), \(d=10.0\), and \(\varepsilon =10^{-6}\), the charge distribution on the spherical conducting surface is compared with the analytic result. Here, the number of random walks per run is \(10^{11}\) and \(\theta \) the polar angle of one of the spheres. We also observe an excellent agreement with previous analytic results by Lekner [18].
4 Algorithm Summary
Charge density and capacitance of spheres are calculated by (24)–(27) through diffusion Monte Carlo methods. With the configuration of Fig. 5, \(N_{\mathrm {total}}\) particles diffuse with the WOS algorithm from a sampled point on the sphere of radius \(R_1\). The departure from the surface (i.e. \(\vec {x} \rightarrow \vec {y}\) in Fig. 4) is done by the random sampling of \(0\le \phi \le 2\pi \) and \(0\le \theta \le \alpha \) with \(\cos \alpha =-r/(2R)\) and an appropriate radius r. The trajectory continues until it hits one of the two spheres or diffuses away to infinity. \(N_1\) is the number of the trajectories that return to the sphere of radius \(R_1\), \(N_2\) the number of the trajectories that arrive at the other sphere of radius \(R_2\), and \(N_3\) the number of the particles that do not touch the two spheres forever. Charge density in (26) and (27) at a specific point on the sphere is related to the ratios of \(N_1\), \(N_2\), and \(N_3\) with a weighting factor:
with
where \(i=1\), 2, and 3.
For the fast simulation of a diffusion to infinity, a hypothetical sphere is used that encloses both spheres with a sufficiently large radius. Once the diffusing particle gets out of this enclosing sphere, whether the particle goes to infinity or not is tested using the probability described in Fig. 1. If the particle does not go to infinity, it goes back to the surface of the hypothetical sphere and continues to diffuse [9].
To compute mutual capacitances, \(\sigma _{11}\) and \(\sigma _{12}\) over a large enough number of points on the sphere should be obtained. Then the integrations of Eqs. (24) and (25) are reduced to the summations as
where \(\vec {x}_n\) is the point on the sphere of radius \(R_1\) and \(N_{\mathrm {sampling}}\) the total number of sampling points for charge densities. The implementation of this process is summarized as a pseudo-code in Algorithm 1.
5 Conclusion
Compared with finite difference methods or finite element methods, Monte Carlo algorithms can be very efficient in some elliptic boundary value problems, where we have singular points or very complicated geometries; Monte Carlo samplings can avoid the singular points very naturally and Monte Carlo methods do not need any discretization. One kind of Monte Carlo algorithms, first- and last-passage algorithms have been used for obtaining charge density distribution on a conducting surface. Especially first-passage algorithms are used for the overall charge distribution and the total capacitance of the electrostatic conductor, and last-passage algorithms are good for obtaining charge density at a specific point on the conducting surface.
In our previous research [5], we obtained the last-passage Laplace Green’s function on a flat surface. In this paper, we obtained the Laplace Green’s function on a spherical surface and further developed the previous last-passage algorithm. We validate the last-passage algorithm on the spherical surface by obtaining charge density on a sphere held at unit potential. Then we demonstrate the last-passage algorithm by computing the mutual capacitance and charge distribution of two spheres and compare them with the previous research by J. Lekner. We find that the results are in good agreement.
It should be mentioned that the last-passage algorithm can be used for computing mutual capacitance in the semiconductor industry. As semiconducting devices have shrunk and the chip sizes have increased, in high-speed chip designs, computing interconnect capacitance has become very important. This electrical parasitic mutual capacitance is now a dominant factor in the high-speed chip design. To extract the capacitance of complex three-dimensional interconnect geometry in the ever-decreasing semiconductor device fabrication, the floating random-walk method [19] has been developed and commercialized. The newly developed last-passage algorithm can be used for mutual capacitance when the geometry of the mutual capacitance includes spheres. That means, when the geometry of the capacitors is a combination of flat and spherical surfaces, we can use the last-passage algorithms now.
The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.
References
Freidlin, M.: Functional Integration and Partial Differential Equations. Princeton University Press, Princeton, New Jersey (1985)
Chung, K.L., Zhao, Z.: From Brownian Motion to Schrödinger’s Equation. Springer-Verlag, Berlin (1995)
Douglas, J.F., Zhou, H.-X., Hubbard, J.B.: Intrinsic viscosity and the electrical polarizability of arbitrarily shaped objects. Phys. Rev. E 49, 5319 (1994)
Hwang, C.-O., Given, J.A., Kim, Y., Lee, S., Lee, S.: First- and Last-Passage Monte Carlo Algorithms for Charge Density on a Conducting Surface Proceedings of the 2016 International Conference on Innovative Material Science and Technology, pp. 139–147 (2016)
Given, J.A., Hwang, C.-O., Mascagni, M.: First- and last-passage Monte Carlo algorithms for the charge density distribution on a conducting surface. Phys. Rev. E 66, 056704 (2002)
Gould, H., Tobochnik, J.: An Introduction to Computer Simulation Methods. Addison-Wesley Publishing Company (1996)
Müller, M.E.: An introduction to computer simulation methods. Ann. Math. Stat. 27, 569 (1956)
Yan, C., Cai, W., Zeng, X.: A parallel method for solving laplace equations with Dirichlet data using local boundary integral equations and random walks. SIAM J. Sci. Comput. 35–4, B868 (2013)
Hwang, C.-O., Hong, S., Kim, J.: Off-centered “Walk-on-Spheres” (WOS) algorithm. J. Comp. Phys. 303, 331 (2015)
Given, J.A., Hubbard, J.B., Douglas, J.F.: Last-passage Monte Carlo algorithm for the mutual capacitance. J. Chem. Phys. 106, 3721 (1997)
Hwang, C.-O., Given, J.A., Mascagni, M.: The simulation-tabulation method for classical diffusion Monte Carlo. Phys. Fluids A 12, 1699 (2000)
Mansfield, M.L., Douglas, J.F., Garboczi, E.J.: Intrinsic viscosity and the electrical polarizability of arbitrarily shaped objects. Phys. Rev. E 64, 061401 (2001)
Hwang, C.-O., Mascagni, M.: Electrical capacitance of the unit cube. J. Appl. Phys. 95, 3798 (2004)
Hwang, C.-O., Kim, M.: Parallel plates algorithms. Math. Comput. Simul. (submitted)
Hwang, C.-O., Given, J.A., Kim, Y., Lee, S., Lee, S.: First- and last-passage Monte Carlo algorithms for charge density on a conducting surface. In: Proceedings of the 2016 International Conference on Innovative Material Science and Technology (IMST 2016), vol. 139, p. 139 (2016)
Chung, K.L.: Green, Brown, and Probability. World Scientific, Singapore (1995)
Hwang, C.-O., Given, J.A.: Last-passage Monte Carlo algorithm for the mutual capacitance. Phys. Rev. E 74, 027701 (2006)
Lekner, J.: Electrostatics of two charged conducting spheres. Proc. R. Soc. A 468, 2829 (2012)
Yu, W., Xu, Z., Li, B., Zhuo, C.: Floating random walk-based capacitance simulation considering general floating metals. IEEE Trans. CAD 37(8), 1711 (2018)
Acknowledgements
This work was supported by the National Research Foundation of Korea (NRF) grant (No. 2017R1E1A1A03070543) funded by the Korea government (Ministry of Science and Information & Communication Technology (ICT)). In addition, this work was supported by the GIST Research Institute (GRI) grant funded by GIST in 2021.
Author information
Authors and Affiliations
Corresponding author
Additional information
Publisher's Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
About this article
Cite this article
Yu, U., Lee, YM. & Hwang, CO. Last-passage Monte Carlo Algorithm for Charge Density on a Conducting Spherical Surface. J Sci Comput 88, 82 (2021). https://doi.org/10.1007/s10915-021-01594-w
Received:
Revised:
Accepted:
Published:
DOI: https://doi.org/10.1007/s10915-021-01594-w