A method to improve the accuracy of calculating the charge density on a plane boundary
By changing the first type of integral equation in the BIE-WOS method to the second type of integral equation, and combining large-scale parallel computing technology, the problem of low calculation accuracy in the existing technology is solved, and a more accurate calculation of the charge density of interconnected lines and dielectric surfaces in nano-integrated circuits is achieved.
Patent Information
- Application Number
- CN202011312648.4
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2020-11-20
- Publication Date
- 2025-04-29
- Estimated Expiration
- 2040-11-20
AI Technical Summary
In the prior art, the calculation accuracy of the planar surface charge density extraction method based on the boundary integral equation and the random mixing method is not high, and it is difficult to effectively solve the calculation problem of the full-chip-level interconnection line and dielectric surface charge density in nano-integrated circuits.
The first type of integral equation used in the original BIE-WOS method is changed to the second type of integral equation, combined with large-scale parallel computing technology, by constructing hemispheres and plane discs, it is converted into the second type of integral equation with closed surfaces, and numerical integration is performed to improve the calculation accuracy.
It effectively improves the accuracy of the algorithm and can more accurately calculate the surface charge density of interconnected lines and dielectric surfaces in nano-integrated circuits, solving the problem of full-chip-level computing.
Smart Images

Figure CN114519314B_ABST
Abstract
Description
Technical Field
[0001] The present invention belongs to the technical field of integrated circuit computer aided design / electronic design automation (CAD / EDA), and in particular relates to a method for improving the accuracy of calculating the charge density of a plane boundary surface by using a boundary integral equation and a random hybrid method (Boundary Integral Equation-Random Walk on Spheres, BIE-WOS) in the direction of integrated circuit parasitic parameter extraction. Background Art
[0002] The prior art discloses that the performance of interconnects is one of the key factors that determine circuit performance and reliability. As the feature size of integrated circuits continues to shrink, the density of layout patterns continues to increase, the spacing between interconnects continues to shrink, and the surface charge density of interconnects and dielectric surfaces continues to rise. In order to prevent dielectric breakdown, it is necessary to accurately calculate the surface charge density of interconnects and dielectrics. In addition, after obtaining the surface charge density of interconnects, it is also easy to extract its parasitic capacitance. Therefore, the calculation of the surface charge density of interconnects and dielectrics becomes an important goal and difficulty of parasitic parameter extraction tools. There are hundreds of millions of interconnects in the large-scale integrated circuit design layout. A large number of interconnects and dielectric materials constitute a complex electromagnetic environment. Accurately calculating the surface charge density on their surfaces is a considerable challenge for EDA software design.
[0003] At present, the methods for extracting capacitor parasitic parameters based on field solvers can be roughly divided into two categories, namely deterministic methods and stochastic methods. Deterministic methods include finite differential method (FDM), finite element method (FEM) and boundary element method (BEM). For example, Raphael based on finite differential method, FastCap based on indirect BEM[1, 2], QMM-BEM based on direct BEM[3], HiCAP and PhiCap based on hierarchical extraction[4, 5] and ParAFEMCap based on finite element method[6]. Deterministic methods have the advantages of high accuracy and wide adaptability, but they are global methods, that is, even if only the solution of the local field is needed, global discretization must be performed, a complete set of linear equations must be established and solved iteratively to obtain the solution on all boundaries in the region. For a full chip layout with hundreds of millions of geometric shapes, it is unrealistic to use a global deterministic method to solve it under current computer technology conditions.
[0004] The random method is a local solution method. Generally speaking, the random method is a specific application of the Feynman-Kac theorem, which expresses the potential (or surface charge density) as a weighted average of the potential values at the regional boundary [7]. The random method can effectively calculate the solution of the local region in the whole chip without obtaining the field distribution of the entire chip. The industry standard parasitic parameter extraction tools Rapid 3D and QuickCap use the random method. In addition, the random method is easy to implement large-scale parallel calculations. However, since the Dirichlet boundary is an absorbing boundary, the current commercial parasitic parameter extraction tools based on the random method cannot obtain the surface charge density of the interconnect or dielectric surface. A more special last-passage random method [7,8] can calculate the surface charge density of the conductor surface, but it is only applicable to the conductor plane where the Dirichlet boundary condition is a constant, and cannot handle the surface charge density extraction problem of the dielectric plane.
[0005] Based on the current status of the prior art, the inventors of this application intend to address the problem of low calculation accuracy in a method for extracting plane surface charge density that combines boundary integral equations and random mixing methods, and provide a method that changes the integral equation of the first kind (IEFK) used in the original BIE-WOS method to the integral equation of the second kind (IESK), which can effectively improve the accuracy of the algorithm. The method proposed in this invention, after applying large-scale parallel computing technology, will help solve the problem of calculating the surface charge density of interconnects and dielectrics at the full chip level in nano-integrated circuits.
[0006] References:
[0007] [1] K.Nabors and J.White, "FastCap: A multipole accelerated 3-Dcapacitance extraction program," IEEE Trans.Computer-AidedDesign, vol.10, pp.1447–1459, Nov.1991.
[0008] [2]J.Phillips and J.White,“A Precorrected-FFT Method For CapacitanceExtraction of Complicated 3-D Structures,”in Proc.IEEE / ACM Int.Conf.Computer-Aided Design,Nov. 1994,pp.268–271.
[0009] [3]Wenjian Yu,and Zeyi Wang,“Enhanced QMM-BEM Solver for Three-Dimensional Multiple-Dielectric Capacitance Extraction Within the FiniteDomain.”IEEE TRANSACTIONS ON MICROWAVE THEORY AND TECHNIQUES,VOL.52,NO.2,FEBRUARY2004pp.560-566
[0010] [4]W.Shi,J.Liu,N.Kakani,and T.Yu,“A fast hierarchical algorithm forthree-dimensional capacitance extraction,”IEEE Trans.Comput.-AidedDes.Inter.Circuits Syst.,vol.21,no.3, pp.330–336,Mar.2002.
[0011] [5]S.Yan,V.Sarim,andW.Shi,“Sparsetransformation andpreconditionersfor3-D capacitance extraction,”IEEE Trans.Comput.-Aided Des.Integr.CircuitsSyst.,vol.24,no.9,pp. 1420–1426,Sep.2005.
[0012] [6]Chen,G.and Zhu,H.and Cui,T.and Chen,Z.and Zeng,X.and Cai,W.“ParAFEMCap:A Parallel Adaptive Finite-ElementMethodfor 3-D VLSIInterconnectCapacitance Extraction.” IEEE Trans.onMicrowaveTheoryandTechniques,no.9,pp.218–231,2012.
[0013] [7]Chi-Ok Hwang and James A.Given,“Last-passage Monte Carlo algorithmfor mutual capacitance,”PHYSICALREVIEWE 74,0277012006.
[0014] [8]James A.Given,Chi-Ok Hwang and Michael Mascagni,“First-and Last-passage Monte Carlo algorithms for the charge density distribution on aconducting surface.”PHYSICAL REVIEWE 66,056704~2002
[0015] [9]Changhao Yan,Wei Cai and Xuan Zeng,“A Parallel Method For SolvingLaplace Equations With Dirichlet Data Using Local Boundary Integral EquationsAnd Random Walks.”SIAM Journal on Scientific Computing,2013(accepted)
[0016]
[10] M. Guiggiani, G. Krishnasamy, T. J. Rudolphi, and F. J. Rizzo, A General Algorithm for the Numerical Solution of Hypersingular Boundary Integral Equations, ASME. J. Appl. Mech., 59 (1992), 604-614. Summary of the invention
[0017] The purpose of the present invention is to provide a method for improving the accuracy of calculating the charge density of the plane boundary surface by using the boundary integral equation and the random mixing method, aiming at the core problems of extracting the parasitic parameters of the interconnection lines of nano-integrated circuits and extracting the charge density of the boundary surface between the conductor and the dielectric.
[0018] Specifically, the present application aims to solve the problem of low calculation accuracy in a method for extracting plane surface charge density by combining the boundary integral equation and the random mixing method, and proposes to change the integral equation of the first kind (IEFK) used in the original BIE-WOS method to the integral equation of the second kind (IESK), which can effectively improve the accuracy of the algorithm.
[0019] The method proposed by the present invention improves the accuracy of calculating the surface charge density of a plane boundary by using a boundary integral equation and a random mixing method. First, a hemisphere is constructed with the point to be calculated on the boundary as the center of the sphere, and the intersection with the boundary of the region is a plane disk. The surface charge density at the point to be calculated is converted into a second-class integral equation on a closed surface composed of a hemisphere and a plane disk. Then, numerical integration is performed on the hemisphere and the plane disk respectively, wherein the unknown potential at the integral point included in the hemisphere surface integral can be obtained by converting the hemisphere surface at The surface is discretized into a uniform grid, and the potentials of these uniform grid points are obtained using the WOS method and bilinear interpolation. The integral on the plane disk contains supersingular integrals, which can be regularized by numerical techniques to reduce the singularity. Finally, the integral results of the hemispherical surface and the plane disk are summed to obtain the surface charge density at the point to be determined on the final boundary.
[0020] The method proposed in the present invention can be applied to both conductor and dielectric plane boundaries. In the calculation of conductor plane boundaries, since the surface potential of the conductor is a constant, the integral on the plane disk is always 0. Therefore, the method for calculating the charge density on the conductor plane boundary is only a special case of the method for calculating the charge density on the dielectric plane boundary.
[0021] The problems to be solved by the method of the present invention are as follows Figure 1 As shown, the potential u(x) in the region Ω satisfies the Laplace equation, is the region boundary, the boundary material is a conductor, a dielectric, or a combination thereof, and the boundary The potential distribution function u0(x) on is known, that is, it satisfies:
[0022]
[0023] The method proposed by the present invention to improve the accuracy of calculating the charge density of a plane boundary surface by using the boundary integral equation and the random mixing method is shown in the flowchart as follows: Figure 2 As shown;
[0024] Input parameters: Region boundaries Geometric information, region boundaries The potential distribution function u0(x) on the surface, the position coordinates of the point x to be determined, the radius a of the hemisphere, and the number of Gaussian integration points N on the hemisphere g1 , the number of Gaussian integration points on the disk plane N g2 , the number of sampling paths N max , absorbing boundary thickness ε, maximum sphere radius R max , the hemisphere is Surface (θ and The number of discrete points N on the spherical coordinate system is the elevation angle and the circumference angle respectively. θ ,
[0025] Output result: surface charge density at point x.
[0026] Step 1: Construct a hemisphere with radius a at point x;
[0027] Step 2: Place the hemisphere on The surface is discretized into uniform grids, and the potentials of these uniform grid points are obtained by using the WOS method, and the potentials of all integration points on the hemispherical surface are obtained by bilinear interpolation;
[0028] Step 3: Use the second kind integral equation method to calculate the surface charge density σ at point x s (x).
[0029] More specifically, in step 1 of the present invention, the selection of radius a should satisfy:
[0030] (1) The intersection of the hemisphere and the region boundary should be a disk plane; it should be noted that this method does not require the entire region boundary is a plane, and only the local area where the hemisphere intersects the region boundary needs to be a plane;
[0031] (2) Except for the disk plane that intersects the region boundary, the sphere does not intersect any other boundary (interconnection line boundary or medium boundary);
[0032] (3) If the radius a is too small, the accuracy of the final result will be affected.
[0033] In step 2 of the method of the present invention, the potential of the integral point on the hemispherical surface is obtained by using the following sub-steps:
[0034] Step 2.1: On the plane, the hemisphere is discretized into Axis discrete points N θ , Discretize into a uniform grid, the coordinates of each grid point are:
[0035]
[0036] Step 2.2: Use the random method to calculate the potential u (Grid ii,jj ) and save the grid point potential; the following uses the commonly used WOS method as an example to illustrate how to calculate any grid point y0 = Grid ii,jj ,0≤ii≤N θ , The potential u(y0) at the position is shown in the flow chart. Figure 3 As shown, the specific steps are as follows:
[0037] Input parameters: Region boundaries Geometric information, region boundaries Potential distribution on the surface, coordinate value of y0, total number of paths N max , absorbing boundary thickness ε, maximum sphere radius R max ;
[0038] Output result: the potential u(y0) of y0.
[0039] Step 2.2.1: Initialization: current path number n = 0, u(y0) = 0;
[0040] Step 2.2.2: Determine the current number of paths n and N max If n≥N max , then the potential at the output y0 is u(y0) / n; otherwise, go to step 2.2.3;
[0041] Step 2.2.3: Initialize the current point y=y0;
[0042] Step 2.2.4: Calculate the current point y and the region boundary The shortest distance d;
[0043] Step 2.2.5: If distance d>R max , it is considered to have reached infinity, n = n + 1, and jump to step 2.2.2; otherwise, jump to step 2.2.6;
[0044] Step 2.2.6: If the distance d < ε, it is absorbed by the boundary, and the potential u0(x) at the point x closest to the boundary is obtained, and it is accumulated into u(y0), n = n + 1, and jump to step 2.2.2; otherwise, jump to step 2.2.7;
[0045] Step 2.2.7: Randomly select the next point y with uniform distribution on the sphere with the current point y as the center and d as the radius, and use it as the current point y, and jump to step 2.2.4.
[0046] It should be noted that the hemispherical grid point potentials calculated in the above steps 2.1 and 2.2 are only calculated once and saved for subsequent calculations;
[0047] Step 2.3: Determine the coordinates of the numerical integration points on the hemisphere.
[0048] Usually, numerical integration method is used to obtain the integral on the hemispherical surface, so the coordinate value of the integration point needs to be obtained first; taking the commonly used Gaussian integral as an example, using the two-dimensional N g1 ×N g1 Point Gaussian integral, Gaussian integral point y in spherical coordinate system i,j The coordinate values are:
[0049]
[0050] Where a is the radius of the hemisphere, ξ i and j All are one-dimensional in the [-1,1] region N g1 The location of the Gaussian point in the point-Gaussian integral;
[0051] Step 2.4: Use the potential values at the hemispherical grid points to calculate the Gaussian integral point y using bilinear interpolation i,j The potential value of
[0052] The method of the present invention adopts the second kind of integral equation method in step 3 to calculate the surface charge density σ at point x s (x) is calculated using the following method:
[0053] Using the spherical Green's function, the second kind of boundary integral equation on the boundary of the closed region formed by the hemispherical surface Γ and the plane S is established;
[0054] By area Ω S The second boundary integral equation on the boundary can be obtained:
[0055]
[0056] Where u(x) and u(y) are the potentials at point x and point y respectively, n y is the unit external normal vector of point y, G1(x,y) is the spherical Green's function, that is, the Green's function G1(x,y) is the spherical area S with a radius a a The Green's function under the homogeneous boundary satisfies:
[0057]
[0058] The mirror image method can be used to solve:
[0059]
[0060] Where a is the radius of the hemisphere, and the subscript 1 indicates that an image charge is introduced into the Green function, i.e., the first term is the original charge point inside the sphere, and the second term is the mirror point outside the sphere.
[0061] Using the spherical Green's function G1(x,y), we can make Therefore, there is no need to discretize the sphere.
[0062] On the boundary of the closed area enclosed by the plane S and the upper hemisphere Γ, only the patch surface S needs to be discretized. If the boundary integral equation (4) is discretized using the boundary element method, it can be transformed into the following linear equation system:
[0063]
[0064] in,
[0065]
[0066]
[0067] where pf is the Hadamard finite partial product of the hypersingular integral, which can be calculated using regularization techniques [9,11].
[0068] Considering that the surface S is a plane, when points x and y are both on the plane S, the second-type integral equation (7) can be further simplified. The first-order partial derivative of the spherical Green's function is:
[0069]
[0070]
[0071] Where x={x1,x2,x3},y={y1,y2,y3},x·y=x1y1+x2y2+x3y3, r=|r|, r=xy.
[0072] So we can get:
[0073]
[0074] The first term of the Cauchy principal integral is zero, because the jump term 1 / 2 on the smooth boundary has been considered in formula (7), and the second term is zero because it has no singularities.
[0075] Substituting formula (10) into formula (7), we can obtain:
[0076]
[0077] in,
[0078]
[0079] Where x = {x1, x2, x3}, y = {y1, y2, y3}, y I =ky, r=|yx|,r I =|y I -x|,
[0080] If y∈Γ, Equation (12) has no singularity, because x lies on the plane S, so r is always greater than zero. However, when x lies near the edge Γ∩S, it becomes nearly singular.
[0081] When y∈S, formula (12) includes and Part I It has supersingularity when x and y coincide, but it can be regularized using the method in
[10] . Part II There is no singularity, but when x is near the edge Γ∩S, it becomes a nearly supersingular integral.
[0082] Using the second kind integral equation method, calculate the surface charge density σ at point x s(x) can be directly obtained by numerical integration of formula (11) without surface discretization and solving linear equations, which greatly simplifies the calculation process. Compared with the method using the first kind of integral equation, although the integral kernel singularity of the second kind of integral equation is higher and a more complex numerical calculation method is required to correctly handle the singular integral, the second kind of integral equation method proposed in the present invention can effectively improve the accuracy of the algorithm.
[0083] The advantages of the method of the present invention are:
[0084] (1) Compared with the rough surface charge density extraction method using the second kind boundary integral equation, combining the boundary integral equation and the random mixing method, the present invention improves the particularity of the plane problem. The proposed method does not need to discretize the interconnection line and the dielectric boundary surface, does not need to establish a complex linear equation system, and does not need to adopt the boundary integral equation algorithm framework. It only needs to implement the random method and numerical integration. The programming is simple and easy to implement.
[0085] (2) Compared with the method using the first kind of integral equation, although the integral kernel singularity of the second kind of integral equation is higher and a more complex numerical calculation method is required to correctly handle the singular integral, the method using the second kind of integral equation proposed in the present invention can effectively improve the accuracy of the algorithm;
[0086] (3) The method of the present invention still has the advantage of the natural parallelism of the random method and is easy to implement large-scale parallel computing.
[0087] The method proposed in the present invention aims at the core problems of extracting parasitic parameters of interconnects of nano-integrated circuits and extracting charge density of conductor and dielectric boundary surfaces. After applying large-scale parallel computing technology, it will help solve the problem of calculating charge density of interconnects and dielectric surfaces at the full chip level in nano-integrated circuits. BRIEF DESCRIPTION OF THE DRAWINGS
[0088] Figure 1 , schematic diagram of the method used to calculate the surface charge density at a point x on a plane boundary.
[0089] Figure 2 , which is a flow chart of the method for calculating the charge density of a plane boundary surface in the present invention.
[0090] Figure 3 , which is a flow chart for calculating the integration point potential using the WOS method.
[0091] Figure 4 a, Schematic diagram of a half-plane example separated by an infinite plane.
[0092] Figure 4 b, Comparison of the accuracy of the first type (IEFK) and the second type (IESK) integral equations,
[0093] Figure 4 c, right-hand side of integral equations of the first and second kind b S , b Γ Compare with b. DETAILED DESCRIPTION
[0094] The method of the present invention is described below through the implementation process of a specific example.
[0095] Implementation of the case study
[0096] Figure 4 A simple half-plane test case is given in Figure a, where the plane z = 0 divides free space into two regions with dielectric constants ò0 and ò1, and a unit point charge q is located at r s =(0,0,-h). The Neumann value on the plane z=0 needs to be calculated.
[0097] The potential u(r) in the region z>0 satisfies the Laplace equation ▽ 2 u(r) = 0, the potential u(r) on the boundary z = 0 is known. The potential u(r) and charge density on the boundary z = 0 The analytical solution is:
[0098]
[0099] According to the BIE-WOS method, a hemisphere with radius a=1 is constructed with o=(0,0,0) as the center, and the first and second boundary integral equations are used to solve them. If the distance h of the charge q is changed, it is h=1, 10, 20, 50 and 100 respectively. Figure 4 b gives the relative error of all discrete surface elements with a radius of 0.7a. At the same time, Table 1 lists the condition number and accuracy comparison results of the first and second kind of integral equations.
[0100] It can be seen that as the charge distance h increases from 1 to 100, the maximum relative error of IEFK increases from 0.57% to 13.82%, while the error of IESK only increases from 0.47% to 1.41%. cond(A) in Table 1 verifies that the condition number of IEFK (358.5) is worse than the condition number of IESK (1), but this is not the main reason for the decrease in accuracy of IEFK. Because the condition numbers of the matrix A of IEFK remain unchanged when the charge distance h increases from 1 to 100, they cannot explain the large degradation of the accuracy of the final result of the first kind of integral from 0.57% to 13.82%. However, for different h, the right side of the matrix system will change, and the results show that this is the real reason for the different accuracy of the BIE-WOS method of IEFK and IESK.
[0101] Table 1: Comparison of condition numbers and accuracy for integral equations of the first kind (IEFK) and the second kind (IESK)
[0102]
[0103]
[0104] For different h, Figure 4 c gives the results of the components of the right-hand side term b of the first and second integral equations, where the x-axis is the index number of the discrete surface element and the y-axis is b S (x), and the right side b(x). Figure 4 c It can be seen that when h changes from 1 to 50, for IEFK, b S (x)>0,b Γ (x)<0, when h=50, b(x) is close to zero, and due to numerical cancellation, b(x) has a huge error. However, for IESK, b S (x)>0, b Γ (x)>0 and b(x)>0, no numerical cancellation occurs on b(x), so the final high calculation accuracy can be guaranteed.
Claims
1. A method for improving the accuracy of calculating the charge density of a plane boundary surface by using a boundary integral equation and a random mixing method, characterized in that: The method comprises the following steps: constructing a hemisphere with the point to be determined on the boundary as the center, and forming a plane disk intersecting with the boundary of the region; converting the surface charge density at the point to be determined into a second-order integral equation on a closed surface composed of the hemisphere and the plane disk; then numerically integrating the hemisphere and the plane disk respectively; and finally summing the integral results of the hemisphere and the plane disk to obtain the surface charge density at the point to be determined on the final boundary; Includes steps: Input parameters: Region boundaries Geometric information, region boundaries The potential distribution function u0(x) on the surface, the position coordinates of the point x to be determined, the radius a of the hemisphere, and the number of Gaussian integration points N on the hemisphere g1 , the number of Gaussian integration points on the disk plane N g2 , the number of sampling paths N max , absorbing boundary thickness ε, maximum sphere radius R max , the hemisphere is at θ- The number of discrete points N on the surface θ , Among them, θ and are the elevation angle and circular angle in the spherical coordinate system respectively; Output result: surface charge density at point x; Step 1: Construct a hemisphere with radius a at point x; Step 2: Place the hemisphere at θ- The surface is discretized into uniform grids, and the potentials of these uniform grid points are obtained by using the WOS method, and the potentials of all integration points on the hemispherical surface are obtained by bilinear interpolation; Step 3: Use the second kind integral equation method to calculate the surface charge density σ at point x s (x); In step 3, the second kind of integral equation method is used to calculate the surface charge density σ at point x. s (x); Using the spherical Green's function, the second kind of boundary integral equation on the boundary of the closed area formed by the hemispherical surface I' and the plane S is established; By area Ω S The second boundary integral equation on the boundary is obtained: Where u(x) and u(y) are the potentials at point x and point y respectively, n y is the unit external normal vector of point y, G1(x, y) is the spherical Green's function, that is, the Green's function G1(x, y) is the spherical area S with a radius a a The Green's function under the homogeneous boundary satisfies: Using the mirror method, we can get: Where a is the radius of the hemisphere, and the subscript 1 indicates that an image charge is introduced into the Green function, i.e., the first term is the original charge point inside the sphere, and the second term is the mirror point outside the sphere; Using the spherical Green's function G1(x, y), There is no need to discretize the sphere; On the boundary of the closed area enclosed by the plane S and the upper hemisphere Γ, only the patch surface S needs to be discretized; the boundary integral equation (4) is discretized using the boundary element method and transformed into the following linear equation system: in, where pf is the Hadamard finite partial product of the supersingular integral, calculated using the regularization technique; Considering that the surface S is a plane, when points x and y are both on the plane S, the second-type integral equation (7) is simplified; the first-order partial derivative of the spherical Green's function is: where \(x = \{x_1, x_2, x_3\}\), \(y = \{y_1, y_2, y_3\}\), \(x\cdot y=x_1y_1 + x_2y_2+x_3y_3\), \(r = |r|\), \(r=x - y\); So we can get: The first term of the Cauchy principal integral is zero. The jump term 1 / 2 on the smooth boundary has been considered in formula (7). The second term is zero, because it has no singularity; Substituting formula (10) into formula (7), we obtain: in, For example,x={x1,x2,x3},y={y1,y2,y3},y I =ky, r=|yx|,r I |y I -x|, If y∈Γ, formula (12) has no singularity, because x lies on the plane S, so r is always greater than zero; when x lies near the edge Γ∩S, it will be a nearly singular integral; When y∈S, formula (12) includes and Part I It has supersingularity when x and y coincide, and can be regularized by the method in (10); Part II There is no singular point. When x is near the edge Γ∩S, it becomes a nearly supersingular integral; Using the second kind integral equation method, calculate the surface charge density σ at point x s (x) can be directly obtained by numerical integration of formula (10), without the need for surface discretization and solving linear equations, thus simplifying the calculation process and improving the accuracy of the algorithm.
2. The method according to claim 1, characterized in that: In step 1, the radius a should be selected to satisfy: (1) The part where the hemisphere intersects with the region boundary should be a disk plane; the method described does not require the entire region boundary to be a plane, but only requires the local area where the hemisphere intersects with the region boundary to be a plane; (2) The sphere does not intersect any other boundary except the disk plane that intersects the region boundary. Other boundaries are interconnection line boundaries or dielectric boundaries. (3) If the radius a is too small, the accuracy of the final result will be affected.
3. The method according to claim 1, characterized in that: In step 2, the potential of the integral point on the hemispherical surface is obtained by using the W0S method and bilinear interpolation, and the potential of the integral point is obtained by the following sub-steps: Step 2.1: In θ- On the plane, the hemisphere is discretized into Axis discrete points N θ , Discretize into a uniform grid, the coordinates of each grid point are: Step 2.2: Use the random method to calculate the potential u (Grid ii,jj ) and save the grid point potential; Use the WOS method to calculate any grid point y0 = Grid ii,jj ,0≤ii≤N θ , The specific steps are as follows: Input parameters: Region boundaries Geometric information, region boundaries Potential distribution on the surface, coordinate value of y0, total number of paths N max , absorbing boundary thickness ε, maximum sphere radius R max ; Output result: the potential u(y0) of y0; Step 2.2.1: Initialization: current path number n = 0, u(y0) = 0; Step 2.2.2: Determine the current number of paths n and N max If n≥N max , then the potential at the output y0 is u(y0) / n; Otherwise go to step 2.2.3; Step 2.2.3: Initialize the current point y=y0; Step 2.2.4: Calculate the current point y and the region boundary The shortest distance d; Step 2.2.5: If distance d>R max , it is considered to have reached infinity, n = n + 1, and jump to step 2.2.2; otherwise, jump to step 2.2.6; Step 2.2.6: If the distance d < ε, it is absorbed by the boundary, and the potential u0(x) at the point x closest to the boundary is obtained, and it is accumulated into u(y0), n = n + 1, and jump to step 2.2.2; otherwise, jump to step 2.2.7; Step 2.2.7: On the sphere with the current point y as the center and d as the radius, randomly select the next point y with uniform distribution, and use it as the current point y, and jump to step 2.2.4; The hemispherical grid point potentials calculated in steps 2.1 and 2.2 above are calculated only once and saved for subsequent calculations. Step 2.3: Determine the coordinates of the numerical integration points on the hemisphere. Using numerical integration method to obtain the integral on the hemispherical surface, we first need to get the coordinates of the integration point; using Gaussian integration, using two-dimensional N g1 ×N g1 Point Gaussian integral, Gaussian integral point y in spherical coordinate system i,j The coordinate values are: Where a is the radius of the hemisphere, ξ i and j All are one-dimensional in the [-1, 1] region N g1 The location of the Gaussian point in the point-Gaussian integral; Step 2.4: Use the potential values at the hemispherical grid points to calculate the Gaussian integral point y using bilinear interpolation i,j The potential value.
4. The method according to claim 1, characterized in that: The method is applicable to the plane boundaries of conductors, dielectrics and mixed conductor-dielectric planes.
5. The method according to claim 1, characterized in that: In the calculation of the plane boundary of the conductor, since the surface potential of the conductor is a constant, the integral on the plane disk is always 0; therefore, the method for calculating the charge density on the plane boundary of the conductor is a method for calculating the charge density on the plane boundary of the dielectric.