A transform method for Laplace’s equation in multiply connected circular domains

A new general transform method for the solution of mixed boundary value problems (BVPs) for Laplace’s equation in multiply connected circular domains is formulated. Circular domains are those with boundaries that are a union of circular arcs and / or straight line segments. Constructive schemes based on the method are showcased by application to some concrete problems in ﬂuid mechanics. In particular, a BVP of topical interest in the study of superhydrophobic surfaces is solved and compared with an existing numerical solution.


Introduction
Consider the 2D domain D in an (x, y)-plane shown in Fig. 1: the region exterior to a circular disc in a channel of uniform width running parallel to the x-axis (a 'disc-in-channel' geometry). This geometry arises in many physical science applications and it is often required to solve mixed boundary value problems (BVPs) for various partial differential equations in such a domain. One such application, of current interest in the field of microfluidics, is treated in Section 8 where the so-called slip lengths over a class of superhydrophobic 'bubble mattresses' are found using the new techniques to be presented in this paper. Despite the simplicity and ubiquity of this domain, it is difficult to find in the literature any constructive analytical techniques to solve basic BVPs defined in it, although a number of attempts have been made. Richmond (1923), for example, derived a conformal mapping, in terms of elementary functions, to a domain that is close to the disc-in-channel geometry of Fig. 1 (the boundary of the 'disc' in his study is only approximately circular) and used it to solve the BVPs of interest to him (see also Martin & Dalrymple, 1988). Poritsky (1960) proposed a quasi-analytical approach for a problem of electrostatics in the disc-in-channel geometry based on the idea of expanding one natural basis set in terms of another. Often, purely numerical schemes are used such as boundary integral, domain decomposition or finite element methods (e.g. Teo & Khoo, 2010;Ng & Wang, 2011).
It is not even obvious how one might best represent a simple analytic function f (z) defined in D. If the circular disc is absent, Fourier transform representations of f (z) in a channel geometry come to mind; but then possible singularities of f (z) inside the disc must be accounted for. On the other hand, if just the disc is present (with no confining channel walls), one might think to represent f (z) using a Laurent series about the disc centre; but then its radius of convergence will be limited by the channel walls and will not give a uniformly valid representation in the whole of D. Intuitively, one ideally needs some amalgam, or hybrid, of these two ideas. The purpose of this paper is to present a new transform technique that is ideally suited to D, and to other, more general, circular geometries like it. The new method leads to uniformly valid (spectral) representations of analytic functions in the relevant domains as well as to concomitant global relations that can be analysed to solve BVPs posed in them. Circular domains are defined to be domains, possibly multiply connected, whose boundaries are a mixture of circular arcs and/or straight lines. The domain D in Fig. 1 qualifies as a doubly connected circular domain.
It is not usual to associate 'transform methods' with circular domains, let alone bounded ones. Traditionally, Fourier and Mellin transforms, and their variants (e.g. sine/cosine transforms), are applied in unbounded channel or wedge geometries having straight line boundaries with special techniques, such as Wiener-Hopf methods, invoked to tailor those transform methods to other geometries (e.g. semi-strips, or half-wedges).
However, in a recent paper Crowdy (2015) the author has presented what he argues are the natural analogues of Fourier-Mellin transforms for general circular domains. There, the focus is on simply connected circular domains of general type. It was also shown in that work how the new transforms are a natural generalization of the transform pairs for simply connected, convex polygonal domains first written down by Fokas & Kapaev (2003). As such, the new transforms of Crowdy (2015) fit neatly into the armoury of new transform techniques that have become known collectively as the unified transform method pioneered in recent years by Fokas (2008) (see also Fokas & Spence, 2012). The formulation in Crowdy (2015) appears to be the first to extend the Fokas unified transform technique to a PDE in general circular domains.
Going further, the present paper gives an extension of the new transform approach of Crowdy (2015) to the case of multiply connected circular domains (again, providing what appears to be the first extension of the unified transform method to domains with such non-trivial topology). Here we restrict attention to BVPs for Laplace's equation defined in such domains. The construction of the transform pairs is given in general in Sections 3 and 4, but the possibilities of the method are best showcased by dint of application to concrete problems as given in Sections 7 and 8.

Transform pairs for polygons
This section presents an elementary derivation of the transform method for convex polygons first introduced by Fokas & Kapaev (2003), who used other techniques. The material in this section is borrowed from Crowdy (2015), but this background is needed to formulate the multiply connected generalization.
We start with a basic geometrical observation. If z lies on some finite length slit on the real axis and z is in the upper-half plane, see  It follows trivially that It is easy to check that the contribution from the upper limit of integration vanishes for the particular choices of z and z to which we have restricted consideration. On the other hand, suppose that z lies on some other finite length slit making angle χ with the positive real axis and suppose that z is in the half-plane shown in Fig. 3 (the half-plane 'to the left' of the slit as one follows its tangent with uniform inclination angle χ ). Now the affine map (2.4) Fig. 4. A convex polygon P as an intersection of N = 3 half-planes with N angles {χ j | j = 1, 2, 3}. Formula (2.7) can be used in the Cauchy integral formula with χ = χ j when z is on side S j (for j = 1, 2, 3).
for example, where the (unimportant) constant α is shown in Fig. 3, takes the slit to the real axis, and z to the upper-half plane, and Hence, on use of (2.3) with the substitutions (2.4), we can write (2.6) or, on cancellation of α and rearrangement, This expression is valid uniformly for all z and z having the geometrical positioning depicted in Fig. 3. Now consider a bounded convex polygon P with N sides {S j | j = 1, . . . , N}. Figure 4 shows an example with N = 3. For a function f (z) analytic in P, Cauchy's integral formula provides that for z ∈ P, We can separate the boundary integral into a sum over the N sides: dz .
(2.9) 1906 D. CROWDY But if side S j has inclination χ j , then (2.7) can be used, with χ → χ j , to re-express the Cauchy kernel uniformly for all z ∈ P and z on the respective sides: (2.10) On reversing the order of integration, we can write where, for integers m, n between 1 and N, we define the spectral matrix to be (2.12) and where L = [0, ∞) is the fundamental contour for straight line edges.

The global relations
Observe that, for any k ∈ C, and for any m = 1, . . . , N, where we have used Cauchy's theorem and the fact that f (z ) e −i e −iχm kz (for m = 1, . . . , N) is analytic inside P. There are N such global relations relating different elements of the spectral matrix, but each is an equivalent statement of the analyticity of f (z) in the domain P.
While the global relations relate all elements of the spectral matrix, it is worth emphasizing that only the diagonal elements of the spectral matrix appear in the integral representation (2.11).
The formulation above differs slightly from that given in Fokas & Kapaev (2003) where the notions of spectral matrix and fundamental contour are not introduced. Note that, by defining, for each j = 1, . . . , N, the change of spectral variable given by in both the spectral functions and the integral representation then, rather than a sum of N contributions over the single fundamental contour in the spectral plane, it is easy to check that the transform pair (2.11) and (2.12) can be written as a sum of contributions from N (generally distinct) rays in the spectral plane: where L j = [0, ∞e −iχ j ) is the ray defined by arg[λ] = −χ j (these coincide with the definitions given in Fokas & Kapaev, 2003). Figure 5 shows the rays for the polygon shown in Fig. 4.  Fokas & Kapaev (2003) for the N = 3 polygon P shown in Fig. 4 (left). The latter rays are images of L under the N transformations (2.14).
Transform pairs for unbounded polygons, such as semi-strips, can also be derived (Fokas & Kapaev, 2003) with minor modifications to the derivation above (the only difference is that the global relations are now valid only in restricted parts of the spectral k-plane).

Transform pairs for circular domains
The following integral representation, established by the author in Crowdy (2015), is uniformly valid for |z| < 1: (3.1) The integrals are taken around the so-called fundamental contour for circular edges. The contour L 1 is the union of the negative imaginary axis (−i∞, −ir], where 0 < r < 1, and the arc of the quarter circle |k| = r in the third quadrant traversed in a clockwise sense; the contour L 2 is the real interval [−r, ∞); the contour L 3 is the arc of the quarter circle |k| = r in the second quadrant traversed in a clockwise sense together with the portion of the positive imaginary axis [ir, i∞). These contours are illustrated in Fig. 6. All integrals in (3.1) are non-singular and it is easy to verify directly that all integrands are exponentially decaying as |k| → ∞ uniformly for all |z| < 1.

Transform pair for interior of unit disc
Let z be a point on the unit circle and let |z| < 1 be a point inside the unit disc as depicted in Fig. 7. Then |z/z | < 1 and it follows on letting z → z/z in (3.1) that The Cauchy kernel is therefore given by the spectral representation  The Cauchy integral formula for some function f (z) analytic in the unit disc D is On substitution of the representation (3.3) for the Cauchy kernel, which is uniformly valid for the chosen z and z, we find (3.5) On swapping the order of integration, we arrive at the transform pair since, for this discrete set of k-values, the integrand of the integral defining ρ(k) is analytic in the unit disc.

Transform pair for exterior of unit disc
The transform pair for the exterior of the unit disc follows similarly. For z on the unit circle with |z| > 1 then |z /z| < 1 and setting z → z /z in (3.1) gives the identity (3.8) Hence, we have the spectral representation of the Cauchy kernel given by (3.9) For |z| > 1, the Cauchy integral formula holds for a function analytic outside the unit disc and is O(1/z) as z → ∞. Use of (the uniformly valid) expression (3.9) in the Cauchy integral formula then produces the following transform pair for such a function: (3.11) The global relation is since, for this discrete set of k-values, the integrand of the integral defining ρ(k) is analytic outside the disc and is O(1/z 2 ) as z → ∞.
If, instead, we are interested in the unbounded region exterior to the disc of radius q and centred at z = δ the modified transform pair is easily shown to be For future use, we prefer to write (3.13) with the addition of two minus signs: which clearly does not affect the transform pair. The reason for this is to preserve the convention that the spectral function is the integral around the boundary of the domain with that domain on the left as the boundary is traversed.

Transform pairs for multiply connected domains
What follows is new, and extends the ideas of Crowdy (2015). Consider now the multiply connected circular region D bounded by M + 1 circles {C j | j = 0, . . . , M } as shown schematically in Fig. 8 for where the integral around ∂D has been separated into the M + 1 integrals around the circular arcs The domain D is the intersection of the interior of circle C 0 and the exterior of the circles {C j | j = 1, . . . , M }. By substituting the relevant spectral representations of the Cauchy kernel for z on each boundary circle, and for z ∈ D, it can be deduced that the required integral representation is where the spectral matrix (Crowdy, 2015) is defined for n = 0, 1, , . . . , M by An important remark is in order. In contrast to the situation for simply connected domains treated in Crowdy (2015), the global relations (3.18) for different m = 0, . . . , M are not equivalent since they are associated with disconnected boundary components. This is a crucial distinction from the simply connected case where there is only a single boundary component. This means that all the relations 1912 D. CROWDY (3.18) must be analysed simultaneously to solve a given BVP. We will show explicitly how to go about this in the examples of Sections 7 and 8.

BVPs for Laplace's equation
Suppose that we seek the single-valued harmonic function φ satisfying (4.1) and some given set of boundary conditions on the boundary circles of a typical multiply connected circular domain D similar to that shown in Fig. 8. These conditions might be of Dirichlet, Neumann or Robin type, or some combination thereof, even including the possibility of several conditions of mixed type on the same boundary circle. For now, it is not necessary to specify the precise nature of these boundary conditions. For Laplace's equation in a multiply connected domain D, we can still write for some function g(z) that is analytic in D, but the non-trivial topology of the domain now means that g(z) is not necessarily single-valued in D. But the single-valuedness of φ implies that the derivative g (z) = dg/dz is both analytic and single-valued in D [see Crowdy (2008) for a discussion of these matters with respect to the so-called Schwarz problem in multiply connected domains]. Therefore, from earlier results, we can write the following transform representation for g (z): where the spectral matrix is defined for n = 0, 1, , . . . , M by A transform representation for g(z) now follows on integration with respect to z: (4.6) The corresponding spectral matrix (4.4) can also be re-expressed in terms of g(z) (instead of g (z)) on use of integration by parts.
It is interesting to remark that the integration with respect to z has increased by one the order of the singularity of the integrand at k = 0 in the M spectral integrals associated with the interior holes in formula (4.6) (but not, on the other hand, in the integral associated with the outer boundary). Removing these M additional singularities would clearly require the conditions to hold. On use of the definitions of the spectral functions in (4.4), this is equivalent to the conditions where square brackets denote the change in the contained quantity on an anticlockwise traversal of the subscripted contour. The conditions (4.8), in turn, are clearly a statement that g(z) is single-valued in the multiply connected domain D (in which case we can, from the outset, start with a transform representation of g(z) having the form (4.3) originally assumed only to hold for g (z) which is always known, a priori, to be single-valued). For generic BVPs, however, conditions (4.7) will not pertain and it is easy to verify that the second order pole at k = 0 forces logarithmic singularities of g(z) associated with each hole, as expected from more general considerations of BVPs in multiply connected domains (Gakhov, 1990;Crowdy, 2008).
In the following sections, we will proceed under the assumption that solving for the functions in the spectral matrix will lead, on substitution of that spectral data into the integral representations for the relevant analytic functions, to the solution of the stated BVP. That this is true is not obvious, but it has been established for BVPs for Laplace equations in convex polygons by Fulton et al. (2004) and for general elliptic equations in arbitrary convex domains by Ashton (2012).

Remark on multiply connected polygons
Before thinking of multiply connected circular domains one might ask: what is the generalization of the results of Fokas & Kapaev (2003) to multiply connected polygons (having just straight line edges)?
The answer is as follows. We have just exploited the fact that a multiply connected circular region can be understood as the intersection of the interior of one outer circle (denoted by C 0 ) and the exterior of M smaller interior circles (denoted by {C j | j = 1, . . . , M }). A multiply connected polygonal region can be defined similarly: it is the intersection of the interior of some outer polygonal boundary ∂P 0 , 1914 D. CROWDY  (2015), for the unbounded region exterior to some bounded polygon (e.g. one of the P j 's) and enclosing the point at infinity, it is geometrically impossible for such a domain to be an intersection of a finite collection of half-planes. Rather, such an unbounded polygonal region can be split up into an atlas of convex polygons that are the intersection of half-planes; this was done, for example, by Charalambopoulos et al. (2010) in the case of the exterior of an equilateral triangle. Consequently, to formulate transform pairs for any multiply connected polygonal region, it is similarly necessary to split that multiply connected domain into a finite union of convex polygonal subdomains in each of which a different transform representation is valid (with continuity imposed across common edges). Figure 9 shows an example doubly connected domain P together with one way (there are many) in which it can be decomposed into four convex polygonal regions {P j | j = 1, . . . , 4} in which various transform expressions can be written down. This is easily done, but we will not do so here. Continuity of φ must be imposed over the common edges of distinct polygons in this atlas as done in the analysis, using the unified transform method, of some problems of Wiener-Hopf type in Crowdy & Luca (2014) (although there the splitting of the domain was done for different reasons).
The same feature just described is not necessarily true, however, for multiply connected circular domains, as we have already seen for domains of the kind shown in Fig. 8. Moreover, it is even possible to find circular domains having a mixture of circular and straight line edges for which the need to split up the domain into smaller subdomains is unnecessary. Domain D of Fig. 1 is one such example.
By now, it should be clear to the reader how to combine the ideas of Section 2 (for straight line edges) and Section 3 (for multiply connected circular-arc edges) to produce hybrid transform pairs for multiply connected domains involving both circular and straight line edges. Essentially, starting from Cauchy's integral formula, for each edge of the domain one substitutes the appropriate spectral representations of the Cauchy kernel, swaps the order of integration and arrives at the relevant transform pairs.

Modified Schwarz problem in multiply connected circular domains
In this section, we briefly indicate how to use the new transform approach to solve a very basic BVP: the classic (modified) Schwarz problem in a multiply connected circular domain. This problem has already been considered by the author in Crowdy (2008) where quite different special function theory methods are employed. Here we tackle the same problem within the new transform approach.
Suppose that f (z) is analytic, and single-valued, in the multiply connected circular region D comprising the unit disc with M smaller discs excised. C 0 denotes the unit circle, {C j , j = 1, . . . , M } denote the other interior circular boundaries. The centre of C j is at δ j ∈ C, its radius is q j ∈ R. Suppose that, on C j for j = 0, 1, . . . , M , where r j (z,z) is a given function. It should be noted that, in order for the single-valued function f (z) to exist, the given data r j (z,z) must satisfy certain compatibility conditions; these are precisely the conditions (4.7) discussed earlier. We present the details below for the case M = 2 (triply connected) shown in Fig. 8. The generalization to higher M will be obvious.
Since f (z) is analytic and single-valued in D, the global relations (3.18), which hold for all k ∈ −N, are We will use these global relations to find the unknown real parts {φ j (z,z) | j = 0, 1, . . . , M } of f (z) on the boundaries. On substitution of (6.1) into (6.2), we find, for k ∈ −N, The latter are all known functions. To find the {φ j (z,z) | j = 0, 1, . . . , M }, we decompose each as φ 0 (z,z) = a 0 + n 1 a n z n +ā n z n , These are representations of the most general real functions on the respective boundaries. It is clear from the BVP that f (z) will only be determined up to an additive real constant so, without loss of generality, we can set, say, a 0 = 0. (6.6) On substitution of (6.5) into (6.4), we arrive at a linear system for the unknown coefficients in (6.5); indeed, the elements of the matrix system can be computed explicitly. The resulting system can be solved by evaluating (6.4) at as many values of k ∈ −N as one likes to form an overdetermined system that can be solved by a least-squares method. With f (z) now determined on the boundaries, Cauchy's integral formula can be used to evaluate f (z) at interior points of D.
In the special case where D is a concentric annulus one takes M = 1 and sets δ 1 = 0 and q 1 = ρ for some 0 < ρ < 1. The global relations (6.2) are then which, when combined, are equivalent to (6.8) where ∂D denotes the entire boundary of D. This is clearly just a statement of the analyticity of the single-valued function f (z) in the annulus D.

Uniform flow past two cylinders
A more interesting example is to consider a problem of classical ideal (inviscid) fluid mechanics. We seek the solution for the steady irrotational uniform flow, parallel to the real axis, past two unit radius cylinders centred at (0, ±h) where h > 1 so that the fluid domain is doubly connected. Figure 10 shows a schematic. Let z = x + iy and let the complex potential for the flow be where φ is the velocity potential and ψ is the streamfunction. There are two ways to proceed with this problem. One option is to exploit the transform pair for the doubly connected exterior of the two circular discs. A second option is to exploit the symmetry of the configuration with respect to reflection in the x-axis and to confine attention to the flow in the region D defined as the upper-half plane exterior to the cylinder centred at ih: y > 0 with |z − ih| > 1 as shown in Fig. 10. Since the latter option involves a combination of circular and straight line edges, and is therefore deemed to be more instructive, we adopt the second approach.
As z → ∞, we require that dw/dz = 1 hence The streamline conditions (no penetration of the fluid into the cylinder, or across the symmetry line) are for some real constant c. It is assumed that the circulation around each cylinder is zero, which is required for consistency with the assumed symmetry of the flow. This requires that w(z) is single-valued in D. Now let where we will use a transform method to findŵ(z) which is analytic and single-valued in D. The boundary conditions onŵ(z), following from (

Transform representation
Sinceŵ(z) is single-valued in D, we can write the following integral transform representation for z ∈ D: (7.6) where L = [0, ∞) is the fundamental contour for straight line edges and the set {L j | j = 1, 2, 3} constitutes the fundamental contour for circular edges. The two spectral functions in this formula are Other elements of the spectral matrix are The global relations are then ρ 11 (k) + ρ 12 (k) = 0, k < 0 (7.9) and ρ 21 (k) + ρ 22 (k) = 0, k ∈ −N. (7.10) Owing to the doubly connected nature of the problem domain these global relations are independent and both must be analysed in finding the unknown boundary data.

Analysis of global relations
To analyse the global relations, we will writê 1 1 ) for some as yet unknown φ(x), and The function φ(x), the coefficients a 0 ∈ R, {a m ∈ C} and c ∈ R are to be found. The global relation (7.9) implies that where we define and On taking a complex conjugate of (7.13), and letting k → −k, we find But from the inverse Fourier transform for the upper-half plane (Fokas & Kapaev, 2003) (or see Section 2) we can write (7.17) Now (7.17) can be substituted into the second global relation (7.10) to produce where (7.19) and, for n 2, where we have carried out some standard exercises in residue calculus. On use of (7.12), the right-hand side of (7.18) is where we define n 3.
(7.22) Equation (7.18) therefore leads to the system of equations: (7.23) Sinceŵ(z) is only determined up to a real constant, we choose a 0 = 0 then, because the right-hand side of (7.23) is independent of c, the first equation in (7.23) determines c once the set {a 1 , a 2 , . . .} has been found. From (7.23) with n → n + 1, the latter coefficients are solutions of the linear system 2π ia n = ∞ 0 I(k, n + 1)[B 1 (−k) +R 1 (−k)] dk − R 2 (n + 1), n 1. (7.24) Hence, (7.24) is for n 1. On use of (7.14) and (7.15) this is equivalent to with the latter equality following from a repeated integration by parts.
Finally, the linear system to be solved for {a n | n = 1, 2, . . .} is given explicitly by This is easily solved. With {a n | n = 1, 2, . . .} and c determined from the solution of (7.28) and the first equation in (7.23) thenŵ(z) on |z − ih| = 1 is now known so ρ 22 (k) can be found. The spectral function ρ 11 (k) follows from (7.16).

Alternative solution via conformal mapping
The solution to this problem can also be obtained by completely different methods using a new calculus for ideal flows in multiply connected domains devised by the author (Crowdy, 2010). An implicit parametric form of the solution W (ζ ) ≡ w(z) is given by (7.29) for ζ in the annulus ρ < |ζ | < 1. The special function K(ζ , ρ) is defined by (1 − ρ 2n ζ )(1 − ρ 2n /ζ ), (7.30) where the prime notation denotes differentiation with respect to ζ . The function P(ζ , ρ) is essentially the so-called Schottky-Klein prime function for the annulus ρ < |ζ | < 1. This alternative solution gives a valuable check on the transform solution.

Superhydrophobic bubble mattresses
An important problem, of topical interest to the microfluidics community, is the problem of slow viscous longitudinal shear flow over a superhydrophobic 'bubble mattress' (Crowdy, 2011). In contrast to the inviscid flow problem of Section 7, this problem involves a highly viscous fluid (or, zero Reynolds number flow). We omit a full description of the physical problem and focus only on giving a mathematical statement of the governing BVP. Further details of the physical problem, and related literature, can be found in Crowdy (2011).
In the absence of a pressure gradient, the axial velocity parallel to the grooves of the bubble mattress is harmonic, i.e. The flow takes place in x > 0 but the problem can be extended to the whole channel region exterior to a circular inclusion.
The required velocity field is the imaginary part, the function χ is its harmonic conjugate. The flow takes place in the region x > 0, −l < y < l and for |z| > 1 (outside the unit disc) with l > 1 (Fig. 11). This domain represents a single period window of a periodic arrangement of bubbles. The boundary conditions are that and w = 0, x = 0, 1 |y| l. (8.4) By the symmetry of the flow, we also deduce that The complex potential h(z) can be taken to be single-valued in the flow region. This completes the statement of the problem.
Let This allows us to extend the region of analyticity of H(z) into the entire channel region outside the unit disc. Henceforth, we seek the solution for H(z) in the full channel domain exterior to the circular inclusion shown in Fig. 11. The boundary conditions can now be written as which implies, on substitution of (8.6), that on |z| = 1, (8.14) On side 1 wherez = z + 2li the boundary condition (8.11) implies (8.15) which, by the periodicity of the coth function and use of (8.8), yields On side 3 wherez = z − 2li the boundary condition (8.11) implies  (8.19) where L 1 = [0, ∞) and L 2 = [0, −∞) and {L j | j = 1, 2, 3} are the fundamental contours for circular edges. (Here we have elected to use the spectral rays for straight line edges used by Fokas & Kapaev (2003) rather than the fundamental contour L for straight line edges, this choice is up to the user.) The three spectral functions are defined to be and ρ 31 (k) = ρ 11 (k) and ρ 13 (k) = ρ 33 (k). The global relations are then (8.24) which are equivalent, and We must analyse both global relations (8.24) and (8.25) in order to find the unknown boundary data.

Analysis of global relations
The global relation (8.24) can be written as  (8.30) and where the real coefficients {b n ∈ R} are to be found. The ansatz (8.29) satisfies both (8.8) and (8.10). Equation (8.26) can be written as where is recognized as the standard Fourier transform ofH(x) and It follows that From an inverse Fourier transform, we find (8.36) It is clear that we must require (8.47) The resulting system of equations is (8.49) This linear system (8.48), together with (8.37), must be solved for {b n | n 1} and the slip length λ. This is straightforward, and the matrix A nm is found to have a condition number close to unity. Once the coefficients {b n | n 1} and the slip length λ are known, then H(z) on |z| = 1 is known and ρ 22 (k) can be readily computed from (8.29) and (8.30). The other required spectral functions follow from (8.26).

Comparison with full numerical solution
For convenience, we introduce the parameter a = 1 l , (8.50) which, for l > 1, lies in the interval 0 < a < 1. A graph of λ/(2a) as computed from the transform method just described is shown in Fig. 12. Superposed on this figure, as crosses, are some data points as computed by a numerical method of Teo & Khoo (2010). The agreement is excellent, and provides verification of the transform analysis. Some interesting observations arise from this calculation. The dashed line in Fig. 12 is a graph of an approximation to the slip length λ given by the simple formula λ approx = π a (2 − π 2 a 2 /6) . (8.51) As a → 0, (8.51) gives λ approx ≈ π a 2 , (8.52) which is consistent with the result of Crowdy (2011) (for protrusion angle θ = π/2 in the notation of that paper) who studied the dilute limit of this problem (i.e. the limit l → ∞ or a → 0). This dilute limit is shown as the dotted line in Fig. 12. Formula (8.51) can therefore be viewed as an (approximate) extension of the result of Crowdy (2011) to the non-dilute case (a not necessarily small). The approximate result (8.51), which appears to be very effective, was produced from the following considerations. On inspection of the numerical results it appears to be reasonable, over a large range of l, to set B(0) to zero in (8.37) and approximate it as Figure 12 clearly shows that divergence between the numerically computed λ and λ approx is only significant as a → 1 and, for a ≈ 0.9, it remains a reasonable approximation only slightly underestimating the actual value. Equation (8.51) is therefore a simple and accurate approximation to the longitudinal slip length provided the bubbles are not too close together. Ng & Wang (2011) proposed their own analytical approximation formula for λ in the form of a polynomial fitted to match their numerical data computed using a triangulated-grid finite element solver: λ Ng+Wang = −0.40 + 8.93a − 46.74a 2 + 134.28a 3 − 173.88a 4 + 86.56a 5 ; (8.56) this is Ng & Wang (2011, Equation (6)) divided by a to match the present scalings. On evaluation, it is found that (8.56) gives numerical results that are very close to those given by (8.51) over the range a ∈ [0, 0.85].

Discussion
A new transform method for multiply connected circular domains, including those with boundaries that are a mixture of circular arcs and straight lines, has been presented. We have given two concrete applications of the method to problems in fluid mechanics. The method is very general and applies to a broad class of BVPs for Laplace's equation in multiply connected domains. Figure 13 shows some other common geometries, often arising in applications, where the method also applies. For clarity of exposition, we have restricted here to cases in which each of the M boundary components is either a full circle, or straight line through the point at infinity, but it should be clear how to adapt the approach to cases in which the jth connected boundary component, say, is subdivided into n j circular arc and/or straight line segments. How to do this with N-sided, simply connected polycircular arc domains was explained in Crowdy (2015).
We have also focussed on presenting the method for BVPs in which each disconnected boundary component features a boundary condition of a single type (e.g. Dirichlet, Neumann or Robin). However,