Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences
Open AccessResearch articles

Interaction of scales for a singularly perturbed degenerating nonlinear Robin problem

Paolo Musolino

Paolo Musolino

Dipartimento di Scienze Molecolari e Nanosistemi, Università Ca’ Foscari Venezia, via Torino 155, 30172 Venezia Mestre, Italy

[email protected]

Contribution: Conceptualization, Writing – original draft, Writing – review & editing

Google Scholar

Find this author on PubMed

and
Gennady Mishuris

Gennady Mishuris

Department of Mathematics, Aberystwyth University, Ceredigion, Aberystwyth SY23 3BZ, UK

Contribution: Conceptualization, Writing – original draft, Writing – review & editing

Google Scholar

Find this author on PubMed

Published:https://doi.org/10.1098/rsta.2022.0159

    Abstract

    We study the asymptotic behaviour of solutions of a boundary value problem for the Laplace equation in a perforated domain in Rn, n3, with a (nonlinear) Robin boundary condition on the boundary of the small hole. The problem we wish to consider degenerates in three respects: in the limit case, the Robin boundary condition may degenerate into a Neumann boundary condition, the Robin datum may tend to infinity, and the size ϵ of the small hole where we consider the Robin condition collapses to 0. We study how these three singularities interact and affect the asymptotic behaviour as ϵ tends to 0, and we represent the solution and its energy integral in terms of real analytic maps and known functions of the singular perturbation parameters.

    This article is part of the theme issue ‘Non-smooth variational problems and applications’.

    1. Introduction

    This article is devoted to the study of the asymptotic behaviour of solutions of a boundary value problem for the Laplace equation in a perforated domain in Rn, n3, with a (nonlinear) Robin boundary condition which degenerates into a Neumann condition on the boundary of the small hole. The problem we wish to consider degenerates in three respects. First, in the limit case, the Robin boundary condition may degenerate into a Neumann boundary condition (i.e. the coefficient of the trace of the solution in the boundary condition may vanish). Second, the Robin datum may tend to infinity. Finally, the size ϵ of the small hole where we consider the Robin condition approaches the degenerate value ϵ=0.

    The behaviour of solutions to boundary value problems with degenerating or perturbed boundary conditions has been studied by many authors. A family of Poincaré problems approximating a mixed boundary value problem for the Laplace equation in the plane has been studied by Wendland et al. [1]. A study of the convergence of the solution of the Helmholtz equation with boundary condition of the type ϵ(u/ν)+u=g to the solution with Dirichlet condition u=g as ϵ0 can be found in the work of Kirsch [2]. Costabel and Dauge [3] studied a mixed Neumann–Robin problem for the Laplace operator, where the Robin condition tends to a Dirichlet condition as the perturbation parameter tends to 0. Boundary value problems for Maxwell equations with singularly perturbed boundary conditions have been analysed, for instance, by Ammari and Nédélec [4]. Also, singularly perturbed transmission problems have been investigated by Schmidt and Hiptmair [5] by means of integral equation methods. Dalla Riva and Mishuris [6] have investigated the solvability of a small nonlinear perturbation of a homogeneous linear transmission problem by using potential-theoretical techniques. The present article represents a continuation of the analysis done in [7], where the authors considered the behaviour as δ0 of the solutions to the boundary value problem

    {Δu(x)=0xΩoΩi¯,νΩou(x)=go(x)xΩo,νΩiu(x)=δFδ(u(x))+gi(x)xΩi,1.1
    where Ωo and Ωi are sufficiently regular bounded open sets such that Ωi¯Ωo. In this equation, the superscript ‘o’ stands for ‘outer domain’ and the superscript ‘i’ stands for ‘inner domain.’ The problem generalizes a linear problem that, under suitable assumptions, admits a unique solution uδ for each δ>0. When δ=0, the problem degenerates into the Neumann problem
    {Δu(x)=0xΩoΩi¯,νΩou(x)=go(x)xΩo,νΩiu(x)=gi(x)xΩi.1.2
    As is well known, this Neumann problem may have infinite solutions or no solutions, depending on compatibility conditions on the Neumann datum. In [7] we proved that, under suitable assumptions, solutions to (1.1) exist and that they diverge if the compatibility condition on the Neumann datum for the existence of solutions to (1.2) does not hold. In [7], we considered a Robin problem as a simplified model for the transmission problem for a composite domain with imperfect (non-natural) conditions along the joint boundary. Such nonlinear transmission conditions frequently appear in practical applications for various nonlinear multiphysics problems (e.g. [815]). All such transmission conditions have been derived using formal variational or asymptotic techniques (see e.g. [1618]). However, accurate analysis of their solvability and solution regularity has not been performed. One of the aims of the present paper is to address this need. On the other hand, the problem in question, though simpler than most of those arising in applications, is rich enough as it contains some features influencing the final result. This refers not only to the condition itself but also the surface on which they hold.

    In [7], we considered the case where the surface on which we consider the Robin condition is the boundary of a fixed hole Ωi. Here, we wish to study the case where the hole becomes small and degenerates into a point. Then a natural question arises: if we replace the set Ωi by a small set ϵωi (with ϵ close to 0) and the parameter δ by a function δ(ϵ) possibly tending to zero as ϵ0, what happens? How does the geometric degeneracy (the set ϵωi collapsing to the origin when ϵ=0) interact with the possible degeneracy of the boundary condition if δ(ϵ)0 as ϵ approaches 0? We also observe that even though several techniques are available for the analysis of linear problems, the presence of a nonlinear boundary condition requires a specific type of analysis since, for example, existence and uniqueness of solutions is not immediately ensured.

    The purpose of the present article is to give answers to these questions. Here we consider only the case of dimension n3. Indeed, our technique is based on potential theory, and the two-dimensional case requires a specific analysis due to different aspects of the fundamental solution of the Laplacian. In particular, if n=2 the fundamental solution Sn of the Laplacian equals (log|x|)/(2π), whereas if n3 the fundamental solution Sn is a multiple of 1/|x|n2. This leads to different rescaling behaviour of Sn(ϵx) and to different behaviour at infinity of single-layer potentials, which are among our main tools in the analysis. We note that the set ϵωi when ϵ is close to zero can be seen as a small hole in the set Ωo. The behaviour of the solutions to boundary value problems in domains with small holes has long been investigated by the expansion methods of asymptotic analysis. Such methods are mainly based on elliptic theory and allow the treatment of a large variety of linear problems. As examples, we mention the method of matching outer and inner asymptotic expansions of Il’in [19] and the compound asymptotic expansion method of Maz’ya et al. [20,21], which allows the treatment of general Douglis–Nirenberg elliptic boundary value problems in domains with perforations and corners. More recently, Maz’ya et al. [22] provided asymptotic analysis of Green’s kernels in domains with small cavities by applying the method of mesoscale asymptotic approximations (see also the papers [2327]). Moreover, we refer to Ammari and Kang [28] for several applications to inverse problems and Novotny and Sokołowski [29] for applications to topological optimization.

    Instead of the methods of asymptotic analysis, here we exploit the so-called functional-analytic approach proposed by Lanza de Cristoforis in [30]. The goal of this approach is to represent solutions to problems in perturbed domains in terms of real analytic maps and known functions of the perturbation parameter. For a detailed presentation of the functional-analytic approach, we refer to Dalla Riva et al. [31]. Here, we mention that the functional-analytic approach has been used to analyse a nonlinear Robin problem for the Laplace equation by Lanza de Cristoforis [32] and Lanza de Cristoforis and Musolino [33] and to analyse nonlinear traction problems for Lamé equations by, for example, Dalla Riva and Lanza de Cristoforis [3436] and Falconi et al. [37].

    As a first step, we introduce the geometric setting in which we are going to consider our boundary value problem. As the dimensional parameter, we take a natural number

    nN{0,1,2}.
    Then, to define the perforated domain, we consider a regularity parameter α(0,1) and two subsets ωi and Ωo of Rn satisfying the following condition:
    ωiand Ωoare bounded open connected subsets of Rnof class C1,αsuch that 0Ωoωiand both Rnωi¯and RnΩo¯are connected.
    The set Ωo plays the role of the unperturbed domain, whereas the set ωi represents the shape of the perforation. We refer, for instance, to Gilbarg and Trudinger [38] for the definition of sets and functions of the Schauder class Ck,α (kN). We fix
    ϵ0sup{θ(0,+):ϵωi¯Ωoϵ(θ,θ)}.
    We note that if ϵ(0,ϵ0), the set ϵωi¯ (which we think as a hole) is contained in Ωo and therefore we can remove it from the unperturbed domain. We define the perforated domain Ω(ϵ) by setting
    Ω(ϵ)Ωoϵωi¯
    for all ϵ(0,ϵ0). When ϵ approaches zero, the set Ω(ϵ) degenerates to the punctured domain Ωo{0}. Clearly, the boundary Ω(ϵ) of Ω(ϵ) consists of the two connected components Ωo and (ϵωi)=ϵωi. Therefore we can identify, for example, C0,α(Ω(ϵ)) with the product C0,α(Ωo)×C0,α(ϵωi). Moreover, after a suitable rescaling, we can identify functions in C0,α(ϵωi) with functions in C0,α(ωi). Having introduced the geometric aspects of our problem, we need to define the boundary data. To do this, we fix two functions
    goC0,α(Ωo)andgiC0,α(ωi).
    Then we take a family {Fϵ}ϵ]0,ϵ0[ of functions from R to R and two functions δ() and ρ() from (0,ϵ0) to (0,+). As we shall see, the function go represents the Neumann datum on the exterior boundary Ωo. The family of functions {Fϵ}ϵ]0,ϵ0[ will allow us to define the nonlinear Robin condition on ϵωi, and δ(ϵ) will be the coefficient of a function of the Dirichlet trace in the Robin condition. We will consider a non-homogeneous Robin condition, and thus the corresponding datum will be gi(/ϵ)/ρ(ϵ). Next, for each ϵ(0,ϵ0), we want to consider a nonlinear boundary value problem for the Laplace operator. Namely, we consider a Neumann condition on Ωo and a nonlinear Robin condition on ϵωi. Thus, for each ϵ(0,ϵ0) we consider the following problem:
    {Δu(x)=0xΩ(ϵ),νΩou(x)=go(x)xΩo,νϵωiu(x)=δ(ϵ)Fϵ(u(x))+gi(x/ϵ)ρ(ϵ)xϵωi,1.3
    where νΩo and νϵωi denote the outward unit normals to Ωo and to (ϵωi), respectively. Our aim is to analyse the behaviour of the solutions to problem (1.3) as ϵ0. As we have already mentioned, when ϵ tends to 0 the hole ϵωi degenerates into the origin 0. Moreover, if δ(ϵ) tends to 0 as ϵ0, the Robin condition may degenerate into a Neumann condition. Furthermore, we will also allow the term ρ(ϵ) to tend to 0, which may generate a further singularity. An aspect we wish to highlight in the present article is how all these singularities interact together. Our main results are represented by theorems 4.4 and 4.6, which describe in detail the asymptotic behaviour of the solutions as ϵ0, and theorem 4.7, which concerns the behaviour of the energy integrals of the solutions. These results highlight the interactions of different scales. Moreover, as we will see, it will be crucial to assume that the quantities ϵδ(ϵ) and ϵn1/ρ(ϵ) have limits as ϵ0. Incidentally, we observe that interactions of scales are well known to possibly cause strange phenomena in the limiting behaviour of solutions. As an example, we mention the celebrated works of Cioranescu and Murat [39,40] and of Marčenko and Khruslov [41], as well as the more recent articles of Arrieta and Lamberti [42], Arrieta et al. [43] and Ferraresso and Lamberti [44]. We also mention the work of Bonnetier et al. [45] concerning small perturbations in the type of boundary conditions and of Felli et al. [46,47] on disappearing Neumann or Dirichlet regions in mixed eigenvalue problems.

    We observe that in the present article, the boundary of the hole depends on ϵ simply through a dilation. However, in the literature one can find examples where the geometry changes in a more drastic way, as in the case of oscillating boundaries (see e.g. [4244]). On the other hand, one may also consider the case where the geometry is fixed and the boundary condition is changing, as in [7].

    This article is organized as follows. In §2, we analyse a toy problem in an annular domain. In §3, we transform problem (1.3) into an equivalent system of integral equations. In §4, we analyse this system and prove our main results on the asymptotic behaviour of a family of solutions and the corresponding energy integrals.

    2. A toy problem

    As we have done in [7], we consider problem (1.3) in the annular domain

    Ω(ϵ)Bn(0,1)Bn(0,ϵ)¯=Bn(0,1)ϵBn(0,1)¯,
    i.e. we take ΩoBn(0,1) and ωiBn(0,1), where, for r>0, the symbol Bn(0,r) denotes the open ball in Rn of centre 0 and radius r. We will then set ϵ0=1, Fϵ(τ)=τ for all τR and all ϵ(0,ϵ0), go=a and gi=b, where a,bR. Moreover, we consider two functions δ,ρ:(0,1)(0,+). Then for each ϵ(0,1) we consider the problem
    {Δu(x)=0xBn(0,1)Bn(0,ϵ)¯,νBn(0,1)u(x)=axBn(0,1),νBn(0,ϵ)u(x)=δ(ϵ)u(x)+bρ(ϵ)xBn(0,ϵ).2.1
    As is well known, for each ϵ(0,1) a solution uϵC1,α(Ω(ϵ)¯) to problem (2.1) exists and is unique (see [31], theorem 6.56). On the other hand, if instead we put ϵ=0 in (2.1), the hole disappears and we are led to consider the Neumann problem
    {Δu(x)=0xBn(0,1),νBn(0,1)u(x)=axBn(0,1).2.2
    As with any Neumann problem, the solvability of (2.2) is subject to compatibility conditions on the Neumann datum on Bn(0,1). In this specific case of constant Neumann datum, problem (2.2) has a solution if and only if a=0. Obviously, if a=0, then the Neumann problem (2.2) has a one-dimensional space of solutions, which consists of the space of constant functions in Bn(0,1)¯; if instead a0, problem (2.2) does not have any solution. On the other hand, if one considers the behaviour of the unique solution uϵ of problem (2.1), the earlier remark clearly implies that in general uϵ cannot converge to a solution of (2.2) as ϵ0 if the compatibility condition a=0 does not hold. Also, even if a=0, we shall see that the solutions may diverge as ϵ0, depending on the behaviour of the functions δ(ϵ) and ρ(ϵ) for ϵ close to 0. Moreover, one would like to understand how the behaviour of δ(ϵ) and ρ(ϵ) affects the asymptotic behaviour of uϵ and whether there is a ‘memory’ of the Robin condition. In the specific case of our annular domain and constant data, we can construct explicitly the solution uϵ. Then we try to understand the behaviour of uϵ as ϵ0. To construct explicitly uϵ, we search for a solution of (2.1) in the form
    uϵ(x)Aϵ1(2n)|x|n2+BϵxΩ(ϵ)¯,
    with Aϵ and Bϵ chosen so that the boundary conditions of problem (2.1) are satisfied. By a straightforward computation, we must have
    uϵ(x)a1(2n)|x|n2+1δ(ϵ)(aϵn1bρ(ϵ))+a(n2)ϵn2xΩ(ϵ)¯.2.3
    We now note that we can rewrite equation (2.3) as
    uϵ(x)a1(2n)|x|n2+1δ(ϵ)ϵn1(abϵn1ρ(ϵ)+a(n2)δ(ϵ)ϵ)xΩ(ϵ)¯.2.4
    In particular, if d0limϵ0ϵδ(ϵ)R, r0limϵ0ϵn1/ρ(ϵ)R and abr0+ad0/(n2)0, then uϵ(x) is asymptotic to (abr0+ad0/(n2))/(ϵn1δ(ϵ)) as ϵ tends to 0, when x is fixed in Bn(0,1)¯{0}. In conclusion, under suitable assumptions on the behaviour of δ(ϵ) and ρ(ϵ) as ϵ0, we see that the value of the solution uϵ at a fixed point xBn(0,1)¯{0} behaves like 1/(ϵn1δ(ϵ)) and there is some sort of interaction of scales influencing the limiting behaviour of the solution. Similarly, if one considers the energy integral of uϵ, a direct computation shows that
    Ω(ϵ)|uϵ(x)|2dx=Ω(ϵ)|(a1(2n)|x|n2)|2dx=Ω(ϵ)a21|x|2n2dx=a2snϵ11rn1dr=a2sn(n2)1ϵn2(1ϵn2),2.5
    where the symbol sn denotes the (n1)-dimensional measure of Bn(0,1). In particular, if a0, the energy integral of the solution Ω(ϵ)|uϵ(x)|2dx behaves like 1/ϵn2.

    Our aim is to recover and understand such behaviour of the solution and of its energy integral in a more general situation, for both the geometry and the boundary conditions. Indeed, we will show that the main features discussed above can be identified in the general solution (compare (2.4) with (4.14) and (2.5) with (4.17)). We emphasize that one can derive a uniform asymptotic solution by the methods of [20,21,26]. In particular, one can identify the uniform limit far from the hole as a first approximation. Then one can correct such a limit in order to improve the approximation on rescaled sets and repeat the procedure in order to reduce the error. We will show in remarks 4.2 and 4.5 how one can deduce the above considerations from the results of §4 (which can thus be seen as analogous formulas in more general settings).

    3. An integral equation formulation of the boundary value problem

    As in [7], to analyse problem (1.3) for ϵ close to 0, we exploit the so-called functional-analytic approach (see [31]). This method is based on classical potential theory, which allows one to obtain an integral equation formulation of (1.3). First we need to introduce some notation. We denote by Sn the function from Rn{0} to R defined by

    Sn(x)1(2n)sn|x|n2xRn{0}.
    Since n3, as is well known, Sn is a fundamental solution of the Laplace operator. By means of the fundamental solution Sn, we construct some integral operators (namely, single-layer potentials) that we use to represent harmonic functions (and thus, in particular, the solutions of problem (1.3)). So let Ω be a bounded open connected subset of Rn of class C1,α. If μC0(Ω), we introduce the single-layer potential by setting
    v[Ω,μ](x)ΩSn(xy)μ(y)dσyxRn,
    where dσ denotes the area element of a manifold imbedded in Rn. It is well known that if μC0(Ω), then v[Ω,μ] is continuous in Rn. Moreover, if μC0,α(Ω), then the function v+[Ω,μ]v[Ω,μ]|Ω¯ belongs to C1,α(Ω¯) and the function v[Ω,μ]v[Ω,μ]|RnΩ belongs to Cloc1,α(RnΩ). The normal derivative of the single-layer potential on Ω, on the other hand, exhibits a jump. To describe the jump, we set
    W[Ω,μ](x)ΩνΩ(x)Sn(xy)μ(y)dσyxΩ,
    where νΩ denotes the outward unit normal to Ω. If μC0,α(Ω), the function W[Ω,μ] belongs to C0,α(Ω), and we have
    νΩv±[Ω,μ]=12μ+W[Ω,μ]on Ω.
    As we shall see in lemma 3.1, to represent the functions on Ω(ϵ)¯ which are harmonic and satisfy the boundary conditions, we will exploit single-layer potentials having densities with zero integral mean on Ωo plus constants. Therefore, we find it convenient to set
    C0,α(Ωo)0{fC0,α(Ωo):Ωofdσ=0}.
    More precisely, in lemma 3.1 we represent a function uC1,α(Ω(ϵ)¯) such that Δu=0 in Ω(ϵ) as a single-layer potential plus the ϵ-dependent constant ξ/(δ(ϵ)ϵn1). The reason for the choice of such a constant is that in view of the results of §2, we expect the presence of a constant behaving like 1/(δ(ϵ)ϵn1) as ϵ0 in the representation formula for the solutions of (1.3). The proof of lemma 3.1 can be deduced from classical potential theory (cf. Folland [48], ch. 3, and Dalla Riva et al. [31], proof of proposition 6.49).

    Lemma 3.1.

    Let ϵ(0,ϵ0). Let uC1,α(Ω(ϵ)¯) be such that Δu=0 in Ω(ϵ). Then there exists a unique triple (μo,μi,ξ)C0,α(Ωo)0×C0,α(ωi)×R such that

    u(x)=ΩoSn(xy)μo(y)dσy+ωiSn(xϵs)μi(s)dσs+ξδ(ϵ)ϵn1xΩ(ϵ)¯.

    By exploiting lemma 3.1, we can establish a correspondence between the solutions of problem (1.3) and those of a (nonlinear) system of integral equations.

    Proposition 3.2.

    Let ϵ(0,ϵ0). Then the map from the set of triples (μo,μi,ξ)C0,α(Ωo)0×C0,α(ωi)×R such that

    12μo(x)+ΩoνΩo(x)Sn(xy)μo(y)dσy+ωiνΩo(x)Sn(xϵs)μi(s)dσs=go(x)xΩo3.1
    and
    12μi(t)+ϵn1Ωoνωi(t)Sn(ϵty)μo(y)dσy+ωiνωi(t)Sn(ts)μi(s)dσs=ϵn1δ(ϵ)Fϵ(ΩoSn(ϵty)μo(y)dσy+1ϵn2ωiSn(ts)μi(s)dσs+ξδ(ϵ)ϵn1)+gi(t)ϵn1ρ(ϵ)tωi3.2
    to the set of functions uC1,α(Ω(ϵ)¯) that solve problem (1.3), which takes a triple (μo,μi,ξ) to
    ΩoSn(xy)μo(y)dσy+ωiSn(xϵs)μi(s)dσs+ξδ(ϵ)ϵn1xΩ(ϵ)¯,3.3
    is a bijection.

    Proof.

    If (μo,μi,ξ)C0,α(Ωo)0×C0,α(ωi)×R, then we know that the function in (3.3) belongs to C1,α(Ω(ϵ)¯) and is harmonic in Ω(ϵ). Moreover, if (μo,μi,ξ) satisfies system (3.1)–(3.2), then the jump formula for the normal derivative of the single-layer potential implies the validity of the boundary condition in problem (1.3). Hence, the function in (3.3) solves problem (1.3).

    Conversely, if uC1,α(Ω(ϵ)¯) satisfies (1.3), then lemma 3.1 for harmonic functions ensures that there exists a unique triple (μo,μi,ξ)C0,α(Ωo)0×C0,α(ωi)×R such that

    u(x)=ΩoSn(xy)μo(y)dσy+ωiSn(xϵs)μi(s)dσs+ξδ(ϵ)ϵn1xΩ(ϵ)¯.
    Then the formula for the normal derivative of a single-layer potential and the boundary conditions in (1.3) imply that (3.1) and (3.2) are satisfied. Hence, the map of the statement is a bijection.

    Now that the correspondence between the solutions of boundary value problem (1.3) and those of the system of integral equations (3.1)–(3.2) is established, we wish to study the behaviour of the solutions to system (3.1)–(3.2) as ϵ0. Note that if ϵ(0,ϵ0), we can write

    ϵn1δ(ϵ)Fϵ(ΩoSn(ϵty)μo(y)dσy+1ϵn2ωiSn(ts)μi(s)dσs+ξδ(ϵ)ϵn1)=ϵn1δ(ϵ)Fϵ(1ϵn1δ(ϵ)(ϵn1δ(ϵ)ΩoSn(ϵty)μo(y)dσy+ϵδ(ϵ)ωiSn(ts)μi(s)dσs+ξ))tωi.

    We now wish to analyse equation (3.2) for ϵ small. As in [7], we need to make a structural assumption on the nonlinearity, i.e. on the family of functions Rτϵn1δ(ϵ)Fϵ(τ/(ϵn1δ(ϵ))) for ϵ close to 0. So we assume the following:

    there exist ϵ1(0,ϵ0),mN, a real analytic function F~:Rm+1Rand a function η:(0,ϵ1)Rm such that η0=limϵ0η(ϵ)Rm and ϵn1δ(ϵ)Fϵ(1ϵn1δ(ϵ)τ)=F~(τ,η(ϵ))for all (τ,ϵ)R×(0,ϵ1).3.4

    As a simple example, one can take as Fϵ a small polynomial perturbation of the identity. For example,

    Fϵ(z)=z+h(ϵ)zm,
    where mN{0,1} and h is a certain function from (0,ϵ1) to R. Then we have
    ϵn1δ(ϵ)Fϵ(1ϵn1δ(ϵ)τ)=τ+h(ϵ)(ϵn1δ(ϵ))m1τm.
    If
    limϵ0h(ϵ)(ϵn1δ(ϵ))m1R,
    then one has
    η(ϵ)=h(ϵ)(ϵn1δ(ϵ))m1andF~(τ,η(ϵ))=τ+η(ϵ)τm=τ+h(ϵ)(ϵn1δ(ϵ))m1τm.
    On the other hand, one could also construct Fϵ starting from a given F~ and η(ϵ). This would allow one to generate more involved nonlinearities (if perhaps less natural).

    Here, we observe that different structures of the nonlinearity may be tackled by modifying our approach. Although the type of nonlinearity we consider is quite specific, we emphasize that our techniques are not confined to linear boundary conditions and apply also in some nonlinear cases. At the same time, we are also interested in the linear case, since the degeneracy appears there as well. Therefore, for us, it is enough to include some (nonlinear) perturbations of the linear case.

    4. Analytic representation formulas for the solution of the boundary value problem

    We observe that under the additional assumption (3.4), equations (3.1) and (3.2) take the forms

    12μo(x)+ΩoνΩo(x)Sn(xy)μo(y)dσy+ωiνΩo(x)Sn(xϵs)μi(s)dσs=go(x)xΩo4.1
    and
    12μi(t)+ϵn1Ωoνωi(t)Sn(ϵty)μo(y)dσy+ωiνωi(t)Sn(ts)μi(s)dσs=F~(ϵn1δ(ϵ)ΩoSn(ϵty)μo(y)dσy+ϵδ(ϵ)ωiSn(ts)μi(s)dσs+ξ,η(ϵ))+gi(t)ϵn1ρ(ϵ)tωi,4.2
    for all ϵ(0,ϵ1). We would like to pass to the limit as ϵ0 in equations (4.1) and (4.2). However, to do so, we need to know the asymptotic behaviour for ϵ close to 0 of the quantities ϵδ(ϵ) and ϵn1ρ(ϵ), which appear in (4.2). Accordingly, we now assume that
    d0limϵ0ϵδ(ϵ)Randr0limϵ0ϵn1ρ(ϵ)R.4.3

    Motivated by (4.1) and (4.2), we replace the quantities ϵδ(ϵ), η(ϵ) and ϵn1/ρ(ϵ) by the auxiliary variables γ1, γ2 and γ3, respectively, and we now introduce the operator Λn(Λno,Λni) from (ϵ1,ϵ1)×Rm+2×C0,α(Ωo)0×C0,α(ωi)×R to C0,α(Ωo)×C0,α(ωi) defined by

    Λno[ϵ,γ1,γ2,γ3,μo,μi,ξ](x)12μo(x)+ΩoνΩo(x)Sn(xy)μo(y)dσy+ωiνΩo(x)Sn(xϵs)μi(s)dσsgo(x)xΩo4.4
    and
    Λni[ϵ,γ1,γ2,γ3,μo,μi,ξ](t)12μi(t)+ϵn1Ωoνωi(t)Sn(ϵty)μo(y)dσy+ωiνωi(t)Sn(ts)μi(s)dσsF~(ϵn2γ1ΩoSn(ϵty)μo(y)dσy+γ1ωiSn(ts)μi(s)dσs+ξ,γ2)gi(t)γ3tωi4.5
    for all (ϵ,γ1,γ2,γ3,μo,μi,ξ)(ϵ1,ϵ1)×Rm+2×C0,α(Ωo)0×C0,α(ωi)×R.

    Then, if ϵ(0,ϵ1), in view of definitions (4.4) and (4.5), the system of equations

    Λno[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ),μo,μi,ξ](x)=0xΩo4.6
    and
    Λni[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ),μo,μi,ξ](t)=0tωi4.7
    is equivalent to the system (4.1)–(4.2). Then if we let ϵ0 in (4.6) and (4.7), we obtain
    12μo(x)+ΩoνΩo(x)Sn(xy)μo(y)dσy+νΩo(x)Sn(x)ωiμi(s)dσs=go(x)xΩo4.8
    and
    12μi(t)+ωiνωi(t)Sn(ts)μi(s)dσs=F~(d0ωiSn(ts)μi(s)dσs+ξ,η0)+gi(t)r0tωi.4.9

    Now we would like to prove for ϵ(0,ϵ1) the existence of solutions (μo,μi,ξ) to the system (4.6)–(4.7) around a solution of the limiting system (4.8)–(4.9). Therefore, we now further assume that

    the system (4.8)(4.9) in the unknown (μo,μi,ξ) admits a solution (μ~o,μ~i,ξ~) in C0,α(Ωo)0×C0,α(ωi)×R.4.10
    We do not discuss here conditions on F~ ensuring the existence of a solution of system (4.8)–(4.9). However, they can be obtained by arguing as in [32], appendix C, or in [7].

    We note that if (μ~o,μ~i,ξ~) is a solution of the system (4.8)–(4.9), then by integrating (4.8) on Ωo and using the equalities

    ΩoΩoνΩo(x)Sn(xy)μ~o(y)dσydσx=12Ωoμ~o(y)dσy
    (cf. [31], lemma 6.11) and ΩoνΩo(x)Sn(x)dσx=1 (cf. [31], corollary 4.6), we obtain ωiμ~i(s)dσs=Ωogo(x)dσx. In the following proposition, we investigate the system of integral equations (4.1)–(4.2) by applying the implicit function theorem to Λn under suitable assumptions on τF~(d0ωiSn(ts)μ~i(s)dσs+ξ~,η0), where τF~ denotes the partial derivative with respect to the variable τ of the function (τ,η)F~(τ,η).

    Proposition 4.1.

    Let assumptions (3.4) and (4.3) hold. Let (μ~o,μ~i,ξ~) be as in (4.10). Assume that

    ωiτF~(d0ωiSn(ts)μ~i(s)dσs+ξ~,η0)dσt0
    and, if d00, also
    τF~(d0ωiSn(ts)μ~i(s)dσs+ξ~,η0)0tωi.
    Then there exist ϵ2(0,ϵ1), an open neighbourhood U of (d0,η0,r0) in Rm+2, an open neighbourhood V of (μ~o,μ~i,ξ~) in C0,α(Ωo)0×C0,α(ωi)×R and a real analytic map (Mo,Mi,Ξ) from (ϵ2,ϵ2)×U to V such that
    (ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ))Uϵ(0,ϵ2)
    and such that the set of zeros of Λn in (ϵ2,ϵ2)×U×V coincides with the graph of (Mo,Mi,Ξ). In particular, (Mo[0,d0,η0,r0],Mi[0,d0,η0,r0],Ξ[0,d0,η0,r0])=(μ~o,μ~i,ξ~).

    Proof.

    From standard results of classical potential theory (see e.g. Dalla Riva et al. [31], Miranda [49] and Lanza de Cristoforis and Rossi [50]), real analyticity results for integral operators with real analytic kernel (see Lanza de Cristoforis and Musolino [51]), assumption (3.4) and real analyticity results for the composition operator (Böhme and Tomi [52], p. 10, Henry [53] and Valent [54], theorem 5.2), we deduce that Λn is a real analytic operator from (ϵ1,ϵ1)×Rm+2×C0,α(Ωo)0×C0,α(ωi)×R to C0,α(Ωo)×C0,α(ωi). By standard calculus in Banach spaces, we verify that the partial differential (μo,μi,ξ)Λn[0,d0,η0,r0,μ~o,μ~i,ξ~] of Λn at (0,d0,η0,r0,μ~o,μ~i,ξ~) with respect to the variable (μo,μi,ξ) is given by

    (μo,μi,ξ)Λno[0,d0,η0,r0,μ~o,μ~i,ξ~](μ¯o,μ¯i,ξ¯)(x)12μ¯o(x)+ΩoνΩo(x)Sn(xy)μ¯o(y)dσy+νΩo(x)Sn(x)ωiμ¯i(s)dσsxΩo
    and
    (μo,μi,ξ)Λni[0,d0,η0,r0,μ~o,μ~i,ξ~](μ¯o,μ¯i,ξ¯)(t)12μ¯i(t)+ωiνωi(t)Sn(ts)μ¯i(s)dσsτF~(d0ωiSn(ts)μ~i(s)dσs+ξ~,η0)(d0ωiSn(ts)μ¯i(s)dσs+ξ¯)tωi
    for all (μ¯o,μ¯i,ξ¯)C0,α(Ωo)0×C0,α(ωi)×R. Now we want to show that the partial differential (μo,μi,ξ)Λn[0,d0,η0,r0,μ~o,μ~i,ξ~] is a homeomorphism from C0,α(Ωo)0×C0,α(ωi)×R onto C0,α(Ωo)×C0,α(ωi). Since (μo,μi,ξ)Λn[0,d0,η0,r0,μ~o,μ~i,ξ~] is the sum of an invertible operator and a compact operator, one immediately verifies that it is a Fredholm operator of index 0. Therefore, to prove that the operator (μo,μi,ξ)Λn[0,d0,η0,r0,μ~o,μ~i,ξ~] is a homeomorphism, it suffices to prove that it is injective. So let us assume that
    (μo,μi,ξ)Λn[0,d0,η0,r0,μ~o,μ~i,ξ~](μ¯o,μ¯i,ξ¯)=0.
    By integrating on the Ωo equality
    (μo,μi,ξ)Λno[0,d0,η0,r0,μ~o,μ~i,ξ~](μ¯o,μ¯i,ξ¯)(x)=0xΩo
    and using the equalities
    ΩoΩoνΩo(x)Sn(xy)μ¯o(y)dσydσx=12Ωoμ¯o(y)dσy
    (cf. [31], lemma 6.11) and ΩoνΩo(x)Sn(x)dσx=1 (cf. [31], corollary 4.6), we obtain
    ωiμ¯i(s)dσs=0.4.11
    As a consequence,
    12μ¯o(x)+ΩoνΩo(x)Sn(xy)μ¯o(y)dσy=0xΩo,
    and thus by [31], theorem 6.25, since Ωoμ¯odσ=0 we have μ¯o=0. By (4.11) and the same argument as in [33], proof of theorem 4.4, the equality
    (μo,μi,ξ)Λni[0,d0,η0,r0,μ~o,μ~i,ξ~](μ¯o,μ¯i,ξ¯)(t)=0tωi
    implies that (μ¯i,ξ¯)=0. In conclusion, we have shown that (μo,μi,ξ)Λn[0,d0,η0,r0,μ~o,μ~i,ξ~] is injective and thus, being a Fredholm operator of index 0, also a homeomorphism. As a consequence, we can apply the implicit function theorem for real analytic maps in Banach spaces (cf. [55], theorem 15.3) and deduce that there exist ϵ2(0,ϵ1), an open neighbourhood U of (d0,η0,r0) in Rm+2, an open neighbourhood V of (μ~o,μ~i,ξ~) in C0,α(Ωo)0×C0,α(ωi)×R and a real analytic map (Mo,Mi,Ξ) from (ϵ2,ϵ2)×U to V with (ϵδ(ϵ),η(ϵ),(ϵn1/ρ(ϵ)))U for all ϵ(0,ϵ2) such that the set of zeros of Λn in (ϵ2,ϵ2)×U×V coincides with the graph of (Mo,Mi,Ξ) and, in particular, (Mo[0,d0,η0,r0],Mi[0,d0,η0,r0],Ξ[0,d0,η0,r0])=(μ~o,μ~i,ξ~).

    Remark 4.2.

    If F is linear, system (4.8)–(4.9) simplifies to

    12μo(x)+ΩoνΩo(x)Sn(xy)μo(y)dσy+νΩo(x)Sn(x)ωiμi(s)dσs=go(x)xΩo4.12
    and
    12μi(t)+ωiνωi(t)Sn(ts)μi(s)dσs=d0ωiSn(ts)μi(s)dσs+ξ+gi(t)r0tωi.4.13
    Then, by arguing as in the proof of proposition 4.1, one verifies that system (4.12)–(4.13) in the unknown (μo,μi,ξ) admits a unique solution (μ~o,μ~i,ξ~) in C0,α(Ωo)0×C0,α(ωi)×R. By integrating (4.12) and (4.13), we deduce that
    1ωidσ(Ωogodσd0ωiωiSn(ts)μ~i(s)dσsdσtr0ωigidσ)=ξ~.
    If we further assume that
    Ωo=ωi=Bn(0,1),go(x)=axBn(0,1)andgi(t)=btBn(0,1)
    for some constants a,bR, then by the well-known identity
    Bn(0,1)Sn(ts)dσt=12nsBn(0,1),
    one obtains
    1sn(asnd012nasnbsnr0)=ξ~
    and thus
    ξ~=abr0+an2d0.

    Now that we have converted (1.3) into a system of integral equations for which we have exhibited a real analytic family of solutions, we introduce a family of solutions to (1.3).

    Definition 4.3.

    Let the assumptions of proposition 4.1 hold. Then we set

    u(ϵ,x)=ΩoSn(xy)Mo[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](y)dσy+ωiSn(xϵs)Mi[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](s)dσs+Ξ[ϵ,ϵδ(ϵ),η(ϵ),(ϵn1/ρ(ϵ))]δ(ϵ)ϵn1xΩ(ϵ)¯,ϵ(0,ϵ2).

    By propositions 3.2 and 4.1 and definition 4.3, we deduce that for each ϵ(0,ϵ2) the function u(ϵ,)C1,α(Ω(ϵ)¯) is a solution to problem (1.3). In the following theorems, we exploit the analyticity result of proposition 4.1 to prove representation formulas for u(ϵ,) and for its energy integral in terms of real analytic maps. We start with the following theorem, which considers the restriction of the solution u(ϵ,) to a set that is ‘far’ from the hole.

    Theorem 4.4.

    Let the assumptions of proposition 4.1 hold. Let ΩM be a bounded open subset of Ωo such that 0ΩM¯. Then there exist ϵM(0,ϵ2) and a real analytic map UM from (ϵM,ϵM)×U to C1,α(ΩM¯) such that ΩM¯Ω(ϵ)¯ for all ϵ(0,ϵM) and

    u(ϵ,x)=UM[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](x)+Ξ[ϵ,ϵδ(ϵ),η(ϵ),(ϵn1/ρ(ϵ))]δ(ϵ)ϵn1xΩM¯4.14
    for all ϵ(0,ϵM). Moreover, if we set
    u~M(x)ΩoSn(xy)μ~o(y)dσyxΩo¯,
    we have that UM[0,d0,η0,r0]=u~M|ΩM¯+Sn|ΩM¯Ωogodσ, and u~M solves the Neumann problem
    {Δu(x)=0xΩo,νΩou(x)=go(x)νΩoSn(x)ΩogodσxΩo.4.15

    Proof.

    Taking ϵM(0,ϵ2) small enough, we can assume that ΩM¯ϵωi¯= for all ϵ(ϵM,ϵM). In view of definition 4.3, it is natural to set

    UM[ϵ,γ1,γ2,γ3](x)ΩoSn(xy)Mo[ϵ,γ1,γ2,γ3](y)dσy+ωiSn(xϵs)Mi[ϵ,γ1,γ2,γ3](s)dσsxΩM¯,
    for all (ϵ,γ1,γ2,γ3)(ϵM,ϵM)×U. By proposition 4.1 and real analyticity results for integral operators with real analytic kernel (cf. [51]), we verify that UM is a real analytic map from (ϵM,ϵM)×U to C1,α(ΩM¯) and that equality (4.14) holds. By proposition 4.1, we also deduce that UM[0,d0,η0,r0]=u~M|ΩM¯+Sn|ΩM¯Ωogodσ and, by standard properties of the single-layer potential (cf. [31], §4.4), that u~M is a solution of problem (4.15). The proof is complete.

    Remark 4.5.

    By proposition 4.1, remark 4.2 and theorem 4.4, if F is linear and

    ξ~=1ωidσ(Ωogodσd0ωiωiSn(ts)μ~i(s)dσsdσtr0ωigidσ)0,
    we deduce that the value of the solution at a fixed point x¯Ωo¯{0} is asymptotic to
    1/(ωidσ)(Ωogodσd0ωiωiSn(ts)μ~i(s)dσsdσtr0ωigidσ)δ(ϵ)ϵn1asϵ0.
    If we further assume that
    Ωo=ωi=Bn(0,1),go(x)=axBn(0,1)andgi(t)=btBn(0,1)
    for some constants a,bR, then if
    abr0+an2d00,
    we deduce that the value of the solution at a fixed point x¯Bn(0,1)¯{0} is asymptotic to
    abr0+ad0/(n2)δ(ϵ)ϵn1asϵ0.
    Thus, we recover the result of §2 on the toy problem.

    We now consider in theorem 4.6 the behaviour of the rescaled solution u(ϵ,ϵt).

    Theorem 4.6.

    Let the assumptions of proposition 4.1 hold. Let Ωm be a bounded open subset of Rnωi¯. Then there exist ϵm(0,ϵ2) and a real analytic map Um from (ϵm,ϵm)×U to C1,α(Ωm¯) such that ϵΩm¯Ω(ϵ)¯ for all ϵ(0,ϵm) and

    u(ϵ,ϵt)=1ϵn2Um[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](t)+Ξ[ϵ,ϵδ(ϵ),η(ϵ),(ϵn1/ρ(ϵ))]δ(ϵ)ϵn1xΩm¯
    for all ϵ(0,ϵm). Moreover, if we set
    u~m(t)ωiSn(ts)μ~i(s)dσstRnωi,
    we have that Um[0,d0,η0,r0]=u~m|Ωm¯ and u~m solves the (nonlinear) Robin problem
    {Δu(t)=0tRnωi¯,νωiu(t)=F~(d0u(t)+ξ~,η0)+gi(t)r0tωi,limtu(t)=0.4.16

    Proof.

    Taking ϵm(0,ϵ2) small enough, we can assume that ϵΩm¯Ωo¯ for all ϵ(ϵm,ϵm). By definition 4.3, we note that if ϵ(0,ϵm), then

    u(ϵ,ϵt)=ΩoSn(ϵty)Mo[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](y)dσy+ωiSn(ϵtϵs)Mi[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](s)dσs+Ξ[ϵ,ϵδ(ϵ),η(ϵ),(ϵn1/ρ(ϵ))]δ(ϵ)ϵn1=1ϵn2(ϵn2ΩoSn(ϵty)Mo[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](y)dσy+ωiSn(ts)Mi[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](s)dσs)+Ξ[ϵ,ϵδ(ϵ),η(ϵ),(ϵn1/ρ(ϵ))]δ(ϵ)ϵn1tΩm¯.
    Accordingly, we set
    Um[ϵ,γ1,γ2,γ3](t)ϵn2ΩoSn(ϵty)Mo[ϵ,γ1,γ2,γ3](y)dσy+ωiSn(ts)Mi[ϵ,γ1,γ2,γ3](s)dσstΩm¯,
    for all (ϵ,γ1,γ2,γ3)(ϵm,ϵm)×U. By proposition 4.1 and real analyticity results for integral operators with real analytic kernel (cf. [51]), we verify that Um is a real analytic map from (ϵm,ϵm)×U to C1,α(Ωm¯) and that equality (4.14) holds. By proposition 4.1, we also deduce that Um[0,d0,η0,r0]=u~m|Ωm¯ and, by standard properties of the single-layer potential (cf. [31], §4.4), that u~m is a solution of the (nonlinear) Robin problem (4.16). The proof is complete.

    Finally, we consider the energy integral Ω(ϵ)|u(ϵ,x)|2dx when ϵ is close to 0.

    Theorem 4.7.

    Let the assumptions of proposition 4.1 hold. Let u~m be as in theorem 4.6. Then there exist ϵe(0,ϵ2) and a real analytic map E from (ϵe,ϵe)×U to R such that

    Ω(ϵ)|u(ϵ,x)|2dx=1ϵn2E[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)]4.17
    for all ϵ(0,ϵe). Moreover,
    E[0,d0,η0,r0]=Rnωi|u~m(t)|2dt.4.18

    Proof.

    Let ϵ(0,ϵ1). By the divergence theorem, we have that

    Ω(ϵ)|u(ϵ,x)|2dx=Ωou(ϵ,x)νΩou(ϵ,x)dσxϵωiu(ϵ,x)νϵωiu(ϵ,x)dσx=Ωou(ϵ,x)νΩou(ϵ,x)dσxϵn1ωiu(ϵ,ϵt)νωi(t)u(ϵ,ϵt)dσt.
    Then, let UM and ϵM be as in theorem 4.4, with ΩMΩoBn(0,rM)¯ for some rM>0 such that Bn(0,rM)¯Ωo. Then one verifies that if ϵ(0,ϵM),
    Ωou(ϵ,x)νΩou(ϵ,x)dσx=ΩoUM[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](x)νΩo(x)UM[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](x)dσx.
    Then, let Um and ϵm be as in theorem 4.6, with ΩmBn(0,rm)ωi¯ for some rm>0 such that Bn(0,rm)ωi¯. Then one verifies that if ϵ(0,ϵm),
    ϵn1ωiu(ϵ,ϵt)νωi(t)u(ϵ,ϵt)dσt=1ϵn2ωiUm[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](t)νωi(t)Um[ϵ,ϵδ(ϵ),η(ϵ),ϵn1ρ(ϵ)](t)dσt.
    As a consequence, we set ϵemin{ϵm,ϵM} and
    E[ϵ,γ1,γ2,γ3]ϵn2ΩoUM[ϵ,γ1,γ2,γ3](x)νΩo(x)UM[ϵ,γ1,γ2,γ3](x)dσxωiUm[ϵ,γ1,γ2,γ3](t)νωi(t)Um[ϵ,γ1,γ2,γ3](t)dσt
    for all (ϵ,γ1,γ2,γ3)(ϵe,ϵe)×U. We verify that E is a real analytic map from (ϵe,ϵe)×U to R and that equality (4.17) holds. Moreover, by the behaviour at infinity of u~m and the divergence theorem on exterior domains (cf. [31], §§3.4 and 4.2), we verify that
    E[0,d0,η0,r0]=ωiu~m(t)νωi(t)u~m(t)dσt=Rnωi|u~m(t)|2dt,
    and accordingly equality (4.18) holds.

    5. Conclusion

    We have studied the asymptotic behaviour of solutions of a boundary value problem for the Laplace equation in a perforated domain of Rn, n3, with a (nonlinear) Robin boundary condition which may degenerate into a Neumann condition on the boundary of a small hole of size ϵ. Under suitable assumptions, for ϵ close to 0, the value of the solution at a fixed point far from the origin behaves as 1/(δ(ϵ)ϵn1), where δ(ϵ) is a coefficient of a nonlinear function of the trace of the solution in the Robin boundary condition. We have also investigated the behaviour of the energy integral of the solutions as ϵ tends to 0: the energy integral behaves as 1/ϵn2 multiplied by the energy integral of a solution of an exterior nonlinear Robin problem. In particular, if δ(ϵ)=ϵr, then to satisfy assumption (4.3) we need to have r1, and we have that the value of the solution at a fixed point behaves as 1/ϵn1+r, whereas the energy integral behaves as 1/ϵn2 (and such behaviour is not affected by the specific power δ(ϵ)=ϵr). As we have seen, our study is confined to the case of dimension n3. We plan to investigate the two-dimensional case (which requires a different analysis due to the logarithmic behaviour of the fundamental solution) in a forthcoming paper. Moreover, together with the study of the planar case, we wish to include numerical examples.

    Data accessibility

    This article has no additional data.

    Authors' contributions

    P.M.: conceptualization, writing—original draft, writing—review and editing, final version of the manuscript, proving of theorems; G.M.: conceptualization, writing—original draft, writing—review and editing, final version of the manuscript, formulation of the problem and discussion of the modelling problem.

    All authors gave final approval for publication and agreed to be held accountable for the work performed therein.

    Conflict of interest declaration

    We declare we have no competing interests.

    Funding

    The authors acknowledge support from the EU through the H2020-MSCA-RISE-2020 project EffectFact, grant agreement ID 101008140. P.M. also acknowledges the support of the grant ‘Challenges in Asymptotic and Shape Analysis - CASA’ from the Ca’ Foscari University of Venice. G.M. also acknowledges support from Ser Cymru Future Generation Industrial Fellowship no. AU224 – 80761.

    Acknowledgements

    The authors thank the referees for several valuable comments. Part of the work was done while P.M. was visiting Martin Dutko at Rockfield Software Limited. P.M. wishes to thank Martin Dutko for useful discussions and the kind hospitality. P.M. is a member of the Gruppo Nazionale per l’Analisi Matematica, la Probabilità e le loro Applicazioni (GNAMPA) of the Istituto Nazionale di Alta Matematica (INdAM). G.M. thanks the Royal Society for the Wolfson Research Merit Award.

    Footnotes

    One contribution of 14 to a theme issue ‘Non-smooth variational problems and applications’.

    Published by the Royal Society under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/4.0/, which permits unrestricted use, provided the original author and source are credited.

    References