1 Introduction

This note is the result of investigations into an open uniqueness question for quadrature domains in the complex plane \(\mathbf {C}\), which appears in papers such as [21, 24, 44]. After describing the problem and reviewing some known results, we will suggest and explore an approach based on methods from real algebraic geometry and symbolic computation.

To get started, it is convenient to fix some notation that will be used throughout.

General notation. By a “domain” \(\Omega \) we mean an open and connected subset of \(\mathbf {C}\); we write \({\overline{\Omega }}\) for its closure, \(\partial \Omega \) for its boundary, and \(\Omega ^e=\mathbf {C}\setminus {\overline{\Omega }}\) for its exterior. A bounded domain \(\Omega \) is said to be “solid” if \(\Omega ^e\) is connected and \(\partial \Omega =\partial (\Omega ^e)\). Note that what we call solid domains are also referred to as “Carathéodory domains” in some sources.

We write “dA” for the normalized Lebesgue measure in the plane \(dA(z)=\frac{1}{\pi }\, dx\, dy\) (so the unit disc has area 1).

We denote by \(L^1(\Omega )\) the usual \(L^1\)-space of functions on \(\Omega \) that are integrable with respect to dA, and we write \(AL^1(\Omega )\), \(HL^1(\Omega )\), and \(SL^1(\Omega )\) for the subsets of \(L^1(\Omega )\) consisting of all analytic, harmonic, and subharmonic functions, respectively.

  • Standard sets: \(D(a,r)=\{z\in \mathbf {C};\,|z-a|<r\}\), \(\mathbf {D}=D(0,1)\), \(\mathbf {T}=\partial \mathbf {D}\).

  • Differential operators: \(\partial =\frac{1}{2}(\partial _x-i\partial _y)\), \(\bar{\partial }=\frac{1}{2} (\partial _x+i\partial _y)\), \(\Delta =4\partial \bar{\partial }=\partial _x^2+\partial _y^2.\)

Given an open set \(\Omega \subset \mathbf {C}\), a subspace (or cone) \({\mathcal {F}}\subset L^1(\Omega )\), and a linear functional \(\mu :{\mathcal {F}}\rightarrow \mathbf {C}\), we consider quadrature identities of the form

$$\begin{aligned} \int _\Omega f\, dA=\mu (f),\qquad f\in {\mathcal {F}}. \end{aligned}$$
(1)

The linear functional \(\mu \) will always be a fixed measure or distribution of appropriate type with compact support in \(\Omega \) (and defined on the appropriate test-class \({\mathcal {F}}\)). Since we will frequently take \(\mu \) to be a combination of point-evaluations, we stress that statements like \({\mathcal {F}}\subset L^1(\Omega )\) should not be taken literally, but rather in terms of the natural injective maps \({\mathcal {F}}\rightarrow L^1(\Omega )\).

If the above conditions hold, we say that \(\Omega \) is a quadrature domain (or “q.d.”) with data \((\mu ,{\mathcal {F}})\), and we write \(\Omega \in Q(\mu ,{\mathcal {F}})\). We are mainly interested in the case \({\mathcal {F}}=AL^1(\Omega )\), but also \({\mathcal {F}}=HL^1(\Omega )\), \({\mathcal {F}}=SL^1(\Omega )\) will play a role. In the last case, (1) must be replaced by the inequality

$$\begin{aligned} \mu (f)\le \int _\Omega f\, dA,\qquad f\in SL^1(\Omega ). \end{aligned}$$

It is easy to see that \(Q(\mu ,SL^1)\subset Q(\mu ,HL^1)\subset Q(\mu ,AL^1)\) and that a solid domain belongs simultaneously to the classes \(Q(\mu ,HL^1)\) and \(Q(\mu ,AL^1)\).

The above classes are conveniently interpreted in terms of the logarithmic potentials

$$\begin{aligned} U^\mu :=\ell *\mu ,\qquad U^\Omega :=U^{\mathbf {1}_\Omega \, dA},\qquad \text {where}\quad \ell (z):=\frac{1}{2}\log \frac{1}{|z|}. \end{aligned}$$

For example, we have that \(\Omega \in Q(\mu ,AL^1)\) if and only if \(\partial U^\mu =\partial U^\Omega \) on \(\Omega ^e\) and \(\Omega \in Q(\mu ,SL^1)\) if and only if \(U^\mu =U^\Omega \) on \(\Omega ^e\) and \(U^\mu \ge U^\Omega \) on \(\mathbf {C}\).

Given these proviso, we can formulate our basic problem in a succinct way (cf. [21]).

(Q). Determine whether or not there exists a functional \(\mu \) such that the class \(Q(\mu ,AL^1)\) contains two distinct, solid domains.

Theorem 1.1

The following uniqueness results are known.

  1. (i)

    If \(\Omega _1,\Omega _2\in Q(\mu ,AL^1)\) are star-shaped with respect to a common point, then \(\Omega _1=\Omega _2\).

  2. (ii)

    If there exists a solid domain \(\Omega \in Q(\mu ,SL^1)\), then this \(\Omega \) is the unique solid quadrature domain, even within the class \(Q(\mu ,AL^1)\).

  3. (iii)

    If \(\mu \) is a positive measure of total mass m and if \({\text {supp}}\mu \) is contained in a disc of radius r where \(r^2< m\), then each solid domain \(\Omega \) of class \(Q(\mu ,AL^1)\) is obtainable from \(\mu \) by partial balayage, and so it belongs to \(Q(\mu ,SL^1)\).

Remark on the proof

Part (i) is due to Novikov [33], cf. also [17, 29]. Part (ii) was proved in Sakai’s book [36], using the technique of partial balayage. An alternative proof is found in the paper [21] by Gustafsson. The statement (iii) was likewise proved by Sakai using partial balayage, see the papers [27, 35, 38]. \(\square \)

The “classical” setting corresponds to point-functionals, i.e., functionals \(\mu \) of the form

$$\begin{aligned} \mu (f)=\sum _{i=1}^m\sum _{j=0}^{n_i-1} c_{ij}f^{(j)}(a_i),\qquad f\in AL^1(\Omega ) \end{aligned}$$

where \(a_i\) are some points in \(\Omega \) and \(c_{ij}\) some complex numbers. A quadrature domain of this type is said to be of order \(n_1+\cdots +n_m\). When \(\mu \) contains no derivatives, i.e., when \(\mu (f)=\sum _{i=1}^n c_if(a_i)\) we speak of a pure point-functional. (Cf. [15].)

Example 1.1

Theorem 1.1(ii) completely settles the uniqueness problem for subharmonic quadrature domains. The following example due to Gustafsson shows that the question (Q) for analytic test functions is of a different kind.

It is shown in [20, Sect.  4] that there exists a quadrature domain \(\Omega \) having the appearance in Fig. 1, satisfying a three-point identity \(\int _\Omega f\, dA=c_1f(a_1)+c_2f(a_2)+c_3f(a_3)\) where \(c_1,c_2,c_3>0\).

It is known (see e.g. [27, Thm. 2.1]) that if \(\Omega \in Q(\mu ,SL^1)\) and if \(\partial \Omega \) has cusps, then those cusps must be contained in the convex hull of the nodes \(a_i\). Hence the quadrature domain in Fig. 1 is not subharmonic.

Fig. 1
figure 1

Gustafsson’s example of a domain in \(Q(\mu ,HL^1)\setminus Q(\mu ,SL^1)\). The triangle depicts the convex hull of the nodes \(a_1,a_2,a_3\)

Now fix a solid quadrature domain \(\Omega \in Q(\mu ,AL^1)\) containing the origin 0. Let \({\varphi }:\mathbf {D}\rightarrow \Omega \) be the conformal map normalized by \({\varphi }(0)=0\) and \({\varphi }'(0)>0\). Recall that \(\Omega \) is uniquely determined by \({\varphi }\) via Riemann’s mapping theorem.

The following theorem, which gives a non-trivial relation for \({\varphi }\), will be the main tool in our subsequent investigations.

Theorem 1.2

Suppose that \(\mu \) has compact support in \(\Omega \) and let \(\nu \) be the pullback of \(\mu \) to \(\mathbf {D}\), i.e., \(\nu (g)=\mu (g\circ {\varphi }^{-1})\) for \(g\in AL^1(\mathbf {D})\). Then \(\varphi \) satisfies the relation

$$\begin{aligned} \varphi (z)= \nu _\lambda ^*\left( \frac{z}{\overline{\varphi '(\lambda )}(1-z{\bar{\lambda }})}\right) ,\qquad z\in \mathbf {D}, \end{aligned}$$
(2)

where \(\nu _\lambda ^*(g):=\overline{\nu _\lambda ({\bar{g}})}\) acts on integrable anti-analytic functions \(g(\lambda )\).

Conversely, if \({\varphi }\) is any univalent solution to (2), then the domain \(\Omega ={\varphi }(\mathbf {D})\) is of class \(Q(\mu ,AL^1)\).

This result is not very easy to spot in the literature, but it has in fact been noticed earlier in somewhat different guises. The first proof might be due to Davis, see [15, Ch. 14], cf. also [14, Sect. 5]. Since the result will be central for what follows, we include an alternative proof (that we have found independently) in Sect. 2.

In the special case of a pure point-functional \(\mu (f)=\sum _1^n c_i f(a_i)\), the relation (2) takes the form

$$\begin{aligned} {\varphi }(z)=\sum _{i=1}^n \frac{{\bar{c}}_i}{{\bar{w}}_i}\frac{z}{1-{\bar{\lambda }}_i z},\qquad w_i={\varphi }'(\lambda _i),\, {\varphi }(\lambda _i)=a_i, \end{aligned}$$
(3)

which appears implicitly in e.g. the books [25, 40, 43].

The main idea behind our approach is to “solve” functional relations such as (3) by using techniques from algebraic geometry. To set up a suitable system of polynomial equations we differentiate (3) and substitute \(z=\lambda _j\), giving

$$\begin{aligned} w_j=\sum _{i=1}^n\frac{{\bar{c}}_i}{{\bar{w}}_i}\frac{1}{(1-\lambda _j{\bar{\lambda }}_i)^2},\qquad a_j=\sum _{i=1}^n\frac{{\bar{c}}_i}{{\bar{w}}_i}\frac{1}{1-\lambda _j{\bar{\lambda }}_i},\qquad j=1,\ldots ,n, \end{aligned}$$
(4)

where the unknown complex numbers \(\lambda _i\) and \(w_i\) are subject to the constraints

$$\begin{aligned} \lambda _1=0,\quad |\lambda _2|<1,\quad \ldots ,\quad |\lambda _n|<1,\quad w_1>0. \end{aligned}$$
(5)

The appearance of inequalities and complex-conjugates means that we are considering the real semi-algebraic geometry of a particular system of rational functions. Such semi-algebraic systems, i.e. those having the special structure of (4), (5) have, to the best of our knowledge, not been systematically studied before.

It is of course possible that no univalent solution to (4), (5) exists; for example the quadrature identity \(\int _\Omega f\, dA=f(a_1)+f(a_2)\) implies that \(\Omega \) is the disjoint union \(D(a_1,1)\cup D(a_2,1)\) if \(|a_1-a_2|\ge 2\). However, after having studied exact solutions for many examples of lower order quadrature domains, we find “empirically” the pattern that the system tends to have at most one solution which may or may not give rise to a univalent mapping \({\varphi }\). The non-univalent solutions fail to be locally univalent, i.e., they (still, empirically) satisfy \({\varphi }'=0\) somewhere in the disc \(\mathbf {D}\). For a different type of quadrature domains, not described by (4) and (5), a similar observation was made by Ullemar in [42].

Remark

The non-univalent solutions \({\varphi }\) are believed to represent quadrature domains on Riemann surfaces with branch points. Such q.d.’s are studied in the references [23, 37, 41].

Remark

Note that our method relies on knowledge of all solutions to (4), (5). To find one or several approximate solutions, one can of course try to apply numerical methods, such as Newton’s iterative method. By appropriately choosing different initial data, we may indeed obtain solutions by such methods in a relatively short time for n up to 10. However, since the number of solutions to the system is unknown, it is impossible to know when one has found all solutions, so this kind of information is of no use when studying the uniqueness question for quadrature domains. Another problem with a numerical approach is that systems such as (4), (5) tend to be quite sensitive to small perturbations of the quadrature data \(\{c_i,a_i\}_1^n\).

Remark

In connection with uniqueness problems in the gravi-equivalent sense, it is pertinent to point (besides the sources already mentioned) to the early works of Zidarov and Zhelev in the context of geophysics, see e.g. [45, 46]. (We thank one of the referees for this remark.)

The literature on quadrature domains is vast, and we have at this stage omitted to mention several important aspects. A somewhat fuller picture is given in Sect. 7, where we briefly compare a few of the more well-known techniques that have been developed over the years, such as Laplacian growth and Schottky-Klein functions.

To illustrate the challenges involved in studying the uniqueness of quadrature domains, we now give an example demonstrating the subtlety of the problem even for a q.d of order 2.

Example 1.2

Let \(\Omega _1\) be the solid quadrature domain obtained from a monopole with charge 1/2 and a dipole with strength \(\sqrt{3}/18\) placed at the origin, i.e.

$$\begin{aligned} \int _{\Omega _1}f\ dA=\frac{1}{2}f(0)+\frac{\sqrt{3}}{18}f'(0),\qquad f\in AL^1(\Omega _1). \end{aligned}$$

For quadrature domains of this type, Aharonov and Shapiro have proved uniqueness in [2]; in fact \(\Omega _1\) is determined as the image of \(\mathbf {D}\) under the conformal map \(p(z)={\sqrt{3}}(2z+z^2)/6\). The boundary \(\partial \Omega _1\) is a cardioid with a cusp at \(p(-1)={-1}/{\sqrt{12}}\), see Fig. 2.

Let us now construct a similar q.d. \({\tilde{\Omega }}_2\) but only using point charges,

$$\begin{aligned} \int _{{\tilde{\Omega }}_2}f\ dA=-\frac{1}{2}f(0)+f(a_2),\qquad f\in AL^1({\tilde{\Omega }}_2). \end{aligned}$$

Clearly both \(\Omega _1\) and \({\tilde{\Omega }}_2\) have area 1/2. We shall choose the parameter \(a_2\) real, such that the boundary of \({\tilde{\Omega }}_2\) has a cusp.

Since the parameters \(c_i,w_i,\lambda _i\) in (3) must be real in this case, the conformal map \({\tilde{\varphi }}:\mathbf {D}\rightarrow {\tilde{\Omega }}_2\) takes the form

$$\begin{aligned} {\tilde{\varphi }}:z\mapsto \frac{c_1z}{w_1}+\frac{c_2}{w_2}\frac{z}{1-\lambda _2z} \end{aligned}$$
(6)

where \(c_1=-1/2\) and \(c_2=1\), \(\lambda _2\) the pre-image of \(a_2\), and \(w_1={\tilde{\varphi }}'(0),\ w_2={\tilde{\varphi }}'(\lambda _1)\).

Writing \(\nu =-9+18\,i\sqrt{2}\), we find that there is a unique choice of \(a_2\) producing a cusp, namely

$$\begin{aligned} a_2=&\frac{1}{2\,\nu ^{1/3}}\,\sqrt{-\nu ^{1/3} \left( 3\,i\sqrt{3}\,\nu ^{2/3}-27\,i\sqrt{3}+3\, \nu ^{2/3}+4\,\nu ^{1/3}+27\right) } \approx 0.1316. \end{aligned}$$

The remaining parameters in the mapping function may be obtained by solving the system corresponding to (6) using the method of “RCTD” described in Sect. 4. The result is

$$\begin{aligned} \lambda _2&=\frac{3a_2}{\sqrt{3a_2^2+3}} \approx 0.2261,\\ w_{{2}}&={\frac{-c_{{2}}{\lambda _{{2}}}^{3} \left( -{\lambda _{{2}}}^{2 }c_{{1}}+{\lambda _{{2}}}^{2}c_{{2}}+{a_{{2}}}^{2} \right) }{a_{{2}} \left( -c_{{1}}{\lambda _{{2}}}^{6}+c_{{2}}{\lambda _{{2}}}^{6}+{ \lambda _{{2}}}^{4}{a_{{2}}}^{2}+2\,c_{{1}}{\lambda _{{2}}}^{4}-2\,{ \lambda _{{2}}}^{2}{a_{{2}}}^{2}-{\lambda _{{2}}}^{2}c_{{1}}-{\lambda _{{2}}}^{2}c_{{2}}+{a_{{2}}}^{2} \right) }} \approx 0.6674, \\ w_1&={\frac{c_{{1}}\lambda _{{2}} \left( {\lambda _{{2}}}^{2}-1 \right) w_{{ 2}}}{w_{{2}}{\lambda _{{2}}}^{2}a_{{2}}-w_{{2}}a_{{2}}+\lambda _{{2}}c_{ {2}}}} \approx 0.5016. \end{aligned}$$

This gives a cusp at

$$\begin{aligned} {\tilde{\varphi }}(-1)=-\frac{c_1}{w_1}-\frac{c_2}{w_2}\frac{1}{1+\lambda _2} \approx -0.2253. \end{aligned}$$
Fig. 2
figure 2

The domains \(\Omega _1\) and \(\Omega _2\)

A translate of \({\tilde{\Omega }}_2\) by \(\alpha :=p(-1)-{\tilde{\varphi }}(-1)\approx -0.06333\) leads to a domain \(\Omega _2={\tilde{\Omega }}_2+\alpha \) having a cusp at the point \(p(-1)\) and satisfying the quadrature identity

$$\begin{aligned} \int _{\Omega _2}f\ dA=-\frac{1}{2}f(\alpha )+f(a_2+\alpha ),\qquad f\in AL^1(\Omega _2). \end{aligned}$$

The resemblance between \(\Omega _1\) and \(\Omega _2\) (Fig. 2) is striking, even though they admit completely different quadrature identities. The similarity between \(\Omega _1\) and \(\Omega _2\) indicates that their potentials should be similar, and in terms of numerical values they are. But there is one essential difference between the two: the potential of \(\Omega _1\) is exactly determined by two terms in its multipole expansion while the potential for \(\Omega _2\) needs the entire infinite series. In detail we have

$$\begin{aligned} U^{\Omega _1}(z)&=\frac{\frac{1}{2}}{2\pi }\log \frac{1}{|z|}+\frac{\frac{\sqrt{3}}{18}}{2\pi }\mathrm {Re}\left( \frac{1}{z}\right) \\ U^{\Omega _2}(z)&\approx \frac{\frac{1}{2}}{2\pi }\log \frac{1}{|z|}+\frac{0.09998}{2\pi }\mathrm {Re}\left( \frac{1}{z}\right) +\cdots \end{aligned}$$

and for comparison note \(\sqrt{3}/18\approx 0.09623\). From this example we see that two very similar domains may have fundamentally different potentials.

2 The Master Formula

As previously stated, Theorem 1.2 appears (in equivalent form) in [15, Eq. (14.12)]. We shall here give a different derivation.

Consider the univalent map \({\varphi }:\mathbf {D}\rightarrow \Omega \) normalized by \({\varphi }(0)=0\) and \({\varphi }'(0)>0\) where \(\Omega \in Q(\mu ,AL^1)\) and where the distribution \(\mu \) is assumed to be of compact support in \(\Omega \). Our point of departure is Poisson’s equation

$$\begin{aligned} \Delta U^\mu =-\mu . \end{aligned}$$

Taking the distributional \(\partial \)-derivative of \(-4U^\mu \) we obtain the Cauchy transform

$$\begin{aligned} \mathcal {C}\mu (z):=\mu *k(z)=\mu (k_z), \end{aligned}$$

where \(k(\lambda ):=\lambda ^{-1}\) denotes the Cauchy kernel and \(k_z(\lambda ):=k(z-\lambda )\). Since \(-4\partial U^\mu ={\mathcal {C}}\mu \), (2) says that \(\bar{\partial }{\mathcal {C}}\mu =\mu \). In particular \({\mathcal {C}}\mu \) is holomorphic on \({\mathbb {C}}\setminus {\text {supp}}\mu \). Taking \(f=k_z\) \((z\not \in \Omega )\) in the quadrature identity (1) we see that

$$\begin{aligned} \partial U^\Omega =\partial U^\mu \qquad \text {on}\qquad \mathbf {C}\setminus \Omega . \end{aligned}$$
(7)

In fact, an application of Bers’ approximation theorem from [4] shows that (7) is equivalent to (1). Now consider the “Schwarz potential” u defined by

$$\begin{aligned} u=\mathcal {C}(\mathbf {1}_\Omega -\mu )=-4\partial (U^\Omega -U^\mu ), \end{aligned}$$

which is zero for \(z\in \Omega ^e\). Using the continuity of \(\partial U^\Omega \) we get \(u=0\) also on \(\partial \Omega \). Moreover, Poisson’s equation gives \(\bar{\partial }u=\mathbf {1}_\Omega \) on \(\mathbf {C}\setminus {\text {supp}}\mu \). Hence the function

$$\begin{aligned} S(z):={\bar{z}}-u(z) \end{aligned}$$

is holomorphic on \(\Omega \setminus \mathrm {supp}\ \mu \) and continuous up to the boundary \(\partial \Omega \), while satisfying \(S(z)={\bar{z}}\) for \(z\in \partial \Omega \). This determines S as the Schwarz function for the boundary curve \(\partial \Omega \), cf. [15, 40].

The following lemma is well-known, see e.g. [15, 40].

Lemma 2.1

The conformal mapping \(\varphi :\mathbf {D}\rightarrow \Omega \) extends holomorphically across \(\mathbf {T}\) to an analytic function on the disk D(0, R) for some \(R>1\).

Proof

As we saw above, the function \(S\circ \varphi \) is defined and holomorphic in some annulus \(1-\epsilon<|z|<1\), continuous up to the boundary and satisfies \(S(\varphi (z))=\overline{\varphi (z)}\) when \(z\in \mathbf {T}\). Likewise, the function \(\varphi ^*(z)=\overline{\varphi ({\overline{z}})}\) is holomorphic in \(\mathbf {D}\) and continuous up to the boundary, and we have the relation

$$\begin{aligned} S(\varphi (z))=\varphi ^*(1/z),\qquad z\in \mathbf {T}. \end{aligned}$$

Now, \(\varphi ^*(1/z)\) is holomorphic in the exterior of \(\mathbf {D}\), so the above formula shows that the functions \(S\circ \varphi (z)\) and \(\varphi ^*(1/z)\) are analytic continuations of each other across the circle \(\mathbf {T}\). In particular, \(\varphi ^*(1/z)\) is analytically continuable inwards across \(\mathbf {T}\), which means that \(\varphi ^*\) as well as \(\varphi \) are analytically continuable outwards across \(\mathbf {T}\), to some disc D(0, R) with \(R>1\). \(\square \)

Lemma 2.2

We have that \( \mathcal {C}\left[ \overline{\varphi '}\cdot \mathbf {1}_\mathbf {D}\right] (z)=\overline{\varphi (z)} ,\, z\in \mathbf {T}.\)

Remark

For an absolutely continuous measure \(\nu =f\, dA\), we prefer to denote its Cauchy transform by \(\mathcal {C}f\) rather than \(\mathcal {C}\nu \).

Proof of Lemma 2.2

Fix a point \(z=e^{i\theta }\in \mathbf {T}\) and a positive number \(r<1\) and put

$$\begin{aligned} I_r=\,&\mathcal {C}\left[ \overline{\varphi '}\cdot \mathbf {1}_\mathbf {D}\right] (rz)=\int _\mathbf {D}\frac{\overline{\varphi '(\lambda )}}{re^{i\theta }-\lambda }\ dA(\lambda ) = re^{-i\theta }\int _{r\mathbf {D}}\frac{\overline{\varphi '(r\zeta e^{i\theta })}}{1-\zeta }\ dA(\zeta ). \end{aligned}$$

Set \( \varphi (\lambda )=\sum _{j=0}^\infty c_j\lambda ^j\), where \(\sum _{j=0}^\infty |c_j|<\infty \) by Lemma 2.1. Inserting the expansion \(1/(1-\zeta )=\sum \zeta ^j\) we find that

$$\begin{aligned} I_r =\sum _{j=1}^\infty {\bar{c}}_j\left( re^{-i\theta }\right) ^jr^{2j}. \end{aligned}$$

Since \(\sum |c_j|<\infty \) we may pass to the limit as \(r\nearrow 1\), leading to

$$\begin{aligned} \lim _{r\rightarrow 1}I_r&=\lim _{r\rightarrow 1}\sum _{j=1}^\infty {\bar{c}}_j\left( re^{-i\theta }\right) ^jr^{2j}=\sum _{j=1}^\infty {\bar{c}}_j\left( e^{-i\theta }\right) ^j=\overline{\varphi (z)}. \end{aligned}$$

The proof of the lemma is complete. \(\square \)

Proof of Theorem 1.2

The quadrature identity (1) pulls back to

$$\begin{aligned} \int _\Omega f\ dA=\int _\mathbf {D}(f\circ \varphi )\cdot |\varphi '|^2\ dA=\nu (f\circ \varphi ) \end{aligned}$$
(8)

where \(\nu (g)=\mu (g\circ \varphi ^{-1})\).

Given an arbitrary \(f\in AL^1(\Omega )\) we define a function \(g\in AL^1(\mathbf {D})\) by

$$\begin{aligned} g=(f\circ \varphi )\cdot \varphi '. \end{aligned}$$

The identity (8) can be written as

$$\begin{aligned} \int _\mathbf {D}g\overline{\varphi '}\ dA=\nu \left( \frac{g}{\varphi '}\right) . \end{aligned}$$
(9)

Now fix a point \(z\in \mathbf {T}\) and choose g to be the Cauchy-kernel \(g=k_z\). With this choice, (9) takes the form

$$\begin{aligned} \mathcal {C}\left[ \overline{\varphi '}\cdot \mathbf {1}_\mathbf {D}\right] (z)=\mathcal {C}\left[ \frac{1}{{\varphi }'}\cdot \nu \right] (z) ,\qquad z\in \mathbf {T}. \end{aligned}$$

By Lemma 2.2 this is equivalent to

$$\begin{aligned} \overline{\varphi (z)}=\mathcal {C}\left[ \frac{1}{\varphi '}\cdot \nu \right] (z)=\nu _\lambda \left( \frac{1}{\varphi '(\lambda )(1/{\bar{z}}-\lambda )}\right) ,\qquad z\in \mathbf {T}. \end{aligned}$$

Taking complex-conjugates and considering the analytic continuation to \(\mathbf {D}\) we obtain

$$\begin{aligned} {\varphi }(z)=\nu _\lambda ^*\left[ \frac{z}{\overline{\varphi '(\lambda )}(1-z{\bar{\lambda }})}\right] ,\qquad z\in \mathbf {D},\end{aligned}$$
(10)

as desired.

Conversely, if \({\varphi }\) is univalent (and normalized) in \(\mathbf {D}\) and satisfies (10), we may read backwards and deduce (8), so \(\Omega ={\varphi }(\mathbf {D})\) belongs to \(Q(\mu ,AL^1)\) where \(\mu \) is the push-forward of \(\nu \). \(\square \)

We conclude this section with three examples of applications of Theorem 1.2, which are known from the literature on quadrature domains.

Example 2.1

Let us compare Theorem 1.2 with the computations in Shapiro’s book [40, Prop. 3.2]. For this purpose we fix a pure point functional \(\mu =\sum _{i=1}^n c_i\delta _{a_i}\), which pulls back to

$$\begin{aligned} \nu (g)=\mu \left( g\circ \varphi ^{-1}\right) =\sum _{i=1}^n c_ig\left( \varphi ^{-1}(a_i)\right) =\sum _{i=1}^n c_ig(\lambda _i),\qquad \lambda _i={\varphi }^{-1}(a_i). \end{aligned}$$

By Theorem 1.2 we know that an arbitrary solid domain \(\Omega \in Q(\mu ,AL^1)\) is of the form \(\Omega ={\varphi }(\mathbf {D})\) where \({\varphi }\) is univalent and normalized and satisfies

$$\begin{aligned} \varphi (z) =\sum _{i=1}^n\frac{{\overline{c}}_i}{\overline{\varphi '(\lambda _i)}}\frac{z}{1-z{\overline{\lambda }}_i} \end{aligned}$$
(11)

Given such a \({\varphi }\), we put \(S(\varphi (z)):=\varphi ^*(1/z)\) and note that S is the Schwarz function for \(\partial \Omega \). In view of (11), S is a meromorphic function with simple poles at \(z=a_j\). As \(z\rightarrow \lambda _j\) the dominant term in \(S\circ \varphi \) satisfies

$$\begin{aligned} S(\varphi (z))\sim \frac{c_j}{\varphi '(\lambda _j)}\frac{1}{z-\lambda _j}. \end{aligned}$$

From this we get that (as \(z\rightarrow a_j\))

$$\begin{aligned} S(z)\sim \frac{c_j}{\varphi '(\lambda _j)}\frac{1}{\varphi ^{-1}(z)-\lambda _j}=\frac{c_j}{\varphi '(\lambda _j)}\frac{z-a_j}{\varphi ^{-1}(z)-\varphi ^{-1}(a_j)}\frac{1}{z-a_j}\sim \frac{c_j}{z-a_j}. \end{aligned}$$

The residues of S are thus just Res\((S;a_j)=c_j\) for all j. Since \(S(z)={\bar{z}}\) on \(\partial \Omega \), an application of Green’s theorem and the Residue theorem now gives

$$\begin{aligned} \int _\Omega f\ dA&=\frac{1}{2\pi i}\int _\Gamma f(z){\overline{z}}\ dz=\frac{1}{2\pi i}\int _\Gamma f(z)S(z)\ dz = \sum c_jf(a_j) \end{aligned}$$

where \(\Gamma \) is the positively oriented boundary of \(\Omega \). We have shown again that \(\Omega \in Q(\mu ,AL^1)\).

We remark that a similar proof applied to a more general point functional \(\mu (f)=\sum c_{ij}f^{(j)}(a_i)\) gives the well known result (see [26, Thm. 3.3.1]) that a solid domain \(\Omega \) is a quadrature domain of finite order if and only if each conformal map \({\varphi }:\mathbf {D}\rightarrow \Omega \) is a rational function, if and only if the Schwarz function of \(\partial \Omega \) extends to a meromorphic function in \(\Omega \).

Example 2.2

Let \(\mu \) be a linear combination of a monopole and a dipole at the origin, i.e. \(\mu (f)=M_0f(0)+M_1f'(0)\). The action of the pullback is then given by

$$\begin{aligned} \nu (g) =\mu (g\circ \varphi ^{-1}) =M_0g(0)+\frac{M_1}{\varphi '(0)}g'(0),\qquad g\in AL^1(\mathbf {D}). \end{aligned}$$

Applying Theorem 1.2, we find that a normalized conformal map \({\varphi }:\mathbf {D}\rightarrow \Omega \), where \(\Omega \in Q(\mu ,AL^1)\), must satisfy

$$\begin{aligned} \overline{\varphi (z)}=\left( \frac{1}{\varphi '(0)}M_0- \frac{\varphi ''(0)}{\varphi '(0)^3}M_1\right) {\bar{z}} +\frac{M_1}{\varphi '(0)^2}{\bar{z}}^2. \end{aligned}$$

Hence \({\varphi }(z)\) is a polynomial of degree two. To determine this polynomial, we need to determine the derivatives \({\varphi }'(0)\), \({\varphi }''(0)\). The computation is postponed to Subsect. 5.1, after we have discussed some algebraic prerequisites.

Example 2.3

Following Davis [15, pp. 162-166] we now take \(\mu \) be a line charge with linear density h on the segment \([-a,a]\) of the x-axis, i.e.,

$$\begin{aligned} \displaystyle \mu (f)=\int _{-a}^af(x)h(x)\ dx. \end{aligned}$$

Suppose that \(\Omega ={\varphi }(\mathbf {D})\in Q(\mu ,AL^1)\); the pullback \(\nu \) by \({\varphi }\) is then given by

$$\begin{aligned} \nu (g)= & {} \mu (g\circ \varphi ^{-1})=\int _{-a}^{a}g(\varphi ^{-1}(x))h(x)\ dx=\int _{-\lambda }^\lambda g(w)h({\varphi }(w)){\varphi }'(w)\, dw,\\&\qquad \lambda =\varphi ^{-1}(a). \end{aligned}$$

Applying Theorem 1.2 we see that \({\varphi }\) must satisfy the functional equation (cf. [15, Eq. (14.25)])

$$\begin{aligned} \overline{\varphi (z)}={\bar{z}}\int _{-\lambda }^\lambda \frac{h({\varphi }(w))}{1-{\bar{z}}w}\, dw. \end{aligned}$$

In particular, if we specialize to a uniform line charge \(h\equiv 1\), we obtain

$$\begin{aligned} {\varphi }(z)=\log \left( \frac{1+z\lambda }{1-z\lambda }\right) . \end{aligned}$$

This relation was found by Davis (see Eq. (14.13)), where it is also shown that if \(a=1\) then \(\lambda =\sqrt{\tanh (1/2)}=.6798\cdots <1\). It is then easy to see that the above map \({\varphi }\) is well-defined by choosing the standard branch of the logarithm. We have shown that there exists a unique solid domain \(\Omega \) in the class \(Q(\mathbf {1}_{[-1,1]}(x)\, dx,AL^1)\); a picture is given in Fig. 3.

Fig. 3
figure 3

The q.d. \(\Omega \) generated by the linear density \(d\mu (x)=\mathbf {1}_{[-1,1]}(x)\, dx\)

3 Schur–Cohn’s Test

As stated earlier, the mapping problem (3) may have non-univalent solutions \({\varphi }\). The Schur–Cohn test provides a convenient way of discarding such \({\varphi }\) that fail to be locally univalent. One might hope that this procedure would leave us with at most one univalent \({\varphi }\), thus settling the uniqueness problem. As we will see later, this is indeed the case for a large class of examples.

The Schur–Cohn test is well known and quite elementary, see [28]. For reasons of completeness, we have found it convenient to briefly review the main ideas behind it here. In what follows, we let \(\mathcal {P}_n\) denote the space of all polynomials p of degree at most n,

$$\begin{aligned} p(z)=a_0+a_1z+\cdots +a_nz^n,\qquad a_0,\ldots ,a_n\in {\mathbb {C}}. \end{aligned}$$

Lemma 3.1

If \(n\ge 1\) then

  1. (i)

    if \(|a_0|>|a_1|+\cdots +|a_n|\) then \( p(z)\ne 0\) for all z with \(|z|\le 1\),

  2. (ii)

    if \( p(z)\ne 0\ \mathrm {for\ all}\ z\ \mathrm {with}\ |z|\le 1\) then \(|a_0|>|a_n|\).

Proof

  1. (i)

    If \(|a_0|>|a_1|+\cdots +|a_n|\) then \(|a_1z+\cdots +a_nz^n|\le |a_1|+\cdots +|a_n|<|a_0|\) for \(|z|\le 1\), so \(|p(z)|\ge |a_0|-|a_1z+\cdots +a_nz^n|>0\).

  2. (ii)

    Assume \(p(z)\ne 0\) for all z with \(|z|\le 1\). Then \(p(0)=a_0\ne 0\) which leads to two cases: (1) if \(a_n=0\) then obviously \(|a_0|>|a_n|\); (2) if \(a_n\ne 0\) we have \(p(z)=a_n(z-z_1)\cdots (z-z_n)\) where \(z_n\) are the zeros of p. But then \(p(0)=a_0=a_n(-z_1)\cdots (-z_n)\) and thus \(|a_0|>|a_n|\) since \(|z_k|>1\) for \(k=1,\ldots ,n\). \(\square \)

For each \(p\in \mathcal {P}_n\) the reciprocal polynomial \(p^\#\in \mathcal {P}_n\) is defined by

$$\begin{aligned} p^\#(z)=z^n\cdot \overline{p(1/{\overline{z}})}={\overline{a}}_n+{\overline{a}}_{n-1}z+\cdots +{\overline{a}}_1z^{n-1}+{\overline{a}}_0z^n. \end{aligned}$$

Let us now define the Schur transform, \(S_n:\ \mathcal {P}_n\rightarrow \mathcal {P}_{n-1}\) by

$$\begin{aligned} (S_np)(z)={\overline{a}}_0p(z)-a_np^\#(z)=\sum _{k=0}^{n-1}({\overline{a}}_0a_k-a_n{\overline{a}}_{n-k})z^k,\qquad p\in \mathcal {P}_n. \end{aligned}$$

We note a few simple facts pertaining to these objects.

First, \(|p^\#|=|p|\) on \(\mathbf {T}\). Moreover, every zero of p on \(\mathbf {T}\) is also a zero of \(p^\#\) and is thus a zero of \(S_np\). Finally,

$$\begin{aligned} (S_np)(0)=|a_0|^2-|a_n|^2\in {\mathbb {R}}. \end{aligned}$$

We now construct a chain of polynomials \(p_0,p_1,\ldots ,p_n\) starting with \(p_0=p\) and then taking successive Schur transforms,

$$\begin{aligned} p_1=S_np_0,\quad p_2=S_{n-1}p_1,\qquad \ldots ,\qquad p_n=S_1p_{n-1}. \end{aligned}$$

The last polynomial \(p_n\) is an element in \(\mathcal {P}_0\) and thus is a constant.

Lemma 3.2

If \(p_k\) has no zeros on the unit circle \(\mathbf {T}\) and \(p_{k+1}(0)>0\), then \(p_k\) and \(p_{k+1}\) have equally many zeros in \(\overline{\mathbf {D}}\).

Proof

Let \(p_k(z)=b_0+b_1z+\cdots +b_{n-k}z^{n-k}\). Then \(p_{k+1}(z)={\overline{b}}_0p_k-b_{n-k}p_k^\#\) and \(p_{k+1}(0)=|b_0|^2-|b_{n-k}|^2>0\) so \(|b_0|>|b_{n-k}|\). Since \(p_k(z)\ne 0\) on \(\mathbf {T}\) we have \(|{\overline{b}}_0p_k|>|b_{n-k}p_k|=|b_{n-k}p_k^\#|\) on \(\mathbf {T}\) and thus Rouché’s theorem implies that \(p_{k+1}\) and \(p_k\) have equally many zeros in \(\overline{\mathbf {D}}\). \(\square \)

We are now ready to formulate Schur–Cohn’s test (e.g. [28]).

Theorem 3.1

A polynomial \(p\in \mathcal {P}_n\) has no zeros in \(\overline{\mathbf {D}}\) if and only if \(p_k(0)>0\) for \(k=1,\ldots ,n.\)

Proof

Assume that \(p=p_0\) has no zeros in \(\overline{\mathbf {D}}\). By Lemma 3.1, \(p_1(0)=|a_0|^2-|a_n|^2>0\). Since \(p_0\) has no zeros on \(\mathbf {T}\), Lemma 3.2 implies that \(p_1\) and \(p_0\) have the same number of zeros in \(\overline{\mathbf {D}}\), i.e., none. If \(n\ge 2\) we can repeat the reasoning with \(p_1\in \mathcal {P}_{n-1}\) instead of \(p_0\) and deduce \(p_2(0)>0\) and so on, and after a finite number of steps we finally get \(p_n(0)>0\).

To prove the reverse implication, assume that \(p_k(0)>0\) for \(k=1,\ldots ,n\). Then \(p_n\) is a non-zero constant, and especially has no zeros on \(\mathbf {T}\). Since \(p_n=S_1p_{n-1}\), \(p_{n-1}\) has no zeros on \(\mathbf {T}\) and by Lemma 3.2, \(p_{n-1}\) and \(p_n\) has the same number of zeros in \(|z|\le 1\), i.e., none. We may now repeat this with \(p_{n-1}\in \mathcal {P}_1\) instead of \(p_n\) and deduce that \(p_{n-2}\) has no zeros in \(\overline{\mathbf {D}}\) and so on, and after a finite number of steps we finally get that \(p=p_0\) has no zeros in \(\overline{\mathbf {D}}\). \(\square \)

4 Real Comprehensive Triangular Decomposition

The method of real triangular decompositions, introduced recently in [7, §4], [5, §10] and [6], provides a suitable framework to deal with the the mapping problem for quadrature domains, in the form of systems such as (4), (5). These methods in turn make use of the idea of a Cylindrical Algebraic Decomposition [3, §5], for an overview of this and other applicable methods from real algebraic geometry see, for example, the book [3].

The objects we consider here are semi-algebraic sets, given a list of polynomial equations and inequalities in \(\mathbf {R}[x_1,\dots , x_m]\) a basic semi-algebraic set is the set of all points in \(\mathbf {R}^m\) which simultaneously satisfy all these equations and inequalities, [3, §3]. In this section we will specifically consider semi-algebraic systems defined by polynomials in the polynomial ring

$$\begin{aligned} \mathbb {Q}[c,x]=\mathbb {Q}[c_1,\dots , c_d][x_1,\dots , x_n], \end{aligned}$$

where we think of the \(c_i\) as parameters and the \(x_j\) as variables.

More precisely, given polynomials \(f_1,\dots , f_r\) and \(p_1,\dots , p_s\) in \(\mathbb {Q}[c,x]\), we define the semi-algebraic system \([f_{=0},p_{>0}]\) to be the following set of equations and inequalities

$$\begin{aligned} f_1(c,x)=0,\quad \ldots ,\quad f_r(c,x)=0,\quad p_1(c,x)>0,\quad \ldots ,\quad p_s(c,x)>0. \end{aligned}$$
(12)

The set of real solutions \((c,x)\in \mathbf {R}^d\times \mathbf {R}^n\) to (12) is called the (parameterized) semi-algebraic set generated by the system, denoted by \(\mathcal {S}([f_{=0},p_{>0}])\). Moreover, for fixed \(c\in \mathbf {R}^d\) we define the specialized semi-algebraic set \( \mathcal {S}_{(c)}([f_{=0},p_{>0}])\) as the set of points \(x\in \mathbf {R}^n\) which satisfy the system (12) for the particular parameter-value c.

Now suppose that \(\mathfrak {T}=\left\{ \mathcal {T}_1,\dots , \mathcal {T}_\ell \right\} \) is a collection of semi-algebraic systems in \(\mathbb {Q}[c,x]\). We extend the definitions of semi-algebraic set and specialization to a parameter value c by

$$\begin{aligned} \mathcal {S}(\mathfrak {T})=\bigcup _{j=1}^\ell \mathcal {S}(\mathcal {T}_j)\subset \mathbf {R}^d\times \mathbf {R}^n, \qquad \mathcal {S}_{(c)}(\mathfrak {T})=\bigcup _{j=1}^\ell \mathcal {S}_{(c)}(\mathcal {T}_j)\subset \mathbf {R}^n. \end{aligned}$$

Moreover, given a semi-algebraic system \(\mathcal {T}=[f_{=0},p_{>0}]\) we define the constructible set of \(\mathcal {T}\) to be the set of complex solutions \((c,x)\in \mathbf {C}^d\times \mathbf {C}^n\) to the system of equations and inequalities

$$\begin{aligned} f_1(c,x)=\ldots =f_r(c,x)=0,\quad p_1(c,x)\ne 0,\quad \ldots ,\quad p_s(c,x)\ne 0. \end{aligned}$$
(13)

We denote by \(\mathbf {C}\mathcal {S}(\mathcal {T})\) the set of solutions \((c,x)\in \mathbf {C}^d\times \mathbf {C}^n\) to (13); given \(c\in \mathbf {C}^d\) we define the associated specialization \(\mathbf {C}\mathcal {S}_{(c)}(\mathcal {T})\) to be the set of points \(x\in \mathbf {C}^n\) such that \((c,x)\in \mathbf {C}\mathcal {S}(\mathcal {T})\).

If \(\mathfrak {T}=\left\{ \mathcal {T}_1,\dots , T_\ell \right\} \) is a finite collection of semi-algebraic systems, it should now be obvious how to extend our definitions of constructible set (\(\mathbf {C}\mathcal {S}(\mathfrak {T})=\bigcup _{j} \mathbf {C}\mathcal {S}(\mathcal {T}_j)\)) and specializations (\(\mathbf {C}\mathcal {S}_{(c)}(\mathfrak {T})=\bigcup _{j}\mathbf {C}\mathcal {S}_{(c)}(\mathcal {T}_j)\subset \mathbf {C}^n\)).

A semi-algebraic system \( \mathcal {T}=[f_{=0},p_{>0}]\) is called square-free if all polynomials \(f_j\) and \(p_i\) occurring in \(\mathcal {T}\) are square-free. (A polynomial \(g\in \mathbb {Q}[c,x]\) is square-free if it has no factor of the form \(w^2\) where \(w\in \mathbb {Q}[c,x]\) is non-constant.)

Definition 4.1

Let \(\mathcal {W}=[f_{=0},p_{>0}]\) be a semi-algebraic system defined by polynomials \(f_1,\dots , f_r\) and \(p_1,\dots , p_s\) in \(\mathbb {Q}[c,x]\) and let \(\mathcal {S}(\mathcal {W})\subset \mathbf {R}^d\times \mathbf {R}^n\) be the associated semi-algebraic set.

A real comprehensive triangular decomposition (RCTD) of \({\mathcal {W}}\) is a pair \((\mathfrak {C},\{\mathfrak {T}_C \;;\; C\in \mathfrak {C}\} )\) where \(\mathfrak {C}\) is a finite partition of \(\mathbf {R}^d\) into non-empty semi-algebraic sets C (called “cells”) and for each \(C\in \mathfrak {C}\), \(\mathfrak {T}_C\) is a finite set of square-free semi-algebraic systems such that exactly one of the following holds:

  1. (1)

    \(\mathfrak {T}_C\) is empty so \(\mathcal {S}(\mathfrak {T}_C)=\mathbf {R}^d\times \mathbf {R}^n\) and \(\mathcal {S}_{(c)}(\mathfrak {T}_C)=\mathbf {R}^n\) for all \(c\in \mathbf {R}^d\),

  2. (2)

    The specialized constructible set \(\mathbf {C}\mathcal {S}_{(c)}(\mathfrak {T}_C)\) is infinite for all \(c\in C\),

  3. (3)

    \(\mathfrak {T}_C=\left\{ \mathcal {T}_1,\dots , \mathcal {T}_\ell \right\} \) is a finite set of semi-algebraic systems satisfying the following conditions:

    • \(\mathbf {C}\mathcal {S}_{(c)}(\mathfrak {T}_C)\) is finite and has fixed cardinality for all \(c\in C\),

    • the specialized semi-algebraic sets \(\mathcal {S}_{(c)}(\mathcal {T}_j)\) are finite and non-empty for all j and further for a fixed \(\mathcal {T}_j\) the specialized semi-algebraic set \(\mathcal {S}_{(c)}(\mathcal {T}_j)\) has fixed cardinality for all \(c\in C\),

    • \(\mathcal {S}_{(c)}(\mathcal {W})=\bigsqcup _{j=1}^\ell \mathcal {S}_{(c)}(\mathcal {T}_j) \) for all \(c\in C\).

The following proposition summarizes the results about RCTD’s that we will apply in the sequel.

Proposition 4.1

(§10 of [5]) Let \(f_1,\dots , f_r\) and \(p_1,\dots , p_s\) be polynomials in the polynomial ring \(\mathbb {Q}[c,x]\) defining a semi-algebraic system \(\mathcal {W}=[f_{=0},p_{>0}]\). Then a real comprehensive triangular decomposition of \(\mathcal {W}\) exists and may be computed by an explicit algorithm which is guaranteed to terminate in finite time. This algorithm is implemented in the RegularChains Maple package [31].

Example 4.1

Define \(f=a\cdot y^2+b\cdot x +3\) and \(g=b\cdot xy-y+2\) where xy denote real variables and ab denote real parameters. Consider the semi-algebraic system

$$\begin{aligned} {\mathcal {W}}=[f=0,g=0,y\le 1]. \end{aligned}$$

Using the RegularChains package, the parameter space \(\mathbf {R}^2\) is partitioned into five cells \(\mathfrak {C}=\{C_0,C_1,C_2,C_3, C_\infty \}\) where

$$\begin{aligned}&C_\infty =\left\{ \left( -\frac{3}{4},0\right) \right\} ,\;\;\;C_0=\left\{ (a,b)\in \mathbf {R}^2 \;;\; a\ne -\frac{3}{4},\;\; b=0 \right\} , \\&C_1=\left\{ (a,b)\in \mathbf {R}^2 \;;\; a<-\frac{64}{27}, \;b\ne 0 \right\} \bigsqcup \left\{ (a,b)\in \mathbf {R}^2 \;|\; a\ge 0,\; b\ne 0 \right\} , \\&C_2=\left\{ (a,b)\in \mathbf {R}^2 \;;\; a=-\frac{64}{27}, \;b\ne 0\right\} \bigsqcup \left\{ (a,b)\in \mathbf {R}^2 \;|\; -2<a<0, \; b\ne 0 \right\} ,\\&C_3=\left\{ (a,b)\in \mathbf {R}^2 \;;\; -\frac{64}{27}<a\le -2, \; b\ne 0 \right\} . \end{aligned}$$

These cells are illustrated in Fig. 4.

Fig. 4
figure 4

The cells \(\mathfrak {C}\) from the RCTD in Example 4.1

It is not hard to show that the specialized semi-algebraic set \(\mathcal {S}_{(a,b)}(\mathcal {W})\) consists of j points for all \((a,b)\in C_j\) where \(j\in \{0,1,2,3\}\); for the cell \(C_\infty \) the RCTD guarantees only that the specialized constructible set associated to the parameter choice \(C_\infty \) has infinitely many points (a corresponding specialized semi-algebraic set may be infinite, empty, or finite). Hence \(\mathfrak {C}\), along with the associated semi-algebraic systems for each parameter cell, gives a RCTD of \(\mathcal {W}\).

In more detail; consider the cells \(C_0\) and \(C_\infty \). The disjoint union of these cells is the line \(b=0\) (i.e. the a-axis in Fig. 4). Along this line the semi-algebraic system \(\mathcal {W}\) simplifies to

$$\begin{aligned} a\cdot y^2+3=2-y=0, \;\; y\le 1. \end{aligned}$$
(14)

Solving for y we obtain \(y=\sqrt{-3/a}\) and \(y=2\), so \(f=g=0\) has a real solution if and only if \(a=-3/4\), so \(\{(a,b)\}=C_\infty \). The corresponding specialized constructible set is

$$\begin{aligned} \left\{ (x,y)\in \mathbf {C}^2\;;\;(y-2)(y+2)=(y+2)=0, \; \; y\ne 1\right\} =\{(x,2)\;;\; x\in \mathbf {C}\}, \end{aligned}$$

which is infinite. However since \(y\le 1\) is never satisfied the semi-algebraic set \({\mathcal {S}}({\mathcal {I}}_{C_\infty })\) is empty. On the other hand in the cell \(C_0\) it is clear that the semi-algebraic system (14) has no real solutions, since \(a\ne \frac{-3}{4}\) when \((a,b)\in C_0\).

Now consider the subset P of the cell \(C_2\) consisting of all (ab) where \(a=-64/27\) and \(b\ne 0\). For \((a,b)\in P\) the system \(\mathcal {W}\) simplifies to

$$\begin{aligned} \left[ x=\frac{7}{3b}, \;\; y=\frac{-3}{2},\;\; y\le 1\right] \;\; \mathrm{or}\;\;\left[ x=\frac{-5}{3b}, \;\; y=\frac{3}{4},\;\; y\le 1\right] . \end{aligned}$$
(15)

Clearly the system (15) has precisely two solutions for any \(b\ne 0\), namely \((x,y)=\left( 7/(3b), -3/2 \right) \) and \((x,y)=\left( -5/(3b), 3/4\right) \). Similarly all other choices of \((a,b)\in C_2\) yield associated specialized semi-algebraic sets with exactly two points.

5 Some Computational Results

In this section we apply the method of real comprehensive triangular decompositions to obtain (new) proofs of uniqueness for certain families of quadrature domains.

5.1 Aharonov–Shapiro 1976

In this subsection, we shall give an alternative proof of a theorem of Aharonov and Shapiro from [2], which states that a solid quadrature domain obeying a quadrature identity of the form

$$\begin{aligned} \int _\Omega f\ dA=M_0f(0)+M_1f'(0),\qquad f\in AL^1(\Omega ) \end{aligned}$$
(16)

is unique. Here \(M_0\) is the area of \(\Omega \), so necessarily \(M_0>0\). The constant \(M_1\) is allowed to be an arbitrary complex number.

Recall that the computations in Example 2.2 show that a normalized mapping function \({\varphi }:\mathbf {D}\rightarrow \Omega \) is necessarily a polynomial of degree at most 2, which solves the system

$$\begin{aligned} {\left\{ \begin{array}{ll} w_1^3{\bar{w}}_1=M_0w_1^2-M_1w_2\\ w_1^2{\bar{w}}_2=2M_1 \end{array}\right. },\qquad w_1:={\varphi }'(0)>0,\quad w_2:={\varphi }''(0). \end{aligned}$$
(17)

It remains to show that this system gives rise to a unique univalent solution. For this, we write

$$\begin{aligned} M_1:=m_1+in_1,\qquad w_2:=u_2+iv_2. \end{aligned}$$

We now obtain the following semi-algebraic system, which is equivalent to (17),

$$\begin{aligned} (\mathcal {W})\qquad {\left\{ \begin{array}{ll} w_1^4=M_0w_1^2-(m_1u_2-n_1v_2)\\ 0=m_1v_2+n_1u_2\\ w_1^2u_2=2m_1\\ -w_1^2v_2=2n_1\\ M_0>0\\ w_1>0\\ (w_1^2M_0-(m_1u_2-n_1v_2))^2-4w_1^2(m_1^2+n_1^2)\ge 0 \end{array}\right. }. \end{aligned}$$

We must verify that the system \(\mathcal {W}\) gives rise to at most one univalent solution for all relevant choices of quadrature data. For this, we first recognize that the last condition in \(\mathcal {W}\) is just the Schur–Cohn constraint (see Theorem 3.1), which ensures that \({\varphi }'(z)\ne 0\) in \(\mathbf {D}\). In general this is only necessary for univalence but for polynomials of degree two it is also sufficient. We shall treat \(c=(m_1,n_1,M_0)\) as parameters. Computing a real comprehensive triangular decomposition of the system \(\mathcal {W}\), using the RegularChains library, we obtain a partition of the parameter space \(\mathbf {R}^3\) into two cells \(C_0,C_1\) having the following properties.

All points in the cell \(C_0\) are such that \(M_0\le 0\); hence no point of \(C_0\) can correspond to a quadrature domain. With \(\omega :=\root 3 \of {(27m_1^2+27n_1^2)/2}\), the cell \(C_1\) is expressed as the disjoint union of the following six subsets of \(\mathbf {R}^3\),

  1. (i)
    $$\begin{aligned} {\left\{ \begin{array}{ll} M_0=\omega \\ m_1=m_1\\ n_1<0 \end{array}\right. },\ \ {\left\{ \begin{array}{ll} M_0=\omega \\ m_1=m_1\\ n_1>0 \end{array}\right. } \end{aligned}$$
  2. (ii)
    $$\begin{aligned} {\left\{ \begin{array}{ll} M_0=\omega \\ m_1<0\\ n_1=0 \end{array}\right. } \end{aligned}$$
  3. (iii)
    $$\begin{aligned} {\left\{ \begin{array}{ll} M_0=\omega \\ m_1>0\\ n_1=0 \end{array}\right. } \end{aligned}$$
  4. (iv)
    $$\begin{aligned} {\left\{ \begin{array}{ll} M_0>\omega \\ m_1=m_1\\ n_1<0 \end{array}\right. },\ \ {\left\{ \begin{array}{ll} M_0>\omega \\ m_1=m_1\\ n_1>0 \end{array}\right. } \end{aligned}$$
  5. (v)
    $$\begin{aligned} {\left\{ \begin{array}{ll} M_0>\omega \\ m_1<0\\ n_1=0 \end{array}\right. },\ \ {\left\{ \begin{array}{ll} M_0>\omega \\ m_1>0\\ n_1=0 \end{array}\right. } \end{aligned}$$
  6. (vi)
    $$\begin{aligned} {\left\{ \begin{array}{ll} M_0>0\\ m_1=0\\ n_1=0 \end{array}\right. }. \end{aligned}$$

For \((M_0,m_1,n_1)\) in each of these sets, we now prove that the system (17) has a unique solution \((w_1,w_2)=(w_1,u_2+iv_2)\) with \(w_1>0\). Indeed, straightforward calculations show that in each of the parameter-domains, (i)–(vi) the semi-algebraic system \(\mathcal {W}\) simplifies to, respectively,

  1. (i)
    $$\begin{aligned} {\left\{ \begin{array}{ll} v_2w_1^2+2n_1=0\\ n_1u_2+v_2m_1=0\\ (9m_1^2+9n_1^2)v_2+2M_0^2n_1=0\\ w_1>0 \end{array}\right. } \end{aligned}$$
  2. (ii)
    $$\begin{aligned} {\left\{ \begin{array}{ll} 2u_2w_1-u_2^2+2M_0=0\\ 9m_1u_2-2M_0^2=0\\ v_2=0 \end{array}\right. } \end{aligned}$$
  3. (iii)
    $$\begin{aligned} {\left\{ \begin{array}{ll} 2u_2w_1+u_2^2-2M_0=0\\ 9m_1u_2-2M_0^2=0\\ v_2=0 \end{array}\right. } \end{aligned}$$
  4. (iv)
    $$\begin{aligned} {\left\{ \begin{array}{ll} v_2w_1^2+2n_1=0\\ n_1u_2+v_2m_1=0\\ (m_1^2+n_1^2)v_3^2-2M_0n_1^2v_2-4n_1^3=0\\ w_1>0\\ -v_{{2}}^{2}M_{{0}}m_{{1}}^{2}n_{{1}}^{2}-v_{{2}}^{2}M_{{0}} n_{{1}}^{4}+2\,M_{{0}}^{2}n_1^{4}+6\,v_{{2}}m_{{1}}^{2}n_ {{1}}^{3}+6\,v_{{2}}n_{{1}}^{5}>0 \end{array}\right. } \end{aligned}$$
  5. (v)
    $$\begin{aligned} {\left\{ \begin{array}{ll} u_2w_1^2-2m_1=0\\ u_2^3-2M_0u_2+4m_1=0\\ v_2=0\\ w_1>0\\ -M_{{0}}m_{{1}}^{2}u_{{2}}^{2}+2\,M_{{0}}^{2}m_{{1}}^{2}-6\, m_{{1}}^{3}u_{{2}}-4\,m_{{1}}n_{{1}}^{2}u_{{2}}>0 \end{array}\right. } \end{aligned}$$
  6. (vi)
    $$\begin{aligned} {\left\{ \begin{array}{ll} w_1^2-M_0=0\\ u_2=0\\ v_2=0\\ w_1>0 \end{array}\right. } \end{aligned}$$

In each case (i)–(vi), we have a unique solution \((w_1,w_2)=(w_1,u_2+iv_2)\) where \(w_1>0\), which concludes our automated proof of uniqueness for solid q.d.’s obeying (16). \(\square \)

Remark

The uniqueness problem for quadrature domains of order 2 has been completely settled, see [19, Cor. 10.1]. More precisely, counting also non-simply connected domains and combining with results in [20], a given point functional \(\mu \) of order 2 can give rise to at most two quadrature domains: one simply connected \(\Omega \), and possibly another one of the form \(\Omega \setminus \{a\}\) where a is a “special point”, i.e., \(S(a)={\bar{a}}\) where S is the Schwarz function. The term “special point” is due to Shapiro in [39], see Theorem 2.9 and the discussion that precedes it.

These methods can also be applied to polynomials of higher degree. (This topic is being considered as part of an ongoing work by the two of the authors.)

5.2 Symmetric Smash Sums

Consider a domain \(\Omega \) that satisfies

$$\begin{aligned} \int _\Omega f\, dA=c(f(a_0)+\cdots +f(a_{n-1})),\qquad a_j=a e^{2\pi i j/n},\, f\in AL^1(\Omega ), \end{aligned}$$
(18)

where \(c>0\) and \(a<0\) are constants and n is a positive integer.

A domain satisfying (18) may be constructed as a potential theoretic sum (or “smash sum”) of discs \(D(a_j;\sqrt{c})\), where excess mass coming from overlapping discs is swept out using the process of partial balayage. When \(n=2\) this process gives rise to the well-known Neumann’s oval.

For \(n\ge 3\) one can surmise that a (connected) domain satisfying (18) should be either simply connected or doubly connected, depending on whether or not \(0\in \Omega \). The doubly connected case has already been treated in [43, Section 5.8] and in [9]. (See Subsect. 7.3 for further remarks in this connection.)

In this subsection, we shall focus on the simply connected case and show how uniqueness of the quadrature domain follows using our general algebraic scheme. (In this particular case, there are, of course, alternative ways to see this, e.g. by Novikov’s theorem, Theorem 1.1, (i).)

Due to the \({\mathbb {Z}}_n\) symmetry of (18), the conformal map \(\varphi :\mathbf {D}\rightarrow \Omega \) attains the simple form

$$\begin{aligned} \varphi (z)=\sum _{j=0}^{n-1}\frac{c_j}{{\bar{w}}_j}\frac{z}{1-{\bar{\lambda }}_jz}=\frac{c}{w}\sum _{j=0}^{n-1}\frac{z}{1-\exp \left( \frac{-2\pi i}{n}j\right) \lambda z}=\frac{cn}{w}\frac{z}{1-\lambda ^nz^n} \end{aligned}$$
(19)

where \(w_j=\varphi '(\lambda _j)=w>0\), \(-1<\lambda =\varphi ^{-1}(a)<0\). This gives

$$\begin{aligned} \varphi '(z)=\frac{cn}{w}\frac{1+(n-1)\lambda ^nz^n}{(1-\lambda ^nz^n)^2}. \end{aligned}$$
(20)

The two unknowns, \(\lambda \) and w (\(-1<\lambda <0\) and \(w>0\)) satisfy

$$\begin{aligned} {\left\{ \begin{array}{ll}\displaystyle a=\frac{cn}{w}\frac{\lambda }{1-\lambda ^{2n}}\\ \displaystyle w=\frac{cn}{w}\frac{1+(n-1)\lambda ^{2n}}{(1-\lambda ^{2n})^2} \end{array}\right. } \quad \iff \quad {\left\{ \begin{array}{ll} aw(1-\lambda ^{2n})=nc\lambda \\ w^2(\lambda ^{4n}-2\lambda ^{2n}+1)=nc((n-1)\lambda ^{2n}+1). \end{array}\right. } \end{aligned}$$

We only need to check solutions \({\varphi }\) which are locally univalent. From (20) (or from Schur–Cohn’s test) one may easily deduce that local univalence holds if \(1\ge (n-1)\lambda ^{n}\). Taking this into consideration we obtain the following semi-algebraic system

$$\begin{aligned} {\left\{ \begin{array}{ll} aw(1-\lambda ^{2n})=nc\lambda \\ w^2(\lambda ^{4n}-2\lambda ^{2n}+1)=nc((n-1)\lambda ^{2n}+1)\\ c>0\\ a<0\\ w>0\\ -1<\lambda <0\\ 1\ge (n-1)^2\lambda ^{2n}. \end{array}\right. } \end{aligned}$$
(21)

When setting up a real comprehensive triangular decomposition, we cannot treat n as a parameter but have to use a fixed value of n (since if n is treated as a parameter the system (21) is no longer semi-algebraic). Below we show the results for \(n=3\) but the corresponding decompositions may be obtained for up to \(n=10\) without any runtime issues.

Thus let \(\mathcal {W}\) be the semi-algebraic system (21) with \(n=3\) fixed. The associated semi-algebraic set \(\mathcal {S}(\mathcal {W})\) is contained in \(\mathbf {R}^2\times \mathbf {R}^2\). We treat \((a,c)\in \mathbf {R}^2\) as parameters and consider the specialized semi-algebraic set \(\mathcal {S}_{(a,c)}(\mathcal {W})\subset \mathbf {R}^2\). Using an RCTD, we have found that the parameter space \(\mathbf {R}^2\) is divided into two cells \(C_0,C_1\) such that for all \((a,c)\in C_1\) the system \(\mathcal {W}\) has precisely one real solution and for all \((a,c)\in C_0\) the system \(\mathcal {W}\) has no real solutions. The cell \(C_1\) is made up of the three regions defined in (i), (i), and (iii) below.

  1. (i)
    $$\begin{aligned} {\left\{ \begin{array}{ll}\displaystyle c=\frac{a^2}{2^{1/3}}\\ a<0 \end{array}\right. } \end{aligned}$$
  2. (ii)
    $$\begin{aligned} {\left\{ \begin{array}{ll}\displaystyle \frac{a^2}{2^{1/3}}<c<a^2\\ a<0 \end{array}\right. }, \ \ \ {\left\{ \begin{array}{ll} a^2<c\\ a<0 \end{array}\right. } \end{aligned}$$
  3. (iii)
    $$\begin{aligned} {\left\{ \begin{array}{ll} c=a^2\\ a<0 \end{array}\right. }. \end{aligned}$$

In each of the domains above, the semi-algebraic system (21) simplifies as follows

  1. (i)
    $$\begin{aligned} {\left\{ \begin{array}{ll} (\lambda ^{6}-1)aw+3\lambda c=0\\ a^2\lambda +c=0 \end{array}\right. } \end{aligned}$$
  2. (ii)
    $$\begin{aligned} {\left\{ \begin{array}{ll} (\lambda ^{6}-1)aw+3\lambda c=0\\ 2a^2\lambda ^6-3c\lambda ^2+a^2=0\\ w>0\\ -4\lambda ^6+1>0\\ -\lambda>0\\ \lambda +1>0 \end{array}\right. } \end{aligned}$$
  3. (iii)
    $$\begin{aligned} {\left\{ \begin{array}{ll} (\lambda ^2-1)w+2\lambda a=0\\ 2\lambda ^4+2\lambda ^2-1=0\\ w>0\\ -4\lambda ^6+1>0\\ -\lambda>0\\ \lambda +1>0 \end{array}\right. } \end{aligned}$$

Since for all choices of the parameters in the cell \(C_1\) we know from the definition of the RCTD that the number of real solutions is fixed, picking a particular a and c in \(C_1\) it is easy to verify that the systems above have unique real solutions which correspond to a quadrature domain. Similarly, one can pick a choice of parameters in \(C_0\) to verify that the system (21) has no real solutions for any parameter choice in \(C_0\). In particular, we have shown that a quadrature domain obeying a quadrature identity (18) is uniquely determined. However, we cannot yet be sure that each solution to the system (21) defines a univalent map since we have only required local univalence. The following proposition shows that this is in fact the case.

Proposition 5.1

The mapping \(\varphi \) in (19), (\(-1<\lambda <0\)), is univalent if and only if

$$\begin{aligned} 1\ge (n-1)^2\lambda ^{2n}. \end{aligned}$$
(22)

Proof

The necessity is obvious since the strict inequality \(1>(n-1)^2\lambda ^{2n}\) is just the Schur–Cohn condition for \(\varphi '\) to have no zeros in \(\overline{\mathbf {D}}\).

To prove sufficiency we shall prove that (22) implies that \(\varphi \) is starlike, and in particular univalent. Recall ( [34, Thm. 2.5]) that starlikeness is equivalent to that \( \mathrm {Re}\left( z\varphi '(z)/\varphi (z)\right) >0\) for all \(z\in \mathbf {D}\). Using the Cauchy-Riemann equations one finds that this is equivalent to

$$\begin{aligned} \frac{d}{dr}|\varphi (z)|^2>0,\qquad z=re^{i\theta }\in \mathbf {D}. \end{aligned}$$

(Cf. [16]). Using Eq. (19) we find that \(\frac{d}{dr}|\varphi (z)|^2\) is a positive multiple of

$$\begin{aligned} 2r\left( \lambda ^{2n}r^{2n}(1-n)-\lambda ^nr^n(2-n)\cos (n\theta )+1\right) . \end{aligned}$$
(23)

This is a quadratic polynomial in \(\lambda ^n\) which is clearly positive when \(\lambda =0\). We need to prove that it remains positive for all r with \(0\le r< 1\) if \(\lambda <0\) and (22) holds. We may assume that \(r=1\) and \(\cos (n\theta )=\pm 1\). For odd n, the non-negativity of (23) then means

$$\begin{aligned} \lambda ^{2n}(1-n)-\lambda ^n(2-n)+1\ge 0 \end{aligned}$$
(24)

while for even n, it means

$$\begin{aligned} \lambda ^{2n}(1-n)+\lambda ^n(2-n)+1\ge 0. \end{aligned}$$
(25)

By assumption \(\lambda <0\), so the inequalities (24), (25) amount to \(\lambda ^n(n-1)\le 1\) and \(-1\ge \lambda ^n(n-1)\) respectively. This is equivalent to \(1\ge \lambda ^{2n}(n-1)^2\) and proves that \(\varphi \) is starlike if and only if (22) holds. \(\square \)

6 Three Examples of Uniqueness

In this section we give several further examples of uniqueness of quadrature domains based on real comprehensive triangular decomposition and other computational methods from algebraic geometry (i.e. Gröbner basis computations).

Example 6.1

Let \(\Omega \) be a quadrature domain satisfying

$$\begin{aligned} \int _\Omega f\ dA=f(0)+f\left( \frac{265}{153}+i\right) +f\left( \frac{265}{153}-i\right) ,\qquad f\in AL^1(\Omega ) \end{aligned}$$
(26)

the nodes are chosen such that the union \(D(0,1)\cup D(265/153+i,1)\cup D(265/153-i,1)\) is just barely simply connected (this is guaranteed by \(265/153<\sqrt{3}\)), see Fig. 5.

To find the mapping \(\varphi :\mathbf {D}\rightarrow \Omega \) we need to solve the six complex equations \(a_i=\varphi (\lambda _i),\ w_i=\varphi '(\lambda _i),\ i=1,2,3\), or equivalently 12 real equations. Without loss of generality we pick \(\lambda _1=0=a_0\) and \(w_1>0\), moreover, using the mirror symmetry in (26) we get \(\lambda _2={\bar{\lambda }}_3\) and \(w_2={\bar{w}}_3\). This means that in total we only have to solve for five real variables. The equations are, of course, rational, meaning that in order to calculate the Gröbner basis we need to clear denominators and algebraically remove the new roots introduced by this. These new roots correspond to poles of the original rational system, and hence are not of interest. The defining equations without roots at the poles are obtained by saturating the resulting ideal by the product of the equations from the denominators.

A Gröbner basis calculation shows that there is a unique solid quadrature domain obeying (26), depicted in Fig. 5. To be explicit, we can write \((\lambda _1,\lambda _2,\lambda _3)=(0,l+ip,l-ip)\) and \((w_1,w_2,w_3)=(x_1,x_2+iy_2,x_2-iy_2)\), where

$$\begin{aligned} l= & {} 0.866022320861578...,\\ p= & {} 0.499994659802733...,\\ x_1= & {} 1.00001067993828...,\\ x_2= & {} 93633.9999439279...,\\ y_2= & {} 0.866068567766777.... \end{aligned}$$
Fig. 5
figure 5

The quadrature domain defined by Eq. (26). (A close inspection shows that the domain is in fact simply connected.)

Note that the values \(l,p,x_1,x_2,y_2\) above can be obtained as exact expressions in a field extension of \(\mathbb {Q}\) from the result of our Gröbner basis computation, however we opt to give the numerical approximations here for simplicity.

Example 6.2

Let \(c_1>0\) and let \(\Omega \) be a quadrature domain satisfying

$$\begin{aligned} \int _\Omega f\ dA=c_1\left[ f\left( -2\right) +f\left( -2e^{2\pi i/3}\right) +f\left( -2e^{4\pi i/3}\right) \right] +f(0),\qquad f\in AL^1(\Omega ). \end{aligned}$$

This is similar to a \({\mathbb {Z}}_3\)-symmetric domain in Subsect. 5.2, but with an additional node at the origin. Assuming \(\Omega \) is simply connected, the mapping \(\varphi :\mathbf {D}\rightarrow \Omega \) is given by

$$\begin{aligned} \varphi (z)=\frac{3c_1}{w_1}\frac{z}{1-\lambda ^3z^3}+\frac{c_2}{w_2}z \end{aligned}$$

where \(\varphi (\lambda )=-2=a,\ \varphi '(\lambda )=w_1,\ \varphi '(0)=w_2\) and \(c_2=1\). Clearing denominators leads to the semi-algebraic system

$$\begin{aligned} {\left\{ \begin{array}{ll} 3c_1\lambda (\lambda ^6-1-w_2)=-2w_1w_2(\lambda ^6-1)\\ 3c_1(\lambda ^12-2\lambda ^6(1-w_2)+1+w_2)=w_1^2w_2(\lambda ^6-1)\\ 3c_1(1+w_2)=w_1w_2^2\\ c_1>0\\ -1<\lambda <0\\ w_1>0,\;w_2>0 \end{array}\right. } \end{aligned}$$
(27)

with variables \(w_1,w_2,\lambda \) and parameter \(c_1\). Let \(\mathcal {W}\) denote the semi-algebraic system (27). Computing a triangular decomposition yields two cells: \(C_0=\{ c_1\in \mathbf {R}\;|\; c_1\le 1 \}\) and \(C_1=\{ c_1\in \mathbf {R}\;|\; c_1> 1 \}\) in the parameter space \(\mathbf {R}\) (with coordinate \(c_1\)). For each \(c_1\in C_0\) the system \(\mathcal {W}\) has no real solutions while for each \(c_1\in C_1\) the system \(\mathcal {W}\) has exactly one real solution \((w_1,w_2,\lambda )\in \mathbf {R}^3\); note that \(\mathbf {R}=C_0 \sqcup C_1\).

Observe that to obtain a solution, we must have \(c_1>1\) in order that \(\Omega \) be simply connected (see Fig. 6(c)). If we interpret \(c_1\) as time, we may think of the quadrature domain \(\Omega (c_1)\) as the result of a Hele–Shaw evolution, starting from \(\Omega (0)=D(0,1)\), after injecting fluid through the nodes \(-2,-2e^{2\pi i/3},-2e^{4\pi i/3}\) at constant rate 1, for \(c_1\) units of time (see Subsect. 7.2 below for more about this). The evolution is depicted in Fig. 6.

Fig. 6
figure 6

Hele–Shaw growth of the unit disc \(\mathbf {D}\) under injection at constant rate 1 at three nodes \(-2,-2e^{2\pi i/3},-2e^{4\pi i/3}\). (a) and (b) lie in the cell \(C_0\) of the RCTD; (a) is the initial domain \(\mathbf {D}\) (time \(t=0\)) with the three injection points plotted, (b) corresponds to \(t=1/2\) and is disconnected. In (c), \(t=c_1>1\) so the domain is simply connected and lies in cell \(C_1\) as do (d), (e), and (f)

Example 6.3

Let \(\Omega \) be a quadrature domain satisfying

$$\begin{aligned} \int _\Omega f\ dA=-\frac{1}{2}(f(a)+f(-a))+2f(0),\qquad f\in AL^1(\Omega ). \end{aligned}$$
(28)

This is a \({\mathbb {Z}}_n\)- symmetric domain with \(n=2\) as in Subsect. 5.2 but with negative weights and an additional node at the origin with strength two. The mapping \(\varphi :\mathbf {D}\rightarrow \Omega \) then satisfies

$$\begin{aligned} \varphi (z)=\frac{2c_1}{w_1}\frac{z}{1-\lambda ^2z^2}+\frac{c_2}{w_2}z \end{aligned}$$

where \(\varphi (\lambda )=a,w_1=\varphi '(\lambda ),w_2=\varphi '(0),c_1=-1/2\) and \(c_2=2\). To obtain unique solutions for this system we have to impose Schur–Cohn constraints; this is as easily done for arbitrary n as for \(n=2\) so we do it for arbitrary n.

Set \(p(z)=b_0+b_1z^{n}+b_2z^{2n}\) where \(b_0=c_1w_2n,\ b_1=\lambda ^n(c_1w_2n(n-1)-2c_2w_2)\) and \(b_2=c_2w_1\lambda ^{2n}\). The derivative of the mapping \(\varphi \) is given by

$$\begin{aligned} \varphi '(z)=\frac{c_1n}{w_1}\frac{1+(n-1)\lambda ^nz^n}{(1-\lambda ^nz^n)^2}+\frac{c_2}{w_2}=\frac{p(z)}{w_1w_2(1-\lambda ^nz^n)^2}. \end{aligned}$$

When \(n=2\) we have that \(b_0=4c_1w_2,\ b_1=\lambda ^2(2c_1w_2-2c_2w_2)\) and \(b_2=c_2w_1\lambda ^{4}\). It follows that the Schur–Cohn constraints for the \(n=2\) case are given by

$$\begin{aligned} {\left\{ \begin{array}{ll} b_0^2-b_2^2\ge 0\\ (b_0^2-b_2^2)^2-(b_0b_1-b_1b_2)^2\ge 0, \end{array}\right. } \end{aligned}$$

where we allow for equality to include the case where \(\partial \Omega \) has cusps, see Fig. 7. Clearing denominators in the equations \(\varphi (\lambda )=a,w_1=\varphi '(\lambda ),w_2=\varphi '(0)\), recalling that \(c_1=-1/2\), \(c_2=2\), and adding the Schur–Cohn constraints we obtain the semi-algebraic system

$$\begin{aligned} {\left\{ \begin{array}{ll} \lambda \, \left( 2\,{\lambda }^{4}{w_1}-2\,{ w_1}+{ w_2} \right) =aw_{{1}}w_{{2}} \left( {\lambda }^{4}-1 \right) \\ w_1^2w_2\left( {\lambda }^{4}-1 \right) ^{2}=\varphi '(\lambda )= \left( 2\,{\lambda }^{8}-4\,{\lambda }^{4}+2 \right) w_{{1}}-{\lambda }^{4}w_{{2}}-w_{{2}}\\ w_1w_2^2=2w_1-w_2\\ b_0^2-b_2^2\ge 0\\ (b_0^2-b_2^2)^2-(b_0b_1-b_1b_2)^2\ge 0\\ a<0\\ -1<\lambda <0\\ w_1,w_2>0 \end{array}\right. } \end{aligned}$$
(29)

where we treat a as a parameter and \((w_1,w_2,\lambda )\) as variables (note \(c_1,c_2\) are fixed). Let \(\mathcal {W}\) denote the semi-algebraic system (29). Let \(\tau =27\sqrt{3}+54\) and set

$$\begin{aligned} \vartheta = -\frac{1}{2} \, \sqrt{\frac{{\tau }^{2/3} + 4 \, {\tau }^{1/3} + 9}{{\tau }^{1/3}}} + \frac{1}{2} \, \sqrt{-{\tau }^{\frac{1}{3}} - \frac{9}{{\tau }^{1/3}} + \frac{16}{\sqrt{\frac{{\tau }^{2/3} + 4 \, {\tau }^{1/3} + 9}{{\tau }^{1/3}}}} + 8}\cong -0.369. \end{aligned}$$

A computation of a real comprehensive triangular decomposition for \(\mathcal {W}\) gives that the parameter space \(\mathbf {R}\) (with coordinate a) is partitioned into two cells: \(C_0=\{a\in \mathbf {R}\;|\; a<\vartheta \} \sqcup \{a\in \mathbf {R}\;|\; a\ge 0 \}\) and \(C_1=\{a\in \mathbf {R}\;|\; \vartheta \le a <0 \}\). For all \(a\in C_0\) the system \(\mathcal {W}\) has no solutions and for all \(a\in C_1\) the system \(\mathcal {W}\) has exactly one real solution \((w_1,w_2,\lambda )\in \mathbf {R}^3\). It follows that for any valid choice of the parameter a there exists a unique quadrature domain. The parameter value \(a=\vartheta \in C_1\) corresponds to the case where \(\partial \Omega \) has two cusps, as shown in Fig. 7.

Fig. 7
figure 7

The unique solid quadrature domain satisfying (28) for the parameter-value \(a=\vartheta \in C_1\)

7 Other Related Topics

In this section, we briefly discuss some of the existing methods and constructions in the theory of quadrature domains with bearing for our above methods. The reader interested in a comprehensive picture of the area should consult one or several of the textbooks [12, 15, 25, 26, 40, 43].

7.1 The Defining Polynomial of a q.d.

Consider a quadrature domain \(\Omega \) of order \(n=n_1+\cdots +n_m\) corresponding to a point-functional

$$\begin{aligned} \mu (f)=\sum _{k=1}^m\sum _{j=0}^{n_k-1} c_{kj}f^{(j)}(a_k). \end{aligned}$$
(30)

It was shown by Aharonov and Shapiro in [2] that the boundary \(\partial \Omega \) is an algebraic curve, and more precisely there exists an irreducible polynomial \(P(z,w)=\sum _{j,k=0}^n a_{jk}z^jw^k\) which is self-conjugate and normalized (i.e. \(a_{jk}={\bar{a}}_{kj}\), \(a_{nn}=1\)) such that

$$\begin{aligned} \partial \Omega =\{z\,;\, P(z,{\bar{z}})=0\}\setminus \{\text {``special points''}\}, \end{aligned}$$
(31)

where a “special point” is an isolated solution to the equation \(P(z,{\bar{z}})=0\). (Cf. [39])

It is convenient to refer to the polynomial P in (31) as the “defining polynomial” for the quadrature domain \(\Omega \). A natural problem, then, is to try to characterize the set of defining polynomials P which are associated with a given point-functional \(\mu \). This problem has been investigated by Gustafsson in [19, 20], and leads to interesting insights with respect to the uniqueness question (Q). To elucidate this, we write a defining polynomial in the form

$$\begin{aligned} P(z,w)=\sum _{j=0}^n w^jp_j(z), \end{aligned}$$
(32)

where each \(p_j(z)\) is a polynomial of degree at most n.

It is shown in [19] (cf.  also [8]) that there is an explicit bijection between point-functionals \(\mu \) (of the form (30)) and the last two polynomials \(p_{n-1},p_n\) which may appear in the expansion (32). As pointed out in [8, 19], the determination of the remaining polynomials \(p_0,\ldots ,p_{n-2}\) is generally a difficult matter. In the simply connected case, this problem is similar to the question of completely characterizing all solutions to the master formula (in Theorem 1.2), for each point-functional \(\mu \).

Remark

The notion of a special point is not just an artefact of the construction. For example, such points appear naturally in the process of forming smash sums of discs (i.e., each time discs “collide” as in Fig. 6 (c)) and they have the physical interpretation of stagnation points for fluid flows. We refer to [8, 20] for further details.

7.2 Connection to Laplacian Growth

Let t be a real parameter and consider a family of quadrature domains \(\Omega (t)\), each obeying the quadrature identity

$$\begin{aligned} \int _{\Omega (t)} f\, dA{=}\mu _t(f){=}(c_1+qt)f(a_1)+c_2f(a_2)+\cdots +c_nf(a_n),\quad f\in AL^1(\Omega ),\nonumber \\ \end{aligned}$$
(33)

where q is a real constant, and where we take \(a_1=0\) for simplicity. If \(q>0\), we may think of \(\Omega (t)\) as an expanding blob of fluid, obtained from \(\Omega (0)\) by injecting fluid at the origin, at constant rate q. (If \(q<0\), the domains contract due to suction.) The resulting evolution \(t\mapsto \Omega (t)\) is known under the names “Hele–Shaw evolution” and “Laplacian growth”, see e.g. [1, 26, 43].

For values of t such that \(\Omega (t)\) is solid, we will write \(z\mapsto {\varphi }(z,t)\) for the Riemann mapping \(\mathbf {D}\rightarrow \Omega (t)\).

Remark

Any domain \(\Omega \in Q(\mu _t,AL^1)\) is by definition of finite order, and hence its boundary is part of an algebraic curve, see Subsect. 7.1. This implies that the boundary curve is smooth everywhere with the possible exception of finitely many singular points, which may be either cusps \(p\in \partial \Omega \) pointing inwards (corresponding to values \(p={\varphi }(z)\), \(z\in \mathbf {T}\) at which \({\varphi }'(z)=0\)), or contact points \(p\in \partial \Omega \) (which satisfy \(p={\varphi }(z_1)={\varphi }(z_2)\) for two distinct points \(z_1,z_2\in \mathbf {T}\)). For solid domains, contact points are excluded.

In the following, we suppose that \(\Omega (t)\in Q(\mu _t,AL^1)\) is a smoothly varying family of solid domains obeying (33) for t in some suitable time-interval. A basic result relates the “time-derivative” \({\dot{{\varphi }}}(z,t):=\frac{\partial }{\partial t} {\varphi }(z,t)\) to the “space-derivative” \({\varphi }'(z,t)=\frac{\partial }{\partial z}{\varphi }(z,t)\).

Theorem 7.1

The conformal map \({\varphi }(z,t)\) onto \(\Omega (t)\) satisfies Polubarinova-Galin’s equation

$$\begin{aligned} {\text {Re}}\left[ {\dot{{\varphi }}}(z,t)\overline{z{\varphi }'(z,t)}\right] =\frac{q}{2},\qquad z\in \mathbf {T}. \end{aligned}$$
(34)

A proof can be found in the book [25, Sect. 1.4.2].

Gustafsson and Lin in [22] have studied the evolution of zeros and poles of the space-derivative \({\varphi }'(z,t)\). We will now briefly indicate a different possible approach based on our basic structure theorem, Theorem 1.2. Denote by \(\lambda _i(t)\) the points in \(\mathbf {D}\) such that \({\varphi }(\lambda _i(t),t)=a_i\), and write \(w_i(t)={\varphi }'(\lambda _i(t),t)\). Also denote \(c_1(t)=c_1+qt\) and \(c_j(t)\equiv c_j\) when \(j\ge 2\). In view of (3), the mapping \({\varphi }(z,t)\) obeys

$$\begin{aligned} {\varphi }(z,t)=\sum _{i=1}^n \frac{\overline{c_i(t)}}{\overline{w_i(t)}}\frac{z}{1-\overline{\lambda _i(t)} z}, \end{aligned}$$
(35)

which gives

$$\begin{aligned} {\varphi }'(z,t)=\sum _{j=1}^n\frac{\overline{c_i(t)}}{\overline{w_i(t)}}\frac{1}{(1-\overline{\lambda _i(t)} z)^2} \end{aligned}$$
(36)

and (denoting complex conjugation by \(a^\dagger ={\bar{a}}\), and abbreviating \(c_i=c_i(t)\), \(\lambda _i=\lambda _i(t)\), etc.)

$$\begin{aligned} {\dot{{\varphi }}}(z,t)=z\left[ \frac{1}{w_1}-\sum _{i=1}^n\left( \frac{c_i{\dot{w}}_i}{w_i^2}\right) ^\dagger \frac{1}{1-z{\bar{\lambda }}_i}+\sum _{i=1}^n\frac{{\bar{c}}_i}{{\bar{w}}_i}\frac{z\overline{\dot{\lambda _i}}}{(1-z{\bar{\lambda }}_i)^2}\right] . \end{aligned}$$
(37)

Substituting (36) and (37) in Polubarinova–Galin’s equation (34), we obtain a non-linear system of ordinary differential equations connecting the functions \(\lambda _i(t)\) and \(w_i(t)\) and their time-derivatives \({\dot{\lambda }}_i(t),{\dot{w}}_i(t)\). We may note from (35) that the poles of the Riemann map \({\varphi }(z,t)\) are given by \(1/\overline{\lambda _i(t)}\).

Similar sets of equations have been studied also in the case when \({\varphi }(z,t)=a_1(t)z+\cdots +a_n(t)z^n\) is a polynomial in z, see [26, Sect. 2.1.1] and [32]. It is clear that our method applies in this case as well, but for reasons of length, we shall not pursue this issue here.

7.3 Multi-connected Quadrature Domains

A theory for multiply connected quadrature domains, using the Schottky double of the domain, is developed in the paper [19]. The mapping problem for such domains has been the subject of numerous investigations. In particular, for a significant class of quadrature domains (e.g. based on forming suitable smash sums of discs), the Riemann map can be constructed using the corresponding “Schottky–Klein prime function”. This approach is found in the works [10, 13] as well as in the recent monograph [12, Ch. 11]. (The Schottky–Klein function is surveyed in the article [11].)

Other relevant works in this connection are [30] on the topology of quadrature domains, and [18], which discusses related uniqueness questions.