From Pythagoras to Fourier and From Geometry to Nature


Chapter 12. Solution of Problems in Gielis Domains

Full-Text Views:

Many applications of Mathematical Physics and Engineering are connected with the Laplacian:

  • The wave equation: utt = a2 Δ2u

  • Heat propagation: ut = κ Δ2u

  • The Laplace equation: Δ2u = 0

  • The Helmholtz equation: Δ2u + k2u = 0

  • The Poisson equation: Δ2u = f

  • The Schrödinger equation: h22mΔ2ψ+Vψ=Eψ

Boundary value problems relevant to the Laplacian are solved in explicit form only for domains with a very special shape, namely intervals, cylinders or domains with special (circular or spherical) symmetries [1]. In what follows, we limit ourselves to consider the extensions of classical problems to 2D normal polar domains of the Gielis type, that is domains 𝒟 which are starlike with respect to the polar coordinate system. Then 𝒟 can be interpreted as an anisotropically stretched unit circle. Other general problems, or relative to more complex shapes, have also been considered in [11, 12, 14, 17, 18]. Further extensions have been made to the case of 3D domains, but the relevant equations are much more involved. A list of such articles can be found in the References section (see [9, 10, 13, 15, 1921, 45]).

12.1 The Laplacian in Stretched Polar Coordinates

We introduce in the x, y plane the polar coordinates:

and the polar equation of 𝒟:
where r(θ) ∈ C2[0, 2π]. We suppose the domain 𝒟 satisfies:
and therefore minθ[0,2π] r(θ) > 0.

We introduce the stretched radius ρ* such that:

and the curvilinear (i.e. stretched) coordinates ρ*, θ in the plane x, y:

Therefore, 𝒟 is obtained assuming 0 ≤ θ ≤ 2π and 0 ≤ ρ* ≤ 1.

We show how to modify some classical formulas and we derive methods to compute the coefficients of Fourier-type expansions representing solutions of some classical problems. Of course, this theory can be easily generalized by considering weakened hypotheses on the boundary or initial data.

The case of the unit circle is recovered assuming ρ* = ρ and r(θ) ≡ 1. We consider a C2(𝒟) function u(x, y) = u(ρ cos θ, ρ sin θ) = U (ρ, θ) and the Laplace operator in polar coordinates:


We start representing this operator in the new stretched coordinate system ρ*, θ. Putting:

the unit circle is recovered by putting R(θ) ≡ 1.

Using this polar equation, the corresponding stretched coordinates ρ*, θ in the plane x, y are given by:

and assuming:
the Laplacian becomes:

For ρ* = ρ and R(θ) ≡ 1 we find the Laplacian in polar coordinates.

12.2 The Dirichlet Problem for the Laplace Equation in Gielis Domains

Consider the Dirichlet problem for the Laplace equation:


In [74] we have proven the result:

Theorem 12.1.



the solution of the internal Dirichlet problem can be represented as:

where a0 = α0/2 and the coefficients am, bm (m = 1, 2, 3, …) are given by solving the infinite system:


As an example, we start from the general Gielis equation [40]:

by choosing particular values of the parameters.

By assuming in Equation (12.12) that c = 22, α = 5, β = 8, m = 10, n1 = n3 = 6 and n2 = 4 we obtain the shape of the relevant domain 𝒟 in Figure 38.

Figure 38

Starfish domain.

Let f(x, y) = cosh(x +y) + 5x2y be the function representing boundary values. Then we obtain the results reported in Table 3. In the first column we show the L2(𝒟) norm of the boundary error fuh (where uh denotes the (2h + 1)th partial sum of the approximating Fourier series) and in the second column the L2(𝒟) norm of the inside error, i.e. the L2(𝒟) norm distance of Δuh from zero.

fu1L2 = 0.000335952 ∥Δu1L2 = 0. × 10−17
fu2L2 = 0.000133587 ∥Δu2L2 = 0. × 10−17
fu3L2 = 0.000101291
fu4L2 = 9.02500 × 10−5
fu5L2 = 5.42434 × 10−5
fu6L2 = 4.75581 × 10−5
fu7L2 = 4.75567 × 10−5
fu8L2 = 4.75565 × 10−5
Table 3

L2 norms of boundary and inside approximation errors.

The obtained results, with P. Natalini as a coauthor (see [74]), show the convergence (in general a.e.) of the approximating sequence of functions to the function f, according to the general results on Fourier series proven by L. Carleson [22].

12.3 The Heat Problem in Gielis Domains

The heat problem for a plate with a general shape is often reduced to the circular case by using the conformal mappings technique (see e.g. [35, 65]), but only very special cases can be treated analytically by using this method since only few explicit equations for the relevant conformal mappings are known. However, it is possible to use the stretched coordinates system in order to obtain a quite general result for a Gielis domain.

Consider a plate with normal polar shape 𝒟 and known diffusivity κ. Suppose the boundary temperature is zero for every t ≥ 0 and the initial temperature is given by the continuous function f(x, y) so that the problem of finding the temperature of the plate for every t > 0 is expressed by:


In [73], with P. Natalini and R. Patrizi as coauthors, the following result was proven:

Theorem 12.2.

The above heat problem admits a classical solution:

such that the following generalized Fourier expansion in terms of Bessel functions holds:

Putting U (ρ, θ, 0) = F (ρ, θ) ≕ G(ρ*, θ) where:

so that:
the coefficients Am,k, Bm,k are given by:

Remark 3.

Note that the above formulas still hold if the function r(θ) is a piecewise continuous function and if the initial data are given by square integrable functions, not necessarily continuous, so that the relevant coefficients αh, βh in Equation (12.15) are finite.


In the following example we consider, for the starlike plate, a Gielis equation of the type:


By assuming in (12.18) that c = 0.015, α = 12, β = 4, m1 = 12, m2 = 6, n1 = 8, n2 = 12 and n3 = 6, we obtain the shape of the relevant domain 𝒟 in Figure 39.

Figure 39

Shape of the domain 𝒟.

Let κ = 1.5 be the constant representing the diffusivity and f(x, y) = sinh(xy)+log(x2y2+1) the function representing the initial temperature. In Table 4, the L2(𝒟) and L2(𝒟) norms of the inside and boundary errors κΔu30tu30 and u30 respectively are shown at the times t = 0, 1, 2, 3, where u30 denotes the 30th partial sum of the expansion in Equation (12.14).

κΔu30tu30L2(𝒟) u30L2(𝒟)
t = 0 0.172694 5.87219 × 10−37
t = 1 101.478 5.70500 × 10−48
t = 2 1.48269 × 10−7 5.09531 × 10−58
t = 3 5.87713 × 10−17 5.77811 × 10−68
Table 4

L2 norms of boundary and inside approximation errors at different times.

In Figure 40 are shown, at time t = 0, the approximating solution u30 and the initial temperature f, both expressed in polar coordinates.

Figure 40

The approximating solution u30 and temperature f at time t = 0.

Remark 4.

We note that when the boundary values have wide oscillations, it is necessary to increase the number N of terms in the relevant Fourier expansion in order to obtain better results.

Remark 5.

The L2 norm of the difference between the exact solution and its approximate values is always vanishing in the interior of the considered domain and generally small on the boundary. Point-wise convergence seems to be true on the whole boundary, with the only exception a set of measure zero, corresponding to cusps or quasi-cusped points (i.e. regular points of the curve such that in a very small neighborhood the tangent makes a rotation of almost 180°). In these points, oscillations of the approximate solution (recalling the classical Gibbs phenomenon) usually appear. Therefore, the theoretical results of L. Carleson [22] are confirmed, even in the considered case.

12.4 The Wave Equation in Gielis Domains

Let us consider a membrane with normal polar shape 𝒟 and made from a material characterized by constant propagation speed a. Moreover, suppose the boundary displacement is zero for every t ≥ 0 and the initial displacement and velocity distributions are given by the continuous functions f (x, y) and g (x, y) respectively, so that the problem of finding the displacement at any location within the body for every t > 0 is expressed by:


In [16], with D. Caratelli and P. Natalini as coauthors, the following result was proven:

Theorem 12.3.


and ϵm is Neumann’s symbol [1]. Then the initial-value problem for the wave equation (12.19) admits a classical solution:
such that the following generalized Fourier expansion in terms of Bessel functions holds:
where ζk(m) denotes the kth positive root of the Bessel function of the first type and order m ∈ ℕ0. Imposing the initial conditions U (ϱ*, ϑ, 0) = F (ϱ*, ϑ) and Ut (ϱ*, ϑ, 0) = G (ϱ*, ϑ), the coefficients Am,k, Bm,k, Cm,k, Dm,k are found to be:
with m ∈ ℕ0 and k ∈ ℕ.


In the following example we assume for the boundary 𝒟 a general polar equation of the type:


By assuming in (12.28) that γ1 = γ2 = 3/4, p = q = 7, ν0 = 10, ν1 = ν2 = 6 and ϑ ∈ [0, 2π], the domain 𝒟 features an equisetum-like shape as can be seen in Figures 41 and 42.

Figure 41

Initial distributions of displacement (top) and velocity (bottom) within the equisetum-shaped domain 𝒟 described by the polar equation (12.28) with parameters γ1 = γ2 = 3/4, p = q = 7, ν0 = 10 and ν1 = ν2 = 6.

Figure 42

Spatial distribution of the displacement v (x, y, t) within an equisetum-shaped domain 𝒟 at different times, as predicted by the Fourier expansion representation (12.25) with orders M = K = 60.

Let f(x,y)=log(1+x2y2)12xycos(x+y) and g (x, y) = x3y2 + 3x2y − 2x be the functions describing the initial distributions of displacement and velocity, respectively, within 𝒟 under the hypothesis of normalized propagation constant a = 1. Then, with regard to the relative boundary error eM,K, the numerical results summarized in Table 5 are obtained. In particular, as it appears from Figure 41, the selection of the expansion orders M = K = 60 leads to a very accurate Fourier representation of the solution of the relevant initial-value problem. Finally, we show in Figure 42 the spatial distribution of the displacement v (x, y, t) within the considered domain 𝒟 at different times, as predicted by Equation (12.25) with the mentioned expansion orders.

eM,K M = 0 M = 30 M = 60
K = 1 99.325% 74.383% 74.382%
K = 30 91.050% 15.745% 15.744%
K = 60 90.612% 4.291% 4.239%
Table 5

Relative boundary error eM,K for different expansion orders of the Fourier-like solution of the initial-value problem for the wave equation (12.19) within the domain 𝒟 described by the polar equation (12.28) with parameters γ1 = γ2 = 3/4, p = q = 7, ν0 = 10 and ν1 = ν2 = 6.