# Improving the Accuracy of the Charge Simulation Method for Numerical Conformal Mapping.

1. Introduction

Conformal mapping appears in complex function theory, which plays important roles for applications in fluid mechanics, image processing, plane elasticity theory, and so on [1, 2]. As is well known, solutions to conformal mapping can be divided into analytical and numerical calculation method. The theoretical basis of analytical method is the famous Riemann Theorem, which proves the existence of conformal mapping function but cannot give the expression for a specific function. This implies that we must use numerical calculation method to solve some related problems when the conformal mapping functions are complex and diverse. Symm [3, 4] proposed the integral equation method for the conformal mapping of a simply connected domain bounded by a closed Jordan curve onto the unit disk, of its exterior onto the exterior of the unit disk and of a doubly connected domain onto a concentric circular annulus. In [5], a numerical integral method is proposed to find the specific length-to-width ratio, which is able to improve the accuracy of the numerical conformal mapping on the boundaries and makes the mapping method more reliable when solving problems that are sensitive to accuracy. Nasser [6-9] proposed a method which can be used to compute the mapping function onto the five canonical regions in a unified way. This method is based on a boundary integral equation with the generalized Neumann kernel. In the year 1969, Steinberger firstly proposed the charge simulation method to study the electric field calculation [10]. Following the work of Steinberger, Amano [11-13] developed the charge simulation methods and obtained many nice and important results on the charge simulation method and numerical conformal mapping. In fact, compared to traditional methods, the charge simulation method has a higher accuracy evaluation.

The matrix of Nasser method is well-conditioned, but the coefficient matrix of the constraint equation of Amano method is ill-conditioned [8]. It is very important to solve the constraint equation in order to obtain the high accuracy charges and the conformal mapping radius. Wilkinson's iteration refinement for linear systems [14] is well known for solution in which Cholesky decomposition is repeatedly used. It has been shown that if the linear equations are not too ill-conditioned, the iterative sequence of the approximate solution will be convergent to the solution. But this method may not be useful for the too ill-conditioned system of linear equations. Reference [15] established a new iteration improvement of the solution based on the dynamic system of ordinary differential equations (ODEs), and the new method is more effective than Wilkinson's method. However, the number of iterations required for this method is often more. This paper applies the Runge-Kutta method for the constraint equations of doubly connected domains; the new approach reduced the number of iterations and improved the accuracy.

This paper is organized as follows. In Section 2, this paper describes the numerical conformal mapping of doubly connected domains based on the charge simulation method. In Section 3, a novel algorithm based on Runge-Kutta method for numerical conformal mapping of doubly connected domains is proposed. In Section 4, it gives some numerical examples and shows the superiority of developed method stated. Section 5 is devoted to the statement of conclusions.

2. Numerical Conformal Mapping of Doubly Connected Domains Based on the Charge Simulation Method

This section will describe the conformal mapping of a finite doubly connected domain D bounded by two closed Jordan curves [C.sub.1] and [C.sub.2] given in the z-plane onto a circular annulus [mu] < [absolute value of [omega]] < 1 in the w-plane (see Figure 1) [13, 16-19]. The origin z = 0 lies inside the domain bounded by [C.sub.2]. Then, the mapping function f(z) can be expressed as

f (z) = [ze.sup.g(z)+ih(z)], z [member of] [bar.D] = D + [C.sub.1] +[C.sub.2]. (1)

Here g(z) and h(z) are conjugate harmonic functions in D and g(z) is the solution of the Dirichlet problem

[mathematical expression not reproducible], (2)

Based on the charge simulation method, g(z) and h(z) can be approximated by

[mathematical expression not reproducible], (3)

[mathematical expression not reproducible], (4)

respectively. The charge points [zeta] = 1, 2, ..., N) are arranged outside the given domain D. More precisely, N/2 charge points are arranged outside the domain bounded by [C.sub.1] and N/2 inside the domain bounded by [C.sub.2]. The charges [Q.sub.i] (i = 1, 2, ..., N) can be determined to satisfy the Dirichlet boundary condition of collection points [z.sub.j] (j = 1, 2, ..., N) arranged on the boundary [C.sub.1] and [C.sub.2]; that is,

[mathematical expression not reproducible]. (5)

According to the condition of existence of conjugate harmonic function, we have

[N.summation over (i=N/2+1)] [Q.sub.i] = 0. (6)

From (5) and (6), we obtain the constraint equations of doubly connected domains based on the charge simulation method. Here [a.sub.ij] = log [absolute value of [z.sub.j] - [[zeta].sub.i]].

3. Runge-Kutta Method for Numerical Conformal Mapping of Doubly Connected Domain

It is important to find an effective way to solve the constraint equation (7). We may consider linear equations of the form

Ax = b, (8)

where A [member of] [R.sub.(N+1)x(N+1)] and b [member of] [R.sup.N+1] and N represents the number of simulation charge points.

Normally, the higher precision of conformal mapping, the more points of simulation charge points. In this case, the condition number cond(A) = [parallel]A[parallel] [parallel][A.sup.-1][parallel] of coefficient matrix A is large. Hence (8) is ill-conditioned linear systems.

Consider the following preconditioned linear equations:

[A.sup.T]Ax = [A.sup.T]b. (9)

Eq. (9) is on both sides plus [micro]x, where [mu] [greater than or equal to] 0. Thus, we can construct iterative formula

([micro]I + [A.sup.T]A)[x.sub.n+1] = [mu][x.sub.n] + [A.sup.T]b, (n = 0, 1, ...). (10)

For the simplicity of expression, we define [mathematical expression not reproducible], and [??] = [micro]I + [A.sup.T]A. Assume that [x.sub.n+1] = [x.sub.n] + [y.sub.n], then

[mathematical expression not reproducible], (11)

This can be viewed as explicit Euler method [15, 20] with step h = 1 for solving the system of ordinary differential equations

[mathematical expression not reproducible], (12)

x(0) = [x.sub.0]. (13)

Runge-Kutta method is a high-precision single-step algorithm for numerical solutions of ordinary differential equations [21-24]. For integral differential equation (12), from [t.sub.n] to [t.sub.n+1], we obtain

[mathematical expression not reproducible]. (14)

To improve the formula order, it is necessary to increase the accuracy of the quadrature formula of the right-hand. Obviously, it is necessary to increase the quadrature node, for which the right-hand quadrature formula of (14) can be expressed as

[mathematical expression not reproducible]. (15)

In general, the more the nodes r, the higher the precision. In order to obtain the explicit method that facilitates calculating, formula (15) is expressed as

[x.sub.k+1] = [x.sub.k] + h[k.summation over (i=1)][c.sub.i][K.sub.i]. (16)

[K.sub.i] in (16) were given by

[mathematical expression not reproducible], (17)

where [mathematical expression not reproducible], and [[alpha].sub.ij] are constants given by [22].

Here we prove the convergence of the algorithm described above, taking fourth-order Runge-Kutta method as an example, and it is easy to prove other Runge-Kutta methods converge.

Theorem 1. For the initial value problems for systems of ordinary differential equations

[mathematical expression not reproducible]. (18)

When the function g(t, x) on its domain about x satisfies Lipschitz condition, the increment function [psi](t, x, h) = 1/6([K.sub.1] + 2[K.sub.2] + 2[K.sub.3] + [K.sub.4]) of the fourth-order Runge-Kutta method satisfies Lipschitz condition on x, where [t.sub.0] [less than or equal to] t [less than or equal to] [??], 0 [less than or equal to] h [less than or equal to] [??], and [absolute value of x] < [infinity].

Proof. For any x and x we obtain

[mathematical expression not reproducible]. (19)

Since g(t, x) on its domain about x satisfies Lipschitz condition, we have

[mathematical expression not reproducible]. (20)

From this we obtain

[mathematical expression not reproducible], (21)

where [mathematical expression not reproducible]. Since x and [??] are arbitrary, therefore [psi](t, x, h) satisfies Lipschitz condition on x.

Theorem 2. For the initial value problems for systems of ordinary differential equations

[dx/dt] = g(t, x), [t.sub.0] [less than or equal to] t [less than or equal to] [??], x(0) = [x.sub.0]. (22)

When the function g(t, x) on its domain about x satisfies Lipschitz condition, the fourth-order Runge-Kutta method convergence.

Proof. Note that [psi](t, x, 0) = g(t, x); then Theorem 1 implies fourth-order Runge-Kutta method convergence.

At this point, we can obtain the Runge-Kutta method for solving the constraint equations. Specific steps are given in Algorithm 1, where SolRD and ItMax represent the relative error and maximum number of iterations, respectively. Letting [x.sub.0] = 0, [mu], SolRD, and ItMax are determined by the constraint equations given in Section 4.

According to the previous analysis, the Runge-Kutta method for numerical conformal mapping based on the charge simulation method is summarized as follows.

(1) Given the number of charge points N, charge points [[zeta].sub.i] (i = 1, 2, ..., N) and constraint points [z.sub.j] (j = 1, 2, ..., N).

(2) By Algorithm 1, we calculate charges [Q.sub.i] (i = 1, 2, ..., N) and the approximation M of transformation radius.

(3) After calculating G(z) and H(z) by (3) and (4) for each point in the closed region [bar.D](D [union] [C.sub.1] [union] [C.sub.2]), including the boundary [C.sub.1] [union] [C.sub.2], substitute into the approximate conformal mapping function F(z) = [ze.sup.G(z)+iH(z)] to calculate the corresponding points.
```ALGORITHM 1: Runge-Kutta method for solving the constraint
equations.

Input: A, b, [x.sub.0], g, h, SolRD, ItMax.
Define g(x) = [([micro]I + [A.sup.T]A).sup.-1]([A.sup.T]b -
[A.sup.T]Ax).
for j = 1: ItMax
[K.sub.1] = g([x.sub.k]);
[K.sub.i] = g([x.sub.k] + h[[summa].sup.i-1.sub.j=1]
[[alpha].sub.ij][K.sub.j]);
[x.sub.k+1] = [x.sub.k] + h[[summa].sup.r.sub.i=1][c.sub.i]
[K.sub.i];
RD = norm([x.sub.k+1] - [x.sub.k])/norm([x.sub.k-1]);
if RD < SolRD, break; end if
end for
return [x.sub.k+1] and k.
```

4. Numerical Examples

In this section, two numerical experiments are made to show the effectiveness of Runge-Kutta method in comparison with Amano's method [13]. In the first example, we consider a doubly connected domain with exact solution given by [25]. In the second example, we consider the conformal mapping of a finite concentration ellipse onto a circular annulus. In the third example, we consider the conformal mapping of a finite doubly connected domain bounded by two squares onto a circular annulus. According to [13], the error formula of numerical conformal mapping of doubly connected domain is

[mathematical expression not reproducible]. (23)

The Runge-Kutta method for numerical conformal mapping of doubly connected domain has been implemented in MATLAB.

Example 1. Region D is the doubly connected domain with the boundaries [C.sub.1] : 10[e.sup.it] + 3[e.sup.i2t] and [C.sub.2] : 5[e.sup.-it] + 0.75[e.sup.-i2t], where 0 [less than or equal to] t [less than or equal to] 2[pi] (see Figure 2).

In this example, the analytic solution for domain D is known, which is given by

f(z) = [-10 + [square root of [10.sup.2] + 4 x 3z]/2x3] (24)

The values of parameters [mu], SolRD, and ItMax in this example are taken as 0.00002, [10.sup.-5], and 5000, respectively. Figure 2 shows the distribution of charge points for N = 240, half arranged outside the domain bounded by [C.sub.1] and half inside the domain bounded by [C.sub.2]. A comparison between the exact values of the conformal mapping f(z) and the approximate values, obtained by Amano's method and by our method with the same number of collection points, at various points on the boundaries [C.sub.1] and [C.sub.2], is given in Table 1. From Figure 3 and Table 1, we can see that our method can achieve higher numerical conformal mapping accuracy. The error of numerical conformal mapping is computed by (23). The approximate conformal mapping function F(z) maps Figures 4 and 5, from which we can see the effectiveness of our method. From Table 2, we can observe that the CPU time of our method is shorter than Amano's method.

Example 2. Region D is the doubly connected domain with the boundaries C: : [x.sup.2]/[8.sup.2] + [y.sup.2]/[4.sup.2] = land [C.sub.2]: [x.sup.2]/[4.sup.2] + [y.sup.2]/[2.sup.2] = 1 (see Figure 6).

For this example, let [mu] = 0.000002; the relative error SolRD and maximum number of iterations ItMax are [10.sup.-8] and 5000. Figure 6 shows the distribution of charge points for N = 200, half arranged outside the domain bounded by [C.sub.1] and half inside the domain bounded by [C.sub.2]. In Figure 3, we consider the numerical conformal mapping error behavior based on the Runge-Kutta methods when changing the number of charge points. Looking at Figure 7, we can observe that numerical conformal mapping error decreases by changing the number of charge points. The accuracy of numerical conformal mapping based on the Runge-Kutta is higher than Amano's method when the number of charge points is more than 110. Some numerical results of CPU time are listed in Table 3; we can see that the CPU time of Amanos method is longer than our method. When the number of charge points N = 170, the accuracy of numerical conformal mapping based on the Runge-Kutta method is 4.152 x [10.sup.-9], whereas for the Amano's method the accuracy is 2.036 x [10.sup.-8] .

By applying the Runge-Kutta method to the constrains equation, the approximate conformal mapping function is constructed by using new charge points and conformal mapping radius. In Figure 8 thick solid lines represent the contour boundaries; the thin solid line represents the contour lines. As can be seen from Figure 9, the approximate conformal mapping function F(z) maps the boundaries to concentric circles, and the inner region still corresponds to the interior region of concentric circles. This proves the effectiveness of the proposed method.

Example 3. Region D is the doubly connected domain with the boundaries [C.sub.1] : [absolute value of x] = 2, [absolute value of y] = 2 and [C.sub.2] : [absolute value of x] = 1, [absolute value of y] = 1 (see Figure 10).

For this example, the values of parameters p, SolRD, and ItMax are taken as 0.00002, [10.sup.-3], and 3000, respectively. Similar to Figure 6, Figure 10 shows the distribution of the charge points for N = 128. For Figure 11, we can find that the accuracy of numerical conformal mapping based on the Runge-Kutta method is better than that of Amano's method. When the number of charge points is more than 96, the error of Amanos method is obviously higher than that of our method. The error of numerical conformal mapping based on the Runge-Kutta method has reached 3.887 x [10.sup.-3] when the number of charge simulation points N = 144. Figures 12 and 13 represent the contour lines and corresponding mapping, respectively. Figure 13 shows that the method used in this paper is effective. As can be seen from Table 4, the CPU time of our method is shorter than that of Amanos method.

5. Conclusion

In this paper, a method has been proposed for improving the accuracy of the charge simulation method for numerical conformal mapping of doubly connected domain onto a circular annulus. The method constructs constraint equation by using the charge simulation method by Amano, and the constraint equations are solved by Runge-Kutta method. Numerical experiments showed that the proposed method can obtain higher numerical conformal mapping accuracy than Amano's method. The results of numerical conformal mapping of doubly connected domain are simulated by using the contour lines. In a forthcoming work, this method is considered for computing conformal mapping of multiply connected domains.

http://dx.doi.org/10.1155/2017/3603965

Competing Interests

The authors declare that they have no competing interests.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (NSFC) (Grant no. 11461037).

References

[1] P. Henrici, Applied and computational complex analysis. Vol. 3, vol. 3 of Pure and Applied Mathematics, John Wiley & Sons, New York, NY, USA, 1986.

[2] R. Schinzinger and P. A. Laura, Conformal Mapping: Methods and Applications, Elsevier, New York, NY, USA, 2003.

[3] G. T. Symm, "An integral equation method in conformal mapping," Numerische Mathematik, vol. 9, no. 3, pp. 250-258, 1966.

[4] G. T. Symm, "Conformal mapping of doubly-connected domains," Numerische Mathematik, vol. 13, pp. 448-457, 1969.

[5] W.-L. Lo, N.-J. Wu, C.-S. Chen, and T.-K. Tsay, "Exact boundary derivative formulation for numerical conformal mapping method," Mathematical Problems in Engineering, vol. 2016, Article ID 5072309, 18 pages, 2016.

[6] A. W. Sangawi, A. H. Murid, and M. M. Nasser, "Linear integral equations for conformal mapping of bounded multiply connected regions onto a disk with circular slits," Applied Mathematics and Computation, vol. 218, no. 5, pp. 2055-2068, 2011.

[7] A. W. Sangawi, A. H. Murid, and M. M. Nasser, "Annulus with circular slit map of bounded multiply connected regions via integral equation method," Bulletin of the Malaysian Mathematical Sciences Society. Second Series, vol. 35, no. 4, pp. 945-959, 2012.

[8] M. M. Nasser, "Numerical conformal mapping via a boundary integral equation with the generalized Neumann kernel," SIAM Journal on Scientific Computing, vol. 31, no. 3, pp. 1695-1715, 2009.

[9] M. M. Nasser, "Numerical conformal mapping of multiply connected regions onto the fifth category of Koebe's canonical slit regions," Journal of Mathematical Analysis and Applications, vol. 398, no. 2, pp. 729-743, 2013.

[10] N. H. Malik, "A review of the charge simulation method and its applications," IEEE Transactions on Electrical Insulation, vol. 24, no. 1, pp. 3-20, 1989.

[11] K. Amano, "Numerical conformal mappings of exterior domains based on the charge simulation method," Transactions-Information Processing Society of Japan, vol. 29, no. 1, pp. 62-72, 1998 (Japanese).

[12] K. Amano, "Numerical conformal mappings of interior domains based on the charge simulation method," Transactions-Information Processing Society of Japan, vol. 29, no. 7, pp. 697-699, 1988 (Japanese).

[13] K. Amano, "Numerical conformal mapping of doubly-connected domains based on the charge simulation method," Transactions-Information Processing Society of Japan, vol. 29, no. 7, pp. 914-919, 1988 (Japanese).

[14] R. S. Martin, G. Peters, and J. H. Wilkinson, "Symmetric decomposition of a positive definite matrix," Numerische Mathematik, vol. 7, no. 5, pp. 362-383, 1965.

[15] X. Wu, R. Shao, and Y. Zhu, "New iterative improvement of a solution for an ill-conditioned system of linear equations based on a linear dynamic system," Computers & Mathematics with Applications, vol. 44, no. 8-9, pp. 1109-1116, 2002.

[16] K. Amano, "A charge simulation method for the numerical conformal mapping of interior, exterior and doubly-connected domains," Journal of Computational and Applied Mathematics, vol. 53, no. 3, pp. 353-370, 1994.

[17] K. Amano, D. Okano, H. Ogata, and M. Sugihara, "Numerical conformal mappings onto the linear slit domain," Japan Journal of Industrial and Applied Mathematics, vol. 29, no. 2,pp. 165-186, 2012.

[18] H. Hakula, T. Quach, and A. Rasila, "Conjugate function method for numerical conformal mappings," Journal of Computational and Applied Mathematics, vol. 237, no. 1, pp. 340-353, 2013.

[19] Y.-B. Lu, D.-A. Wu, Y.-Z. Wang, and S.-S. Zheng, "The accuracy improvement of numerical conformal mapping using the modified gram-schmidt method," in Proceedings of the 19th International Conference on Industrial Engineering and Engineering Management: Assistive Technology of Industrial Engineering, pp. 555-563, Springer, October 2013.

[20] X. Wu and Y. Fang, "Wilkinson's iterative refinement of solution with automatic step-size control for linear system of equations," Applied Mathematics and Computation, vol. 193, no. 2, pp. 506-513, 2007

[21] D. G. Yakubu and A. M. Kwami, "Implicit two-derivative Runge-Kutta collocation methods for systems of initial value problems," Journal of the Nigerian Mathematical Society, vol. 34, no. 2, pp. 128-142, 2015.

[22] W. Tang and Y. Sun, "Construction of Runge-Kutta type methods for solving ordinary differential equations," Applied Mathematics and Computation, vol. 234, pp. 179-191, 2014.

[23] L. Brugnano and D. Trigiante, "Convergence and stability of boundary value methods for ordinary differential equations," Journal of Computational and Applied Mathematics, vol. 66, no. 1-2, pp. 97-109, 1996.

[24] K. J. In't Hout, "Convergence of Runge-Kutta methods for delay differential equations," BIT. Numerical Mathematics, vol. 41, no. 2, pp. 322-344, 2001.

[25] P. K. Kythe, Computational Conformal Mapping, Birkhauser, Boston, Mass, USA, 1998.

Fuming Lai, (1) Yingzi Wang, (2) Yibin Lu, (1) and Jian Wang (1)

(1) Faculty of Science, Kunming University of Science and Technology, Kunming 650500, China (2) Computer Center, Kunming University of Science and Technology, Kunming 650500, China

Correspondence should be addressed to Yibin Lu; luyibin@kmust.edu.cn

Received 8 September 2016; Revised 6 December 2016; Accepted 9 January 2017; Published 31 January 2017

Caption: Figure 1: Numerical conformal mapping of doubly connected domains.

Caption: Figure 2: Charge points for Example 1.

Caption: Figure 3: Error curves for Example 1.

Caption: Figure 4: The contour lines for Example 1.

Caption: Figure 5: The image is mapped for Example 1.

Caption: Figure 6: Charge points for Example 2.

Caption: Figure 7: Error curves for Example 2.

Caption: Figure 8: The contour lines for Example 2.

Caption: Figure 9: The image is mapped for Example 2.

Caption: Figure 10: Charge points for Example 3.

Caption: Figure 11: Error curves for Example 3.

Caption: Figure 12: The contour lines for Example 3.

Caption: Figure 13: The image is mapped for Example 3.
```Table 1: A comparison with Amano's method for Example 1.

z                     Exact /(z)

13.0000 + 0.000001    1.00000 + 0.000001
9.01722 + 8.731021    0.80902 + 0.587791
0.66312 + 11.27391    0.30902 + 0.951071
-5.51722 + 7.747211   -0.30902 + 0.951061
-7.16312 + 3.024681   -0.80902 + 0.587791
-7.00000 + 0.000001   -1.00000 + 0.000001
5.75000 + 0.000001    0.50000 + 0.000001
4.27685 - 3.652221    0.40451 - 0.293891
0.93832 - 5.196121    0.15451 - 0.475531
-4.25000 - 0.000001   -0.50000 - 0.000001

Amano's method        Our method

0.99984 + 0.018061    1.00000 + 0.001951
0.81851 + 0.574501    0.80891 + 0.587921
0.32181 + 0.946801    0.30435 + 0.952561
-0.32155 + 0.946891   -0.30465 + 0.952471
-0.81426 + 0.580511   -0.80848 + 0.588531
-0.99997 + 0.007971   -1.00000 + 0.002371
0.49994 + 0.007621    0.50000 + 0.000961
0.41374 - 0.280751    0.40398 - 0.294621
0.16536 - 0.471861    0.15353 - 0.475851
-0.49964 - 0.019021   -0.50000 - 0.002411

Table 2: Comparison of CUP time (s) for Example 1.

Charge points    N = 120   N = 140   N = 160   N = 180

Amano's method   0.3248    0.3660    0.5155    0.7028
Our method       0.0467    0.0878    0.0882    0.0930

Charge points    N = 200   N = 220   N = 240

Amano's method   1.0124    1.4303    1.6480
Our method       0.1028    0.1359    0.2202

Table 3: Comparison of CUP time (s) for Example 2.

Charge points    N = 110   N = 120   N = 130   N = 140

Amano's method   0.2528    0.3625    0.3071    0.3582
Our method       0.0162    0.0131    0.0215    0.0159

Charge points    N = 150   N = 160   N = 170

Amano's method   0.4211    0.5469    0.6525
Our method       0.0171    0.0191    0.0221

Table 4: Comparison of CUP time (s) for Example 3.

Charge points    N = 48   N = 64   N = 80   N = 96

Amano's method   0.0316   0.0758   0.1238   0.1600
Our method       0.0032   0.0111   0.0069   0.0147

Charge points    N = 112   N = 128   N = 144

Amano's method   0.2389    0.2920    0.3851
Our method       0.0045    0.0096    0.0118
```