Parametric Gevrey asymptotics in two complex time variables through truncated Laplace transforms

This work is devoted to the study of a family of linear initial value problems of partial differential equations in the complex domain, dealing with two complex time variables. The use of a truncated Laplace-like transformation in the construction of the analytic solution allows one to overcome a small divisor phenomenon arising from the geometry of the problem and represents an alternative approach to the one proposed in a recent work (Lastra and Malek in Adv. Differ. Equ. 2020:20, 2020) by the last two authors. The result leans on the application of a fixed point argument and the classical Ramis–Sibuya theorem.


Introduction
This work is devoted to the study of a family of singularly perturbed partial differential equations in the complex domain of the form Q(∂ z )u(t 1 , t 2 , z, ) = P t k 1 +1 1 ∂ t 1 , t 2 , ∂ t 2 , ∂ z , z, + f (t 1 , t 2 , z, ), (1) under initial data u(0, t 2 , z, ) ≡ u(t 1 , 0, z, ) ≡ 0, with Q(X) ∈ C[X] and P(T 1 , T 21 , T 22 , Z, z, ) being a polynomial in (T 1 , T 21 , T 22 , Z) with holomorphic coefficients w.r.t. (z, ) on H β × D(0, 0 ). Here, H β and D(0, 0 ) stand for the horizontal strip {z ∈ C : | Im(z)| < β} and the disc at the origin and radius 0 , for some β > 0, 0 > 0, respectively. The forcing term f (t 1 , t 2 , z, ) is holomorphic on C × D(0, h ) × H β × E, for any open sector E centered at 0 and contained in D(0, 0 ), for some h > 0, and remains close to a polynomial in t 1 , analytic in t 2 on D(0, h ) and in z on H β , as becomes close to the origin in C. The variable acts as a small complex parameter. The concrete assumptions on the elements involved in the main problem (1) are to be described and analyzed in detail throughout the work. The study of a problem of such form is motivated by the recent research [12] of the second and third authors. The main aim in the preceding work was related to the descrip-tion of the asymptotic behavior of the analytic solutions, with respect to the perturbation parameter, near the origin, of singularly perturbed equations Q(∂ z )u(t 1 , t 2 , z, ) =P t k 1 +1 1 ∂ t 1 , t k 2 +1 2 ∂ t 2 , ∂ z , z, +f (t 1 , t 2 , z, ), (2) withP(T 1 , T 2 , Z, ) being a polynomial in (T 1 , T 2 , Z, z, ) with holomorphic coefficients w.r.t. (z, ) on H β × D(0, 0 ). Two main novelties are considered here with respect to it. On the one hand, the irregular singular operators related to the second time variable stay rigid in (2), as a polynomial function of the operator t k 2 +1 2 ∂ t 2 . In the present study, the irregular operators in this variable fit a more general scheme within the problem, under certain technical assumptions (see (5) and (6)). This, at first sight slight, variation on the form of the main problem varies its underlying geometry radically. On the other hand, the appearance of different types of solutions observed in [12], known as inner and outer solutions, which describe boundary layer expansions do not appear in the present situation, since we study local solutions in time t 1 , t 2 near the origin in the complex domain. It is worth mentioning that, despite the fact that the form of the main equation under study resembles that of [12], the nature of the singularities appearing in the problem require one to appeal different approaches and apply novel techniques, to be briefly described below.
This work continues a line of research on the study of the asymptotic behavior of solutions of singularly perturbed PDEs in the complex domain, under the action of two time variables: dealing with a symmetric factorized (resp. asymmetric) leading term [11] (resp. [8]), the mentioned work [12], and the corresponding q-analog [10] in the framework of q-difference-differential equations.
In order to illustrate the origin of the family of equations under study in this work, we provide a simple equation of the form (1). More precisely, the change of variable x 1 = 1/t 1 , x 2 = 1/t 2 applied to the equation 1 + 2 t 2 1 ∂ t 1 t 2 2 ∂ t 2 u(t 1 , t 2 , z, ) = f (t 1 , t 2 , z, ) determines the second order PDE where v(1/t 1 , 1/t 2 , z, ) := u(t 1 , t 2 , z, ). Moreover, if x 1 , x 2 , are assumed to be real variables, with = 0, the previous PDE is of hyperbolic type. These kinds of equations are of great importance in mathematical modeling of wave propagation such as sound or electromagnetic waves. The technique developed in the present work consists on searching for solutions of the main problem (see (8) for its precise expression) in the form of a Fourier, truncated Laplace and Laplace transform of certain function, for every fixed value of the perturbation parameter : The integration path L 1, stands for the segment [0, h 1 ( )e √ -1θ 1 ], for some holomorphic function → h 1 ( ) on the domain of definition of the perturbation parameter, approaching to infinity when tends to 0, and some θ 1 ∈ R which does not depend on . The integration with respect to the path L 2 stands for a usual Laplace transform along certain half line [0, ∞)e √ -1d 2 , for some d 2 ∈ R, whereas the function ω belongs to a certain Banach space which depends on the choice of the perturbation parameter. In previous work [11] and [8], we have searched for solutions of the form of a double Laplace transform, and Fourier integral.
This construction through a truncated Laplace transform is proposed in order that the solutions (3) remain close, as tends to 0, to a double usual Laplace transform in both variables t 1 , t 2 . For such a complete double Laplace representable solution, a direct analysis of the asymptotic behavior w.r.t. is unfortunately not possible (as shown in our previous work [12]). However, such study turns achievable within the new approach regarding truncated Laplace transform solutions.
Given a finite family of finite sectors E = (E p ) 0≤p≤ι-1 which conform a good covering (see Definition 3), the first main result in the work, Theorem 1, states the existence of a solution of the main problem in the form (3) for every 0 ≤ p ≤ ι -1, remaining holomorphic in a domain T 1 × T 2 × H β × E p , where T 1 , T 2 are finite sectors with vertex at the origin. Moreover, the exponential decrease of the difference of two solutions associated to consecutive sectors in E enables the application of the classical Ramis-Sibuya theorem (RS) in order to achieve the second main result of our study, namely the asymptotic relation of the analytic solutions and the formal solution of the main problem in powers of , with coefficients in some complex Banach space (see Theorem 3).
In recent years, several steps have been taken to contribute to the knowledge of the asymptotic behavior of analytic solutions of singularly perturbed partial differential equations in the complex domain. We first refer to the recent work [21,22], by Yamazawa and Yoshino, and Yoshino resp., in which the parametric Borel summability of semilinear systems of PDEs is studied, first in the case of fuchsian operators, and second combining both irregular and fuchsian operators. We refer to [2,16] as introductory texts on the classical theory of summability of formal solutions of differential equations in the complex domain.
The appearance of truncated Laplace transform is closely related to the classical theory of asymptotic approximation of analytic functions (examples of this situation is the classical proof of Ritt's theorem for Gevrey asymptotics; see [2] Proposition 10, and also Lemma 1.3.2 in [16]). Truncated Laplace transform also appears as a recent object of study in the literature, related to differential operators [14,15], but also from the numerical point of view [13]. We also refer to [20], where the authors study fractional impulsive differential equations. In [1], variable-order fractional difference equations are studied by means of results on Laplace transform. Other references are [18] on related results in the thermodynamic theory of porous elastic bodies, and [3] on a numerical study on the effects of heat transfer and Hall current on the sinusoidal motion of solid particles through a planar channel. The choice of an integration path for Laplace transform which depends on each fixed value of the perturbation parameter has been inspired from [5,17].
Throughout the work, we use bold letters to indicate a vector of two variables: we write τ for the pair (τ 1 , τ 2 ), u for (u 1 , u 2 ), T for (T 1 , T 2 ), etc.
The paper is organized as follows.
In Sect. 2.1, we recall some properties on Fourier transform which allow one to transform the main problem, stated in Sect. 2.2, in the form a convolution problem, described in Sect. 2.3. The geometry of the problem is an important matter in this work, which needs to be explained in detail. Section 3 is focused on this issue. The Banach spaces involved in the construction of the analytic solution of the auxiliary problem, and some of their main properties, are stated in Sect. 4. Such function is constructed in Sect. 5. The analytic solution of the main problem is obtained in Sect. 6 (Theorem 1), and the work concludes with the description of the parametric Gevrey asymptotic expansions of the analytic solution, obtained in Sect. 7 (Theorem 3).

Layout of the main and auxiliary problems
In this initial section, we describe in detail the main problem under study (8) (Sect. 2.2), and the conditions on the elements involved in it. The solution of this problem is reduced to a convolution auxiliary problem in the Borel plane (17) when inspecting solutions in the particular form of a triple Fourier, Laplace and truncated Laplace transform (see Sect. 2.3). We first give some words about inverse Fourier transform on certain Banach spaces which act on the transformation of the problem (Sect. 2.1).

Inverse Fourier transform on certain function spaces
The transformation of the main problem with respect to variable z requires recalling some basic facts about inverse Fourier transform when acting on certain Banach spaces of real functions of exponential decrease at infinity. The pair (E β,μ , · (β,μ) ) is a Banach space.
The next result will be needed in our reasoning. We refer to [6], Proposition 7, for its proof.
Remark In Sect. 3 we assume further geometric conditions on these polynomials. In particular, observe that condition (18) implies that deg(Q) ≥ deg(R D 1 D 2 ).
We choose μ ∈ R with μ > max The main aim in this work is to study the following initial value problem: for the initial conditions u(t 1 , 0, z, ) ≡ u(0, t 2 , z, ) ≡ 0. Let us describe the form of the elements involved in the problem. Let 0 > 0 and β > 0. For all 1 ≤ 1 ≤ D 1 -1 and 1 ≤ 2 ≤ D 2 -1, the term c 1 2 (z, ) are holomorphic functions on H β × D(0, 0 ). We recall that H β stands for the horizontal strip (4). The function c 1 2 is defined by where m → C 1 2 (m, ) is continuous for m ∈ R and is subject to uniform exponentially flat upper bounds with respect to ∈ D(0, 0 ), i.e. there exists C c > 0 such that Observe that m → C 1 2 (m, ) belongs to E (β,μ) with sup ∈D(0, 0 ) for any given open sector E centered at 0, and contained in D(0, 0 ) \ {0}, for some positive number h . The forcing term is constructed as follows. Let N 1 ≥ 0 and F n 1 ,n 2 (m, ) ∈ E (β,μ) under uniform bounds with respect to in the disc D(0, 0 ). More precisely, assume that , which turns out to be a holomorphic function on C 2 with respect to the first two variables, with coefficients in E (β,μ) . We write where κh 1 ( ) is a holomorphic function on any open sector centered at 0 in the punctured disc D(0, 0 ) \ {0} (see (20)), θ 1 is a real number to be determined and γ (n, z) stands for the incomplete Gamma function which is an entire function w.r.t. z, when n is a fixed positive real number. Observe that the forcing term F depends in particular on the choice of θ 1 .
The following property related to the lower incomplete Gamma function will be crucial in the construction of the auxiliary equation of the problem. Namely, We recall that the infinite Laplace transform satisfies for every positive natural numbers n, k.
This property will be used with respect to the second time variable, whereas a truncated Laplace transform depending on each value of the perturbation parameter near the origin is applied on the first variable. Both (11) and (12) give rise to adequate algebraic properties which allow one to reduce the main equation in the form of an auxiliary problem.
Regarding (10), F is holomorphic w.r.t. T 1 on C , T 2 on the disc D(0, T 0 /2) and on H β w.r.t. z. Furthermore, according to (11), (12), we observe that as tends to 0, for (well chosen) fixed T 1 . Therefore, F is getting closer to a polynomial in T 1 as tends to 0. The function f defined by for any given open sector E centered at 0 and contained in D(0, 0 ) \ {0}, with h > 0 such that 0 < h 0 < T 0 /2. From the remark above, we check in particular that f becomes close to a polynomial in t 1 as becomes closer to the origin.

Auxiliary problems
We search for solutions of (8) in the form of an inverse Fourier transform We make the change of variable T 1 = t 1 and T 2 = t 2 in (8). The classical properties of the inverse Fourier transform (see Proposition 1), together with (5), lead to an auxiliary functional equation satisfied by the expression U(T 1 , T 2 , m, ), namely Let 0 < κ < 1 and θ 1 , d 2 ∈ R. Let → h 1 ( ) be a holomorphic function defined on the domain of definition of the perturbation parameter, to be detailed afterwards. For every fixed value of the perturbation parameter , we search for solutions of (13) written as the Laplace transform with respect to T 2 along direction d 2 and the truncated Laplace transform with respect to T 1 along direction θ 1 -λk 2δD 2 arg( ) applied to a second auxiliary function. More precisely, we search for solutions of (13) of the form where d 1 = d 1 ( ) := arg(h 1 ( )) + θ 1 and L d 1 , stands for the segment [0, consists of the half-line with endpoint at the origin and direction d 2 . The domain of definition of ω and U d 1 d 2 will be discussed in the subsequent sections.
The assumption (6) guarantees the existence of d 2 k 2 ∈ N such that The following result states a one-to-one correspondence between the solution of (13) and (17). Its proof, which is omitted, can be adapted with minor modifications from [8], Lemma 1. (14). Then we have Lemma 1, (15) together with the shape of the solution in (14), the assumptions on the coefficients c 1 2 and the forcing term f , and Lemma 2 entail ω being a solution of the following auxiliary convolution equation in the Borel plane: So far, the solution is of symbolic nature. The geometry of the problem, detailed in the following section, together with Sect. 4 provide convergence and growth estimates of such solution.

On the geometry of the problem
In this section, we preserve the objects and assumptions detailed in Sect. 2.2 on the elements involved in the construction of the main problem under study (8), giving rise to the auxiliary problem (17). This section is devoted to the study of the geometry of the problem, which is crucial in the asymptotic approximation of the solution.
We define for every m ∈ R the polynomial In the case that τ 1 = 0 one can factorize P m in the form We assume that the polynomials Q and R D 1 D 2 satisfy where S Q,R D 1 D 2 stands for an unbounded sector Let λ be a real number which satisfies Let E be a sector with vertex at the origin which is contained in the disc D(0, 0 ) and for every ∈ E we define and the quantities 11 and with r 11 := λk 2δD 2 , and r 22 := λk 1 δ D 1 . The next result summarizes the main properties of the geometric construction above, which will be used to state the asymptotic behavior of the solutions of the main problem.
• Provided that λ > 0 is small enough, for any couple of directions (θ 1 , d 2 ) which satisfy there exist an unbounded sector S d 2 with bisecting direction d 2 and small opening, and a sector S d 1 , , with for all ∈ E. • Let S d 1 , and S d 2 be as above. We put Then, there exists C P > 0 which does not depend on ∈ E such that Proof Let τ = (τ 1 , τ 2 ) ∈ C 2 such that P m (τ ) = 0. One has τ 2 = q (τ 1 , m), for some 0 ≤ ≤ k 2δD 2 -1. In the case that |q (τ 1 , m)| ≤ 2r 2 ( ), from the definition of P m and (18) we derive that The first statements follows from this.
The second statement is a direct consequence of the fact that for all τ 1 ∈ C and m ∈ R one has for every 0 ≤ ≤ k 2δD 2 -1. Regarding the construction of S d 1 , , for all τ 1 ∈ S d 1 , we have The pair (θ 1 , d 2 ) can be chosen accordingly, provided that λ, η Q,R D 1 D 2 ,δ 1 > 0 are small enough.
It only remains to prove that for some constant C P > 0, all m ∈ R and τ ∈ (Ω 1 ( ) × Ω 2 ( )). Usual estimates guarantee that with c 3 = c 1 /c 2 , and some c 4 > 0. Taking into account that we get the existence of c 5 > 0 such that which concludes the proof.

Banach spaces of functions with exponential growth
In this section, we recall the definition and main properties of certain Banach spaces previously used by the authors in [6], and adapted to the several variable case in [8,11]. The dependence of the domains of definition involved in the norm with respect to the values of the perturbation parameter has previously been consider in [7]. Let E be a sector of finite radius in the complex plane. For every ∈ E we consider the following two domains: a finite sector Ω 1 ( ) with vertex at the origin, bisecting direction d 1 which depends on , and radius r 1 ( ); and the union of an infinite sector S d 2 with vertex at the origin, fixed bisecting direction d 2 and positive opening which do not depend on together with the disc D(0, r 2 ( )) for some r 2 ( ) > 0, say Ω 2 ( ), i.e. Ω 2 ( ) = S d 2 ∪ D(0, r 2 ( )).
In the following we write d = (d 1 , d 2 ).
The first result follows directly from the definition of the norm of the Banach space in Definition 2. ∈ E, and let (τ , m) → a(τ , m) be a bounded continuous function on m) is a holomorphic function, then a(τ , m)h(τ , m) belongs to F d (ν,β,μ,k, ) .

Lemma 4 Let
Some parts of the proof of the following result can be adapted from that of Proposition 2 in [6]. We decided to include it completely for the sake of completeness and for making this work self-contained.
Proof Let f ∈ F d (ν,β,μ,k, ) . We have Therefore, one has with The classical estimates for m 1 ≥ 0 and m 2 > 0, yield At this point, we provide upper bounds for C 2.3 ( ). Let g(y) = ay 1+cy b , for some a, c > 0 and b > 1. The function g attains its maximum at y 0 = (c(b -1)) -1/b . We apply this result to the case a = |τ 2 | k 2 , b = σ 1 /σ 1 and c = |τ 2 | k 2 σ 2 to arrive at for some C(σ 1 ,σ 1 ) > 0. We plug the bound (33) into (32) and make the change of variable h = | | k 2 h at the integral in C 2.3 ( ) to arrive at for some C 24 > 0 only depending onσ 1 , σ 1 . Let x 0 > 0. The previous sup particularized for those τ 2 ∈ Ω 2 ( ) such that x := |τ 2 / | k 2 > x 0 reads Let For x > x 0 , one can perform the change of variable h = xu in the integral of (x), in order to get

du.
A partial fraction decomposition allows one to write F k 2 (x) = F 1,k 2 (x) + F 2,k 2 (x), where We observe that for some positive constants F 1,k 2 , F 2,k 2 , valid for all x ≥ x 0 . Under the second assumption in (29), we find that sup x>x 0 (x) is upper bounded by a constant. We conclude that the expression in (35) is upper bounded by for some C 25 > 0. It only rests to provide upper bounds for C 2.3 ( ) regarding the set of τ 2 ∈ Ω 2 ( ) such that 0 ≤ x ≤ x 0 . We observe that We plug this last expression into (32) to arrive at sup (τ 1 ,τ 2 )∈Ω 1 ( )×Ω 2 ( ) for some C 26 , C 27 , C 28 , C 29 > 0. We conclude that the left-hand side in (38) is upper bounded by C 29 | | -r 11 k 1σ1 +k 2 (1+σ 4 ) .
In view of (36) and (39), we derive that which concludes the result.
The proof of the following result can be reproduced under minor adjustments from that of Proposition 2, [9].

Solution of an auxiliary problem
In this section, we preserve the elements and assumptions made on the main problem under study described in Sect. 2.2. More precisely, we assume the conditions (5), (6), (7) are satisfied by the parameters and elements involved. The coefficients c 1 2 (z, ) and the forcing term f (t, z, ) are constructed accordingly. We also assume the geometry of the problem is set in accordance with the assumptions made in Sect. 3 (see condition (18)), and preserve the values of r 1 ( ), r 2 ( ) for each ∈ E and λ > 0 (see (19) and (21)).
We provide a solution of the auxiliary problem (17) by means of a fixed point method in the Banach spaces introduced in Sect. 4.
Let 0 > 0 be such that In view of the choice of 0 > 0 in (50), applied to (49), one has In view of (48) and (51), we find that the operator H , restricted to B(0, ) ⊆ F d (ν,β,μ,k, ) turns out to be a contractive map in the complete metric space for the distance d(x, y) = · (ν,β,μ,k, ) . The classical contractive mapping theorem guarantees the existence of a unique fixed point, say As a result, and regarding the proof of the previous result, one attains the following statement. Proposition 2, the function ω d  k (τ , m, ) is a solution of the auxiliary equation (17). Moreover, for every ∈ E, it satisfies

Corollary 1 Under the assumptions made in
for every τ ∈ Ω 1 ( ) × Ω 2 ( ) and m ∈ R. The constant C ω > 0 can be chosen uniformly for all ∈ E.

Analytic solutions of the main problem
The main aim in this section is to provide analytic solutions of (8) for each of the elements of a family of sectors with respect to the perturbation parameter in the form of a truncated Laplace, Laplace and Fourier transforms. We first fix the geometric elements in this construction.

Definition 3
Let ι be an integer number, ι ≥ 2. Let E := (E p ) 0≤p≤ι-1 , where E p stands for a finite open sector with vertex at the origin, radius smaller than 0 . We assume the intersection of three different elements in E is empty, and 0≤p≤ι-1 E p = U \ {0}, for some neighborhood of the origin U ⊆ C. For the sake of simplicity, we arrange the sectors in order that nonempty intersections of sectors in E correspond to consecutive indices in the ring of integers modulo ι. Under this configuration, we say that E describes a good covering in C .

Definition 4
Let ι be an integer number, ι ≥ 2, and let E := (E p ) 0≤p≤ι-1 be a good covering in C . Let T j be an open sector with vertex at the origin in C and finite radius r T j > 0, for j = 1, 2. For all 0 ≤ p ≤ ι -1 we consider two bounded sectors S d j,p of bisecting direction d j,p , and small opening. In the following statements, we identify the indices p = ι and p = 0. We say that the set is admissible if there exists δ > 0 such that for j = 1, 2 one has for every 0 ≤ p ≤ ι -1, ∈ E p , t j ∈ T j and ξ j ∈ R (which may depend on t j and ) such that Let ι ≥ 2 be an integer number. Let E = (E p ) 0≤p≤ι-1 be a good covering and consider an admissible set {T 1 , T 2 , E, (S d 1,p ) 0≤p≤ι-1 , (S d 2,p ) 0≤p≤ι-1 }, which is associated to the good covering E . We briefly discuss the feasibility of such a construction. Indeed, let 0 ≤ p ≤ ι-1 be fixed. We can first choose the direction ν 1,p (related to a fixed direction θ 1 depending on p) such that (54) holds for j = 1. Then, select the direction ν 2,p = d 2 in order that (54) holds for j = 2 together with the condition stated in the second item of Lemma 3.
Let 0 ≤ p ≤ ι -1. For each 0 ≤ p ≤ ι -1, we consider the main problem under study (8) under the assumptions (5)- (7), and departing from the coefficients c 1 2 (z, ) and the forcing term f (t, z, ) defined in Sect. 2.2. In virtue of the geometry of the problem described in Sect. 3 and Corollary 1, in particular the assumption of condition (18) and the choice of λ and r j ( ) for j = 1, 2 in (19) and (21) resp., one finds that for every ∈ E p there exist a vector of directions d p = (d p,1 ( ), d p,2 ), a bounded sector with vertex at the origin S d p,1 , and bisecting direction d p,1 , with S d p,1 , ⊆ D(0, r 1 ( )) and an infinite sector S d p,2 of bisecting direction d p,2 such that the problem (17)  In view of Corollary 1, one finds that, for every ∈ E p , the function (τ , m) → ω for every τ ∈ Ω p,1 ( ) × Ω p,2 ( ) and m ∈ R. The constant C ω dp k can be uniformly chosen for all ∈ E p . The application of Fourier, Laplace and truncated Laplace transforms to the function ω d p k (τ , m, ) leads to a solution of the main problem under study: for every 0 ≤ p ≤ ι -1 and ∈ E p , we define the function u p (t, z, ) by The integration path L d p,1 , stands for the segment [0, κh 1 ( )e We observe that the choice of the admissible set, compatible with the good covering, together with the bounds in (55) guarantee that (t, z) → u p (t, z, ) is holomorphic on the domain (T 1 ∩ D(0, h )) × (T 2 ∩ D(0, h )) × H β , for 0 < β < β and some h > 0. We recall that H β stands for the horizontal strip Indeed, the construction of ω d p k (τ , m, ) and the definition of u p (t, z, ) in (56) allow one to affirm that the function is holomorphic on the domain ( The properties of Fourier transform (see Sect. 2.1) and Laplace transform (see Lemma 2), together with the definition of the elements involved in the main equation guarantee that (57) represents a solution of the main problem (8).
From now on, we refer to consecutive solutions of (8) to solutions associated to consecutive sectors in the corresponding good covering, which have nonempty intersection.
The next property on the difference of two consecutive solutions will be crucial in order to provide the asymptotic behavior of the solution at 0 regarding the perturbation parameter.
Proof The first part of the proof is guaranteed from the construction of the function u p (t, z, ) for every 0 ≤ p ≤ ι -1.
Let ∈ E p ∩ E p+1 and assume that L d p,1 , can be transformed into L d p+1,1 , by a path deformation and the same holds for L d p,2 and L d p+1,2 without meeting any (τ 1 , τ 2 ) ∈ D(0, r 1 ( )) × (C \ D(0, r 2 ( ))) with P m (τ 1 , τ 2 ) = 0 for m ∈ R, i.e. the movable singularities in D(0, r 1 ( )) × (C \ D(0, r 2 ( ))) fall apart from the arguments between d p,1 and d p+1,1 with respect to the first component, nor between d p,2 and d p+1,2 with respect to the second component. Whenever this configuration holds, Cauchy theorem ensures that D(0, h )) × H β . The same argument can be applied to all ∈ E p ∩ E p+1 concluding that the sectors E p and E p+1 can merge in the configuration of the good covering.
It is worth mentioning that the following cases state three equivalence classes regarding each element in the good covering. A continuity argument shows that, for all 0 ≤ p ≤ ι -1, if there exists ∈ E p ∩ E p+1 such that one of the following mutually excluded cases holds for such , then the same case holds for every element in E p ∩ E p+1 .
Case 1: Assume that L d p,1 , ≡ L d p+1,1 , and L d p,2 differs from L d p+1,2 . This situation occurs in the case that the first component of every singularity in the Borel plane does not fall between the directions d p,1 and d p+1,1 but at least the second component of one singular point in the Borel plane occurs within angles between d p,2 and d p+1,2 .
Case 2: Assume that L d p,2 ≡ L d p+1,2 and L d p,1 , differs from L d p+1,1 , . We only provide details of the steps which differ from the proof of Case 1. We have u p+1 (t, z, )u p (t, z, ) = I 21 -I 22 , where We for 0 < h < (δ 1 /ν 1 ) 1/k 1 , and some C 22 > 0. This entails the existence of C ω dp k ,5 , for all ∈ E p ∩ E p+1 , t ∈ ((T 1 ∩ D(0, h )) × (T 2 ∩ D(0, h ))) and z ∈ H β , with α defined in (59). Case 3: Assume that L d p,2 does not coincide with L d p+1,2 and L d p,1 , differs from L d p+1,1 , . For a more compact writing, we will only display the integration paths in which the integrals involved are subdivided. Each of them can be reduced to the situation in case 1 or case 2 above. In the following steps, we preserve the notation for , p and p+1 , and consider In the previous expression, we have extended in a natural manner the notation adopted for the integration paths in Case 1 and Case 2. Analogous bounds as those stated for the integral I 25 (resp. I 23 ) are also valid for I 33 (resp. I 34 ), in Case 2. For the expression I 35 (resp. I 36 ) one can consider the estimates used to study I 14 (resp. I 13 ), involved in Case 1. We conclude to the existence of C ω k ,6 , C 24 > 0 such that u p+1 (t, z, )u p (t, z, ) ≤ C ω k ,6 exp -C 24 | | α , for all ∈ E p ∩ E p+1 , t ∈ ((T 1 ∩ D(0, h )) × (T 2 ∩ D(0, h ))) and z ∈ H β , with α defined in (59). Figure 1 illustrates the deformation of the paths involved in the procedure.

Parametric Gevrey asymptotic expansions of the analytic solutions
In this section, we analyze the asymptotic behavior of the analytic solutions of the main problem (8) obtained in the previous section, regarding the perturbation parameter approaching the origin. The classical criterion for k-summability of formal power series with coefficients in a Banach space, known as the Ramis-Sibuya theorem (see [2], p.121, or Lemma XI-2-6 in [4]) will be used to describe the Gevrey asymptotic approximation of the solution. The assumptions made in Sect. 2.2 and construction of the elements related to the main problem under study (8) are maintained in this section.
We first give some words on this classical summability theory for the sake of completeness.

k-summable formal power series and Ramis-Sibuya theorem
Let (E, · E ) be a complex Banach space. for every integer N ≥ 1 and ∈ V 1 . Watson's lemma guarantees the uniqueness of such a function, known as the k-sum of the formal power series.
Theorem 2 (RS) Let ι ≥ 2 and let (E p ) 0≤p≤ι-1 be a good covering in C . For every 0 ≤ p ≤ ι -1 we consider a holomorphic function G p : E p → E, and define the function Θ p ( ) := G p+1 ( ) -G p ( ) holomorphic in Z p := E p ∩ E p+1 . We assume the following statements hold: • G p is a bounded function for ∈ Z p , → 0 for all 0 ≤ p ≤ ι -1.
• Θ p is an exponentially flat function of order k in Z p for all 0 ≤ p ≤ ι -1, i.e. there exist K, M > 0 such that Θ p ( ) E ≤ K exp(-M | | k ), valid for all ∈ Z p , and each 0 ≤ p ≤ ι -1. Then, each of the functions G p ( ), for 0 ≤ p ≤ ι -1 admits a common formal power serieŝ G( ) ∈ E[[ ]] as Gevrey asymptotic expansion of order 1/k on E p . In addition to this, if the opening of E p 0 is larger than π/k for some 0 ≤ p 0 ≤ ι -1, then G p 0 ( ) is unique, being the k-sum ofĜ( ) on E p 0 .