Fast computation of analytic capacity

A boundary integral equation method is presented for fast computation of the analytic capacities of compact sets in the complex plane. The method is based on using the Kerzman--Stein integral equation to compute the Szeg\"o kernel and then the value of the Ahlfors map at the point at infinity. The proposed method can be used for domains with smooth and piecewise smooth boundaries. When combined with conformal mappings, the method can be used for compact slit sets. Several numerical examples are presented to demonstrate the efficiency of the proposed method. We recover some known exact results and corroborate the conjectural subadditivity property of analytic capacity.


Introduction
Capacities -such as analytic, logarithmic, and conformal -are fundamental tools in complex analysis and have several applications to problems in different fields, e.g., in approximation theory, potential theory, electronics, and fluid dynamics [5,6,8,19,31,37].These capacities can be expressed explicitly in only a handful of special cases, and therefore numerical methods are needed to compute these capacities in the majority of instances.
The Riemann mapping theorem states that any unbounded simply connected domain G in the extended complex plane C = C ∪ {∞} with ∞ ∈ G and whose boundary consists of more than one point can be mapped one-to-one onto the unit disk D by a conformal map f .If we assume that f (∞) = 0 and f ′ (∞) > 0, (1) then this mapping f is unique and known as the Riemann mapping function.Here, the derivative of analytic function f at the point at infinity is The so-called Ahlfors map can be regarded as an extension of the Riemann mapping function for multiply connected domains.That is, given an unbounded multiply connected domain G of connectivity m, the associated Ahlfors map is the unique analytic function function f that maps G onto D such that [8,31,37] f (∞) = 0, f ′ (∞) > 0, and f ′ (∞) is maximal. ( The Ahlfors map w = f (z) is then an m-to-one covering of G onto D which maps each boundary component of G one-to-one onto the unit circle.
Let us now introduce the notion of analytic capacity.Let E be a compact subset of the complex plane C and let G = C\E be its complement in the extended complex plane C = C ∪ {∞}.The analytic capacity of E is defined to be [31, p. 15] γ(E) = sup |f ′ (∞)| where the supremum is taken over all analytic functions f : G → C such that |f (z)| ≤ 1 for all z ∈ G, and f ′ (∞) is as in (2).In this paper, we assume that G is an unbounded multiply connected domain of connectivity m with ∞ ∈ G.
It is well-known that analytic capacity is inextricably linked to the Ahlfors map [6,8,19,31,37].If w = f (z) is the unique Ahlfors map from the unbounded multiply connected domain G in the z-plane onto the unit disk D in the w-plane satisfying the conditions (3), then the analytic capacity of the set E is given by [8,31] In particular, when E is compact and connected such that G = C\E is a simply connected domain, the Ahlfors map w = f (z) is a conformal map from G onto the unit disk D and hence the analytic capacity γ(E) is equal to the logarithmic capacity cap(E) of E. However, for a general compact set E, we have γ(E) ≤ cap(E).See [37, p. 13] for details.Closed-form expressions for analytic capacity are special and known only in a handful of cases, and it is informative to survey some of these briefly here.If E is a disk of radius r, then [31, p. 17] γ(E) = r, and if E is a square with sides of length ℓ, then [16] γ where Γ(•) is the gamma function.For a complex line segment E = [a, b] ⊂ C, we have [31, p. 17] γ(E) = 1 4 |b − a|.
For a compact and connected set E ⊂ C, we have where diam(E) denotes the diameter of E. For more details, see [31, p. 17].
When E and F are disjoint connected compact subsets of C, Suita [29] proved that the subadditivity property holds.For general compact sets, the proof of this property is still an open problem.However, Tolsa [30] proved the semi-additivity of analytic capacity: there exists a constant c such that for all compact sets E, F ⊂ C, the analytic capacity satisfies Moreover, proving the conjectural subadditivity of analytic capacity for arbitrary compact sets E, F ⊂ C is equivalent to proving it for all disjoint compact sets E, F ⊂ C that are finite unions of disjoint closed disks, all with the same radius [17,36].Several numerical examples have been considered by Younsi & Ransford in [36] for purely circular compact sets.Numerical results for sets other than circular ones were also presented in [36,35].All of these examples provide convincing evidence to suggest that the conjectural subadditivity property for analytic capacity is true.It should also be pointed out that, from (8), in the case of multiple real slits, equality as opposed to inequality holds in (9).The subadditivity property will be of significant interest in the ensuing discussion and will be corroborated numerically in several cases.Several numerical methods are available in the literature for computing the logarithmic and conformal capacities.One of these methods is based on the boundary integral equation (BIE) with the generalized Neumann kernel [16,24].For the numerical computation of analytic capacity, to the best of our knowledge, the only available numerical method is that given in [36].The method is based on using quadratic minimization for the numerical computation of upper and lower bounds for the analytic capacity, which, in principle, converges to its exact value.This method has been used before to compute the logarithmic capacity [27].In this paper, we present a fast and accurate BIE method for the numerical computation of the analytic capacity.Our presented method will be used to compute the analytic capacity for a wider class of compact sets, including those with smooth boundaries, piecewise smooth boundaries, and sets consisting of only slits.The method is based on using a BIE for the Szegö kernel (refer to [3,7,12,13,19] for the definition and basic properties of the Szegö kernel).The BIE has been used by Kerzman & Trummer [13] to compute the conformal mapping for simply connected domains.Bell [2] proved that the BIE can be used to compute the Ahlfors map of bounded multiply connected domains.In [4], the BIE has been implemented numerically to compute both the Riemann mapping function and the Ahlfors map for bounded simply and multiply connected domains, respectively.See also [3,12,25,32].
Besides this introductory section, our paper is structured in the following way.In Section 2, we introduce a numerical method for computing the analytic capacity for compact sets bounded by smooth or piecewise smooth Jordan curves, and several numerical examples for such sets are presented in Section 3. In Section 4, we consider compact slit sets.We provide our concluding remarks in Section 5. Finally, in Appendix A, for a given multiply connected rectilinear slit domain, we review an iterative method from [23] for the construction of a preimage unbounded multiply connected domain bounded by smooth Jordan curves.

The numerical method 2.1 Analytic capacity and the Szegö kernel
Let E be a compact subset of the complex plane C and let G = C\E be its complement in the extended complex plane C. We assume that G is an unbounded multiply connected domain bounded by m smooth, or piecewise smooth, Jordan curves Γ 1 , . . ., Γ m .Domains bounded by slits will be considered in Section 4 below.
From (4), the analytic capacity γ(E) is calculated by computing the derivative f ′ (∞) of the Ahlfors map f from the domain G onto the unit disk D with the normalization (1).
A BIE method for computing the Ahlfors map for a bounded multiply connected domain is presented in [2,4]; however, the domain G is unbounded, and therefore a preliminary step is required.We first conformally map the unbounded multiply connected domain G onto a bounded multiply connected domain D using the Möbius transformation where α is a point in the interior of any of the curves Γ j .The point at infinity is mapped onto the origin, M (∞) = 0 ∈ D. Let w = F (ζ) be the Ahlfors map from the bounded multiply connected domain D onto the unit disk D such that F (0) = 0, F ′ (0) > 0, and F ′ (0) is maximal.It follows immediately that is the Ahlfors map from the unbounded domain G onto the unit disk D which satisfies the conditions (3).Note that As D is a bounded multiply connected domain, the Ahlfors map w = F (ζ) from D onto the unit disk D can be computed using the method presented by Bell [2] (see also [4]).However, computing the analytic capacity γ(E) requires only computing the derivative F ′ (0) since γ(E) = f ′ (∞) = F ′ (0).In fact, it follows from [2] that where S(ζ, 0) is the Szegö kernel for the bounded multiply connected domain D with respect to the base point 0 ∈ D. Our numerical method is based on computing S(0, 0) using the BIE related to the Szegö kernel in multiply connected domains [2].Then

Integral equation for the Szegö kernel
Assume that each boundary component Γ j is parametrized by a 2π-periodic function ζ j (t), t ∈ J j = [0, 2π], j = 1, . . ., m.For domains with corners, the parametrization ζ j (t) is defined as described in [16].We define the total parameter domain J as the disjoint union of the m intervals J j = [0, 2π], j = 1, . . ., m.The whole boundary Γ is therefore parametrized by See [16,22,23] for more details.Further, the boundary ∂D of the bounded multiply connected domain D is parametrized by The Szegö kernel for the bounded multiply connected domain D with respect to the base point 0 ∈ D can be computed by solving the second-kind Fredholm integral equation [2,4] S(η(t), 0) where The kernel A(η(t), η(s)) is continuous with A(η(t), η(t)) = 0.The integral equation in (14) is known as the Kerzman-Stein BIE [12,13,25,32].It was proved in [2] that this BIE can be used also for the computation of the Szegö kernel for bounded multiply connected domains.Multiplying both sides of ( 14) by η ′ (t)/η(t) and defining the BIE ( 14) can be written as Now, since or equivalently, Taking the conjugate of both sides and then multiplying by i/|ζ ′ (t)|, we obtain which can be written in the concise form This BIE ( 17) is a modification of the Kerzman-Stein BIE (14).By (11), computing the analytic capacity γ(E) requires computing the value of the Szegö kernel S(0, 0).Since the Szegö kernel S(ζ, 0) is an analytic function in the domain D, by the Cauchy integral formula, we have Then, by solving the BIE (17) for ψ(t) and using ( 15) and ( 17), we have It then follows at once from (11) that It is immediate from (18) that the integral J ψ(t)|ζ ′ (t)|dt must be real and hence

Numerical solution of the integral equation
The Kerzman-Stein BIE (14) has been used to compute conformal mappings for bounded and unbounded simply connected domains [13,20,25,32], and in [2,4] to compute the Ahlfors map for bounded multiply connected domains.A combination of the usage of the Kerzman-Stein BIE ( 14) and the Fast Multipole Method (FMM) [9,10] has been presented in [25] for computing conformal mappings for bounded simply connected domains.
In this paper, to compute the analytic capacity γ(E), we will solve the BIE ( 17) which is a modified version of the Kerzman-Stein BIE (14).To accomplish this, we too shall employ the FFM when solving (17).Since the integrand in ( 17) is 2π-periodic, the BIE ( 17) can be best discretized by the Nyström method with the trapezoidal rule [13,20,25,32].
For domains with smooth boundaries, we use the trapezoidal rule with equidistant nodes.We discretize each interval J p = [0, 2π], for p = 1, 2, . . ., m, by n equidistant nodes s 1 , . . ., s n where and n is an even integer.Writing s = [s 1 , . . ., s n ], we discretize the parameter domain J by the vector t = [s, s, . . ., s] which consists of m copies of s, i.e., t = [t 1 , t 2 , . . ., t mn ] where for p = 1, 2, . . ., m and q = 1, . . ., n, For a real or a complex function µ(ζ(t)) defined on the boundary Γ, the trapezoidal rule then yields Discretizing the BIE (17) using the trapezoidal rule ( 21) and substituting t = t i , we obtain the linear system where ψ n is an approximation of ψ.
or equivalently, which can be written in the following concise form: Here, .* and ./are the MATLAB element-by-element multiplication and division, respectively, x = ψ n (t), r is the mn × 1 vector with entries r i = 1 2π , and B is the mn × mn matrix with entries The linear system (22) will be solved using the GMRES iterative method [28] where the matrix-vector product can be computed using the FMM.The FMM allows us to use the method presented to compute the analytic capacity of compact sets consisting of a high number of components.If we define the left-hand side of ( 22) to be a function of the unknown vector x, then the value of the function F(x) can be computed quickly and accurately using the MATLAB function zfmm2dpart in the MATLAB toolbox FMMLIB2D developed by Greengard & Gimbutas [9].This method for computing the solution ψ to the BIE ( 17) is summarized in the following MATLAB function where the tolerances for the FMM and the GMRES method are taken to be 0.5 × 10 −15 and 10 −14 , respectively, and the GMRES method is run without restart: The preceding method assumes that the boundaries of the domains of interest are smooth, i.e. without corners.In the case of domains with corners (excluding cusps), to obtain accurate results, a re-parametrization of the boundary of the domain is performed as described in [16].Assume that the boundary component Γ j has ℓ corner points.We first parametrize each boundary component Γ j by a 2π-periodic function ζj (t) for t ∈ J j = [0, 2π].The function ζ j (t) is assumed to be smooth with ζ ′ j (t) ̸ = 0 for all values of t ∈ J j such that ζ j (t) is not a corner point.We assume that ζ ′ j (t) has only the first kind discontinuity at these corner points.At each corner point, the left tangent vector is taken to be the tangent vector at this point.As above, let J be the disjoint union of the m intervals J j = [0, 2π], j = 1, 2, . . ., m and ζ(t), t ∈ J, be a parametrization of the whole boundary Γ.Then, we parametrize the boundary Γ by ζ(t) = ζ(δ(t)), where the function δ(t) is defined in [16, pp. 696-697].The function δ(t) is chosen such that the singularity in the first derivative of the solution of the BIE in the vicinity of the corner points is removed [15,16].With the new parametrization ζ(t), the BIE ( 17) can be solved accurately using the above MATLAB function.However, for domains with corners, we usually need a larger number of points n (which should be a multiple of the number of corners on each boundary component) for discretizing the BIE compared to domains with smooth boundaries (see [15,16,22] for further details).
Once the solution ψ of the BIE (17) has been found, we can proceed to compute the analytic capacity γ(E) using the formula (19).This can be undertaken using the following MATLAB function: Various numerical examples will be presented in the proceeding two sections.We will take, in turn, domains bounded by Jordan curves and domains bounded by slits.

Domains bounded by Jordan curves
In this section, we will use the method presented in the previous section to compute numerical approximations γ(E) to the analytic capacity γ(E) of compact sets bounded by smooth and piecewise smooth boundaries.
Example 1 Consider the compact set E = E 1 ∪ E 2 where E 1,2 = {z ∈ C | |z ± c| ≤ r} and 0 < r < c.Then [36] γ where is the second Jacobi theta function.The relative error in the computed approximate values γ(E) obtained with n = 2 9 are given in Table 1 for several values of c and r.For c = 2 and r = 1, our obtained value is γ(E) = 1.875595019097120 which is in the interval (1.875595019097112, 1.875595019097164) given in [36].The exact value is Example 2 Let E be the square with corners 1, −i, −1, i.In this case, the domain G = C\E is an unbounded simply connected domain.Thus, by (6), the analytic capacity of E is In this example, the boundary components of the compact set E have corners.We use our numerical method with various values of n to approximate the analytic capacity γ(E) and the results are presented in Table 2.As can clearly be seen in Table 2, the relative error decreases as n increases.Example 3 As a validation of our numerical method, let us also consider the four compact sets shown in Figure 1.These sets were considered in [36,Figures 2,3,5 & 6].We are not aware that the analytic capacities γ(E) of these sets are known analytically.The approximate values of γ(E) computed by our method are presented in Table 3.  10  2.12070613733366 0.57 It is important to point out that for the union of 4 ellipses, the elapsed time of 0.55 seconds suggests that our computations are more than 1000 times faster than the method used for the same problem in [35,36].However, around a decade has passed since the computations in [36] were performed.We also note that our method is much faster than the method presented in [35,36] for non-circular compact sets compared with circular ones.
Example 4 Consider the square [−2, 2] × [−2, 2] and the four sub-squares with centers ±1 ± i.Let ε > 0 be a real parameter.We consider three cases (i)-(iii).In the case (i), consider moving the centers of the sub-squares via the parameter ϵ to the points (1 + ε)(±1 ± i) (see Figure 2 (left)).In case (ii), let us fix the lower two sub-squares and Note that the length of each side of the original square is 4, and hence, by (6), Further, let us label the sub-squares by F 1 , . . ., F 4 , then the length of each of these subsquares is 2, and hence, by (6), In each of the three cases, we compute the analytic capacity γ(E ε ) as a function of ε which are presented in Figure 3.The results indicate that γ(E 0 ) is a lower bound for γ(E ε ) and 4 j=1 γ(F j ) = 2γ(E 0 ) is an upper bound for γ(E ε ).These results collectively provide numerical evidence to corroborate the conjectural subadditivity property of analytic capacity for these compact sets.Our results demonstrate the expected phenomenon that the analytic capacity increases as the sub-squares move further apart from each other, and values of the analytic capacity in the case (i) are larger than those in the cases (ii) and (iii) corresponding to the greater number of sub-squares.
Example 5 Consider the four disks with radius 1 and centers at ±1 ± i.Let ε > 0 be a real parameter.We consider three cases (i)-(iii).In case (i), consider moving the centers  of the four disks via the parameter ϵ to the points (1 + ε)(±1 ± i) (see Figure 4 (left)).In case (ii), let us fix the lower two disks and consider moving the centers of the upper disks to the points (1 + ε)(±1 + i) (see Figure 4 (middle)).Finally, in case (iii), let us fix three of the disks and move the center of the remaining disk to the point (1 + ε)(1 + i) by increasing the value of ϵ (see Figure 4 (right)).We label the union of the compact sets generated with ε in each case by E ε .
Let us label the disks by F 1 , . . ., F 4 , then γ(F j ) = 1, j = 1, . . ., 4. In each of the three cases, we compute the analytic capacity γ(E ε ) as a function of ε which are presented in Figure 5.The results demonstrate the expected phenomenon that the analytic capacity increases as the disks move further apart from each other, and values of the analytic capacity in case (i), corresponding to the case when all disks are moving away from each other, are larger than those in cases (ii) and (iii).
Remark 1 If E is a compact set such that G = C\E is bounded by finitely many analytic curves, then [35,36] where the minimum is taken over all functions g(z) that are analytic on G, continuous up to the boundary, with g(∞) = 0. Thus, taking g(z) ≡ 0 gives the inequality Hence, if E is a union of finitely many disjoint closed disks E 1 , . . ., E m with radii r 1 , . . ., r m , we have This explains the results presented in Figure 5 indicating that 4 j=1 γ(F j ) = 4 is an upper bound for γ(E ε ).

Example 6
We consider m = 100 random non-overlapping disks.For k = 1, 2, . . ., 100, the radius r k of the disk E k is chosen randomly in (0.2, 0.8) and its center c k is chosen in the square [−10, 10] × [−10, 10] such that all disks are non-overlapping.Then we randomly choose an integer ℓ ∈ [1,99].We define See Figure 6 (right) for an example of such compact sets E and F .
For this problem, we run our method 50 times to obtain 50 different locations for these disks as well as different sets E and F .For each run j, we use the above presented method with n = 2 9 to compute approximate values of the quantities γ(E), γ(F ) and γ(E ∪ F ).The values of the ratio are plotted as a function of the run number j as shown in Figure 6 (left).As can be seen from the graphs in this figure, we have verified that the conjectural subadditivity property of analytic capacity holds for each of the 50 random compact sets we considered, and in particular that γ(E ∪ F ) ≤ γ(E) + γ(F ).Again, the conjectural subadditivity property for analytic capacity holds for the compact sets in this example.

Example 7
In this example, we study compact sets consisting of disjoint disks of equal radii.In particular, we will validate numerically the conjectural subadditivity property of analytic capacity in several cases.The consideration of such compact sets is important since proving the conjectural subadditivity property of analytic capacity for arbitrary compact sets E, F ∈ C is equivalent to proving it for all disjoint compact sets that are finite unions of disjoint closed disks, all with the same radius [17,36].Let E m be a union of m disjoint disks and F be a disk with center x + iy such that these m + 1 disks are non-overlapping and of unit radii.Let the real function u(x, y) be defined by u(x, y) = γ(E m ∪ F ), i.e. the function u(x, y) is defined for all points (x, y) ∈ R 2 such that the disk F does not overlap any of the other m disks.Note that the domain of definition of the function u(x, y) is unbounded.Hence, in our numerical computations, we consider only the region −10 ≤ x, y ≤ 10 and we assume that the minimum distance between any two disks is 0.02.We compute approximate values of the function u(x, y) using n = 2 10 and then plot several level curves of u(x, y).We consider four cases of m, namely m = 1, 2, 3, 4. For m = 1, E m consists of only one disk which is assumed to be the unit disk.For m = 2, 3, 4, we assume that the centers of the disks forming The computed level curves of the function u(x, y) are presented in Figure 7.The results presented demonstrate that the values of the analytic capacity γ(E m ∪ F ) increase when F moves away from E m and decrease when F moves towards E m .In particular, when F is far away from E m , we have The results presented also illustrate that which indicates that the conjectural subadditivity property for analytic capacity again holds in this example.Finally, we compute approximate values γ(E m ) of the analytic capacity for the above compact sets E m , this time assuming that the centers of the disks are re 2kπi/m , k = 1, 2, . . ., m, for 1.5 ≤ r ≤ 20.The ratios γ(E m )/m are plotted in Figure 8 as functions of r.As we can see, each of these ratios tends to 1 as r increases; that is, the analytic capacity of E m (the union of m disks) tends to the sum of the individual analytic capacities of these m disks which is equal to m (see Remark 1 above).The results presented in (26) and Figures 7 and 8 provide experimental evidence for a previously known result due to Pommerenke who showed that, roughly speaking, the analytic capacity of a compact set is approximately equal to the sum of the analytic capacities of its components if they are far away from each other (see also [14, p. 267]).

Domains bounded by slits
In this section, we will use our numerical method to compute numerical approximations γ(E) to the analytic capacity γ(E) of sets consisting of slits.We will consider here only rectilinear slits.However, the method presented can be extended to other types of slits.
For a 1 , . . ., a m , b 1 , . . ., b m ∈ C such that E j = [a j , b j ] are non-overlapping for j = 1, 2, . . ., m, let E = m j=1 E j and let Ω = C\E, i.e., Ω is the unbounded multiply connected slit domain obtained by removing the m slits E 1 , . . ., E m from the extended complex plane C. The method presented in Section 2 is not directly applicable to such a domain Ω.However, an iterative method has been presented in [23] for constructing of a preimage unbounded multiply connected domain G bounded by smooth Jordan curves and the unique conformal mapping ζ = Φ(z), Φ : G → Ω, such that Φ is normalized near infinity by the condition The inverse function z = Φ −1 (ζ) is then the conformal mapping from Ω onto G.The method presented in Section 2 is now applicable to the new domain G.This iterative method is reviewed in Appendix A. Let w = f (z) be the Ahlfors map from the unbounded domain G onto the unit disk D such that f (∞) = 0, f ′ (∞) > 0, and f ′ (∞) is maximal.Then the function is an Ahlfors map from the unbounded slit domain Ω onto the unit disk D. We have where we used (27) and the fact that lim z→∞ zf (z) = f ′ (∞).Since f ′ (∞) is maximal, it then follows that g ′ (∞) is maximal too.Thus, by computing the preimage domain G and the conformal mapping Φ : G → Ω normalized by the condition ( 27), we will have where w = f (z) is the Ahlfors map from the unbounded domain G onto the unit disk D with the normalization f (∞) = 0 and f ′ (∞) > 0. Since the domain G is bounded by smooth Jordan curves, the value of f ′ (∞) can be computed as explained in Section 2. We consider three examples.In our first example, we consider rectilinear slits on the real line.We know the exact value of the analytic capacity in this case and hence the error in the computed approximate values can be calculated.Examples with unknown explicit formulae are also presented.

Example 8
We consider several of the compact sets used in the process of generating the middle-thirds Cantor set.Let E 0 = [−1, 1], and let E k for k ≥ 1 be defined recursively by Note that E k consists of m = 2 k sub-intervals of the [0, 1], each of length 2/3 k .We denote these sub-intervals by I j for j = 1, 2, . . ., m.By (8), the exact value of γ(E k ) is known and given by from which it is immediate that γ(E k ) → 0 as k → ∞.Note for this example that The proposed method is used to compute approximate values γ(E k ) to the analytic capacity γ(E k ) for k = 1, 2, . . ., 10 and the obtained results are presented in Table 4.We compute also the relative error in the computed approximate values.As can be seen in Table 4, our numerical method gives accurate results for sets consisting of a very high number of slits.Example 9 We next consider the union of two equal rectilinear slits of unit length: one slit F 1 = [0.1,1.1] is fixed on the real line, and the other is taken to be F 2,ε = [0.1eiεπ , 1.1e iεπ ] where we vary ε between zero and one.A schematic of this configuration is shown in Figure 9 (8).When ε = 0, we have E 0 = F 1 = F 2,0 = [0.1,1.1] and hence γ(E 0 ) = 1/4.For 0 < ε < 1, there is no exact value of γ(E ε ).We use our method to compute γ(E ε ) for 0.05 ≤ ε ≤ 1 and the numerical results are presented in Figure 9 That is, the value of the analytic capacity of E ε is maximum when the two slits are collinear.Furthermore, we always have Example 10 We next consider the union of four equal rectilinear slits of length 2 − 2ε, 0 < ε < 1, such that these four slits make the square [−1, 1] × [−1, 1] when ε = 0.A schematic of this configuration is shown in Figure 10 (right) when ε = 0.1.We denote these four slits by F k,ε with k = 1, 2, 3, 4. We define We use our method to approximate γ(E ε ) for 0.01 ≤ ε ≤ 0.99 and the obtained numerical results are presented in Figure 10 (left).As ε → 0, it is clear that the approximate values γ(E ε ) approach

Concluding remarks
This paper has shown how to use a numerical boundary integral equation method to quickly and accurately compute analytic capacity, an important conformal invariant.This quantity has been widely studied from a mostly theoretical perspective with several deep analytical results having been established [1,2,3,6,8,11,18,19,26,30,31,33,36,35,37].Analytic capacity is intimately connected to the Ahlfors map and Szegö kernel -two fundamental objects in complex analysis -and arises from the generalization of the Riemann map to multiply connected domains.
In our work, two particular classes of configurations were considered over which our calculations of analytic capacity were performed: compact sets bounded by smooth and piecewise smooth Jordan arcs, and domains consisting of a finite number of rectilinear slits.Throughout, we made connections with previous results; in particular, we have been able to corroborate the bounds found by Younsi & Ransford [36] for the analytic capacity for several compact sets they considered, and we were able to validate the conjectural subadditivity property of analytic capacity for numerous configurations.We were also able to validate numerically other exact results of analytic capacity, and illustrate several of its properties.Furthermore, the presented numerical results demonstrate that the analytic capacity γ(E ∪ F ), roughly speaking, increases when the distance between E and F increases and decreases when this distance decreases.
Our work has been mainly numerical in approach, and the key to its success lies in the BIE scheme based on the Kerzman-Stein BIE [2,13,20,32] and the FMM [9,10].The method presented can be used for domains with smooth and piecewise smooth boundaries as well as for domains with many boundary components.We used the method also to approximate the analytic capacity for compact sets consisting of rectilinear slits.However, for the latter case, a preliminary conformal mapping step is required; this has also been shown to be expedient in other works [16,23,24].
The MATLAB codes for the presented computations in this paper can be found at the link https://github.com/mmsnasser/ac.
The application of this method requires that the domain G is known.However, in our case, the slit domain Ω is known and the domain G is unknown and needs to be determined alongside the conformal mapping ζ = Φ(z) from G onto Ω.This preimage domain G as well as the conformal mapping Φ(z) will be computed using the following iterative method from [23] which generates a sequence of multiply connected domains G (0) , G (1) , G (2) , . . ., that converges numerically to the required preimage domain G. Let L j = |E j | be the length of the slit E j and let β j be its center, j = 1, . . ., m.In the iteration step i = 0, 1, 2, . .., we assume that G (i) is an unbounded multiply connected domain bounded by the m ellipses Γ (41) The parameters c (i) j and a (i) j are computed as follows: 1. Set c (0) j = (1 − 0.5r)L j , j = 1, . . ., m.
2. For i = 1, 2, 3, . . ., • Compute the conformal mapping from the preimage domain G (i−1) in the zplane onto the canonical rectilinear slit domain Ω (i) obtained by removing m rectilinear slits E  where ε is a given tolerance and Max is the maximum number of iterations allowed.
The above iterative method generates sequences of parameters c (i) j and a (i) j that converge numerically to c j and a j , respectively, and then the boundary components of the preimage domain G are parametrized by (30).In our numerical implementations, we used ε = 10 −13 and Max = 100.
It is clear that each iteration of the above method requires solving the integral equation with the generalized Neumann kernel (37) and computing the function h in (38) which can be done with the Matlab function fbie presented in [22].In our numerical computations, the value of n as well as the values of the other parameters in fbie are chosen to be the same as those used in the method described in Section 2 for computing the analytic capacity.The given slit domain Ω and the computed preimage domain G for Examples 9 and 10 are presented in Figures 11 and 12

Figure 2 :
Figure 2: The three kinds of compact sets considered in Example 4.

Figure 3 :
Figure 3: Graphs of the approximate values γ(E ε ) of the analytic capacity as a function of ε in each of the three cases in Example 4. The dashed line denotes the upper bound 2γ(E 0 ) and the dotted line denotes the lower bound γ(E 0 ).

Figure 4 :
Figure 4: The three kinds of compact sets considered in Example 4.

Figure 5 :
Figure 5: Graphs of the approximate values γ(E ε ) of the analytic capacity as a function of ε in each of the three cases in Example 5.The dotted line denotes the upper bound 4 j=1 γ(F j ) = 4.

Figure 6 :
Figure 6: The ratio in (25) as functions of the run number j for random radii in (0.2, 0.8) (left), and an example of the compact set E ∪F with 100 random disks, where the elements of E are shown in blue and the elements of F are shown in red (right).

Figure 7 :
Figure 7: The level curves of the function u(x, y) = γ(E m ∪ F ) in Example 7 when m = 1, 2, 4, 5, and the disks forming E m are colored blue.

Figure 8 :
Figure 8: The ratios γ(E m )/m as functions of the distance r measured from the origin to the centers of the m disks.

Figure 9 :
Figure 9: Graph of the approximate values γ(E ε ) of the analytic capacity for the compact set E ε in Example 3 as a function of ε (left).A schematic of the compact set E ε showing the parameter ε (right).

Figure 10 :
Figure 10: Graph of the approximate values γ(E ε ) of the analytic capacity for the union of rectilinear slits in Example 10 as a function of ε (left).A schematic of the compact set E ε for ε = 0.1 (right).The red dot is the exact value of the analytic capacity of the square [−1, 1] × [−1, 1] and the dashed-line is the upper bound 2 − 2ε of γ(E ε ).

Table 3 :
The approximate values γ(E) of the analytic capacity γ(E) for Example 3.
m.The approximate values γ(E m ) of the analytic capacity of E m are presented in the following table.

Table 4 :
(29)approximate values γ(E k ) to the analytic capacity γ(E k ) and their relative errors for Example 8, where γ(E k ) is given by(29).