Theoretical and numerical local null controllability of a quasi-linear parabolic equation in dimensions 2 and 3

Theoretical and numerical local null controllability of a quasi-linear parabolic equation in dimensions 2and 3 E. Fern´andez-Caraa, J. L´ımacob, I. Mar´ın-Gaytec aDpto. EDAN & IMUS, Universidad de Sevilla, Aptdo. 1160, 41080 Sevilla, Spain bInst. Matem´atica - Universidade Federal Fluminense, RJ, Brasil cDpto. EDAN, Universidad de Sevilla, Aptdo. 1160, 41080 Sevilla, Spain Abstract This paper is devoted to the theoretical and numerical analysis of the null controllability of a quasi-linear parabolic equation. First, we establish a local controllability result. The proof relies on an appropriate inverse function argument. Then, we formulate an iterative algorithm for the computation of the null control and we prove a convergence result. Finally, we illustrate the analysis with some numerical experiments. AMS Classifications: 35A35, 35K59, 49K20, 49M15, 65M60. Keywords and phrases: Quasi-linear parabolic PDEs, null controllability, iterative algorithms for control problems. 1. Introduction and preliminaries In this work, we analyze the theoretical and numerical null controllability in spatial dimensions 2 and 3 of the following PDE system        yt−∇·(a(y)∇y) = v˜ 1ω,(x, t)∈Q:= Ω ×(0, T), y= 0,(x, t)∈Σ := ∂Ω×(0, T), y(x, 0) = y0(x), x ∈Ω. (1) Here, Ω ⊂RNis a bounded connected open set (N≤3), ω⊂⊂ Ω is a nonempty open set (the control domain), ˜ 1ω∈C∞ 0(Ω) satisfies 0 <˜ 1ω≤1 in ωand ˜ 1ω= 0 outside ωand a∈C3(R) possesses bounded derivatives of order ≤3 and satisfies 0< m ≤a(r)≤M∀r∈R. It will be said that (1) is (globally) null-controllable at time Tif, for any y0∈H1 0(Ω), there exists a control function v∈L2(ω×(0, T)) and an associated state satisfying y(x, T) = 0 in Ω.(2) Email addresses: [email protected] (E. Fern´andez-Cara), [email protected] (J. L´ımaco), [email protected] ( I. Mar´ın-Gayte) On the other hand, it will be said that (1) is locally null-controllable at time Tif there exists ε > 0 such that, for any y0∈H1 0(Ω) with ky0kH1 0≤ε, there exists a control function v∈L2(ω×(0, T)) and an associated state satisfying (2). Recently, important progress has been made in the controllability analysis of linear and semilinear PDEs. We refer to the works [7, 8, 12, 15] and the references therein. For the controllability of equations with nonlocal terms, see [6]. For systems of the form (1) in spatial dimension 1, see [11]. Consequently, it is natural to try to extend the known results to equations like (1). In this paper, the first main result is the following: Theorem 1.1. Under the previous assumptions on the coefficient a, there exists ε > 0such that, if y0∈H3(Ω) ∩H1 0(Ω) and ky0kH3≤ε, the nonlinear system (1) possesses a solution satisfying (2). Note that, in this result, the initial data is assumed to belong and be small in H3(Ω) ∩H1 0(Ω). This is more restrictive than the condition required for the definition of a locally null-controllable system. As shown below, this is needed in the proof. In order to prove Theorem 1.1, we will employ a technique relying on the so called Liusternik’s Inverse Function Theorem, see [1]. Since the nonlinear term appears in the main part of the partial derivative operator, several nontrivial difficulties appear (among other things, we will have to work with regular and small initial data). Thus, in a first step, we consider the linearized system at zero        yt−a(0)∆y=v˜ 1ω+h(x, t),(x, t)∈Q, y= 0,(x, t)∈Σ, y(x, 0) = y0(x), x ∈Ω. (3) It is well known that, under some appropriate assumptions on h, (3) is null-controllable. More precisely, the adjoint of (3) is given by        −ϕt−a(0)∆ϕ=F(x, t),(x, t)∈Q, ϕ= 0,(x, t)∈Σ, ϕ(x, T) = ϕT(x), x ∈Ω, (4) where ϕT∈L2(Ω); the announced null controllability property is implied by a well known Carleman inequality that can be established for any solution to a system of the form (4). In a second step, we rewrite the null controllability problem for (1) as an equation in a well chosen space of “admissible” state-control pairs: H(y, v) = (0, y0),(y, v)∈Y. (5) 2 Then, we apply Liusternik’s Theorem and we deduce the (local) desired result. To this purpose, we previously have to establish some nontrivial estimates for the null controls and the associated states of (3). This paper is also devoted to the computation of a null control for (1). This is not a simple task; see [5, 14, 3, 4] for some achievements concerning the numerical controllability of linear and nonlinear PDEs. Here, we will argue as in [6, 11], taking advantage of the surjectivity of H0(0,0). Thus, let Ybe the Hilbert space where we can find a solution (y, v) to (5) (see (30a) in Section 3 for the precise definition of Y). We introduce the following iterative algorithm: ALG 1: 1. Choose (y0, v0)∈Y. 2. Then, for given n≥0 and (yn, vn)∈Y, compute (yn+1, vn+1)=(yn, vn)− H0(0,0)−1(H(yn, vn)−(0, y0)).(6) In these iterates, we use H0(0,0)−1, which is by definition an inverse to the left of H0(0,0); the precise definitions of Hand H0(0,0) will be given in Section 4. Note that ALG 1 is an elementary quasi-Newton method and consequently has the following property: the finite dimensional approximations of the iterates lead to a set of algebraic systems whose coefficient matrices are always the same. This is very interesting and convenient from the numerical viewpoint, since it allows to perform just one factorization at the beginning and then compute quickly every (yn+1, vn+1). In our second main result, we prove the convergence of ALG 1 and we furnish some estimates: Theorem 1.2. Let y0∈H3(Ω)∩H1 0(Ω) be given with ky0kH3≤ε(εis furnished by Theorem 1.1). There exists κ∈(0,1) such that, if (y0, v0)∈Yand k(y0, v0)−(y, v)kY≤κ, then the (yn, vn)converge to (y, v)and satisfy k(yn+1, vn+1)−(y, v)kY≤θk(yn, vn)−(y, v)kY(7) for all n≥0for some θ∈(0,1). Remark 1.1. A natural question is whether Theorems 1.1 and 1.2 also hold for similar systems with PDEs of the form yt−∇·(a(x, t;y)∇y) = v˜ 1ω, 3 that is, with a nonlinear diffusion coefficient nonhomogeneous in space and time. The answer is yes, provided a:Q×R7→ Ris regular enough. More precisely, if a=a(x, t;r) satisfies (for instance) 0< m ≤a(x, t;r)≤M∀(x, t, r)∈Q×R and ∃∂n+ma ∂xn i∂rm,∂n+ma ∂tn∂rm∈C0 b(Q×R) for 0 ≤m≤3,1≤i≤Nand n= 0,1, the arguments in Sections 2–4 can be easily adapted to this situation. 2 The paper is organized as follows. Section 2 is devoted to recall some known results concerning the null controllability of the linearized system (3). In Section 3, we prove the first main result in this paper (Theorem 1.1). Finally, Section 4 focuses on the proof of the second main result (Theorem 1.2), the numerical computation of the (yn, vn) and the presentation of the results of some numerical experiments. In the sequel, we will denote by Ca generic positive constant, usually depending on Ω, ω,T and sometimes the function a. Also, k·k and (·,·) will stand for the usual norm and scalar product in L2(Ω), respectively. Finally, we will use ∂iwto denote the partial derivative of wrespect to xi. 2. Carleman inequalities and the null controllability of (3) Let ω0⊂⊂ ωbe a non-empty open set. The following technical result, due to Fursikov and Imanuvilov [12], is fundamental. Lemma 2.1. There exists a function α0∈C2(Ω) satisfying:    α0(x)>0∀x∈Ω, α0(x)=0 ∀x∈∂Ω, |∇α0(x)|>0∀x∈Ω\ω0. (8) Let m=m(t) be a function satisfying m∈C∞([0, T]), m ≥T2 8in [0, T/2], m(t) = t(T−t) in [T/2, T]; here we can consider, for example, m(t) = T2 4−T2 8e−1 (T/2−t)2+4 T2for t∈[0, T/2]. which satisfies the aforementioned properties. Let us set ψ(x, t) := eλα0(x) t(T−t), β(x, t) := β(x) t(T−t):= eRλ −eλα0(x) t(T−t), where α0is as in Lemma 2.1, R > kα0kL∞+ log(2), λ > 0 and φ(x, t) := eλα0(x) m(t), α(x, t) := α(x) m(t):= eRλ −eλα0(x) m(t). The following result from [12] is well known: 4 Lemma 2.2. There exist positive constants λ0,s0and C0only depending on Ω,ωand Tsuch that, for any s≥s0and λ≥λ0, any F∈L2(Q)and any ϕT∈L2(Ω),the associated solution to (4) satisfies the Carleman estimate ZZQ e−2sβ (sψ)−1|ϕt|2+|∆ϕ|2+λ2(sψ)|∇ϕ|2+λ4(sψ)3|ϕ|2dx dt ≤C0 ZZQ e−2sβ |F|2dx dt +ZZω0×(0,T ) e−2sβλ4(sψ)3|ϕ|2dx dt!. (9) From now on, we fix s=s0and λ=λ0. The following result holds: Proposition 2.1. There exist a positive constant Conly depending on Ω,ω,sand Tsuch that, for any F∈L2(Q)and any ϕT∈L2(Ω), the associated solution to (4) satisfies the Carleman estimate ZZQ e−2sα φ−1|ϕt|2+|∆ϕ|2+φ|∇ϕ|2+φ3|ϕ|2dx dt ≤C ZZQ e−2sα |F|2dx dt +ZZω0×(0,T ) e−2sαφ3|ϕ|2dx dt!, (10) Proof: We want to get an estimate of the following integral: ZZQ e−2sα φ−1|ϕt|2+|∆ϕ|2+φ|∇ϕ|2+φ3|ϕ|2dx dt. First, we see that ZZQ e−2sα φ−1|ϕt|2+|∆ϕ|2+φ|∇ϕ|2+φ3|ϕ|2dx dt =ZZΩ×(0,T/2) e−2sα φ−1|ϕt|2+|∆ϕ|2+φ|∇ϕ|2+φ3|ϕ|2dx dt +ZZΩ×(T/2,T ) e−2sα φ−1|ϕt|2+|∆ϕ|2+φ|∇ϕ|2+φ3|ϕ|2dx dt. (11) The last integral is bounded by the right hand side of (10) (for some appropriate C) thanks to Lemma 2.2. Indeed, we have ZZΩ×(T/2,T ) e−2sα φ−1|ϕt|2+|∆ϕ|2+φ|∇ϕ|2+φ3|ϕ|2dx dt =ZZΩ×(T/2,T ) e−2sβ ψ−1|ϕt|2+|∆ϕ|2+ψ|∇ϕ|2+ψ3|ϕ|2dx dt ≤C ZZQ e−2sβ |F|2dx dt +ZZω0×(0,T ) e−2sβψ3|ϕ|2dx dt!. In view of the definitions of αand β, one has e−2sβ =e−2sα for t∈(T/2, T); on the other hand, for t∈(0, T/2), the following holds: e−2sβ ≤1≤K0e−2sα for K0:= e2sK1and K1:= e8Rλ/T 2. In a similar way, we also have e−2sβψ3≤Ce−2sαφ3in Q: one has e−2sβψ3=e−2sαφ3 for t∈(T/2, T); also, there exist constants K2and K3(only depending on Ω, ωand T) such that e−2sβψ3≤K2≤K3e−2sαφ3for t∈(0, T/2). 5 As a consequence, we can bound the last two integrals above respectively by CZQ e−2sα |F|2and CZω0×(0,T ) e−2sαφ3|ϕ|2 and our assertion follows. On the other hand, in order to bound the first integral in the right hand side of (11), we can use energy estimates. Indeed, by multiplying (4) by ϕand integrating in space and then in time in [t, s], with t∈ [0, T/2] and s∈[T/2,3T/4], we can see that 1 2kϕ(·, t)k2+a(0) ZZΩ×(t,s) |∇ϕ|2dx dσ =ZZΩ×(t,s) F·ϕ dx dσ +1 2kϕ(·, s)k2 ≤CZZΩ×(0,3T/4) |F|2dx dσ +a(0) 2ZZΩ×(t,s) |∇ϕ|2dx dσ +1 2kϕ(·, s)k2 and kϕ(·, t)k2+a(0) ZZΩ×(0,T/2) |∇ϕ|2dx dσ ≤CZZΩ×(0,3T/4) |F|2dx dσ +kϕ(·, s)k2. Now, integrating with respect to sin [T/2,3T/4] , one has T 4kϕ(·, t)k2+a(0)T 4ZZΩ×(0,T/2) |∇ϕ|2dx dσ ≤CZZΩ×(0,3T/4) |F|2dx dσ +ZZΩ×(T/2,3T/4) |ϕ|2dx dσ. (12) On the other hand, by multiplying (4) by −∆ϕand integrating again in space and time, the following holds: 1 2k∇ϕ(·, t)k2+a(0) ZZΩ×(t,s) |∆ϕ|2dx dσ =ZZΩ×(t,s) F(−∆ϕ)dx dσ +1 2k∇ϕ(·, s)k2, whence k∇ϕ(·, t)k2+a(0) ZZΩ×(0,T/2) |∆ϕ|2dx dσ ≤CZZΩ×(0,3T/4) |F|2dx dσ +k∇ϕ(·, s)k2 and, as before, integration with respect to sin [T/2,3T/4] yields T 4k∇ϕ(·, t)k2+a(0)T 4ZZΩ×(0,T/2) |∆ϕ|2dx dσ ≤CRRΩ×(0,3T/4) |F|2dx dσ +RRΩ×(T/2,3T/4) |∇ϕ|2dx dσ. (13) Finally, we note from (4) that |ϕt|2≤1 2|∆ϕ|2+1 2|F|2.(14) From (12)- (14), we obtain: ZZΩ×(0,T/2) |ϕt|2+|∆ϕ|2+|∇ϕ|2+|ϕ|2dx dt ≤C ZZΩ×(0,3T/4) |F|2dx dt +ZZΩ×(T/2,3T/4) |∇ϕ|2+|ϕ|2dx dt!. 6 Taking into account that αand φare uniformly bounded from above and from below in Ω × (0,3T/4), we can deduce from this inequality the required estimate. Indeed, let us first note that αis bounded from above in Ω×(0,3T/4) by a constant K4(again depending only on Ω, ωand T). Consequently, ZZΩ×(0,3T/4) |F|2dx dt ≤K5ZZQ e−2sα|F|2dx dt, where K5:= e2sK4. On the other hand, we also have e2sβψ−1≤K6and e2sβψ−3≤K7in Ω × (T/2,3T/4) for similar positive constants. Recalling Lemma 2.2 and arguing as above, we get: ZZΩ×(T/2,3T/4) |∇ϕ|2+|ϕ|2dx dt ≤K6ZZQ e−2sβψ|∇ϕ|2dx dt +K7ZZQ e−2sβψ3|ϕ|2dx dt ≤C ZZQ e−2sβ |F|2dx dt +ZZω0×(0,T ) e−2sβψ3|ϕ|2dx dt! ≤C ZZQ e−2sα |F|2dx dt +ZZω0×(0,T ) e−2sαφ3|ϕ|2dx dt! Hence, the first integral in the right hand side of (11) is also bounded by the right hand side of (10). This ends the proof. 2 Let us introduce the weights ρ:= esα, ρ3:= esαφ−3/2, ρ5:= esαm5/2, ρ7:= esαm7/2, ρ9:= esαm9/2, ρ11 := esαm11/2(15) and the constants α1:= min x∈Ω α(x) = eRλ −eλkα0kL∞, α2:= max x∈Ω α(x) = eRλ −1. It is then clear that 2α1> α2and the following inequalities hold: esα1/m(t)< esα(x,t)< esα2/m(t)< e2sα1/m(t). The main result in this section is the following. It relies on the null controllability of (3), that is of course well known. The estimates in this result are also known in the literature and have been used with various purposes, although we believe that is is useful to recall them here, as well as their proofs. Theorem 2.1. Assume that the function hin (3) satisfies ρ3h∈L2(Q),ρ9ht∈L2(Q)and h(·,0) ∈H1 0(Ω). Then (3) is null-controllable. In fact, for each y0∈H1 0(Ω), there exist null controls vwith ρ7v∈L2(0, T;H2(ω)) ∩C0([0, T]; H1(ω)),(ρ7v)t∈L2(ω×(0, T)), v(·,0) ∈H1(ω) (16) 7 and associated states ysatisfying (2). Furthermore, if y0∈H3(Ω) ∩H1 0(Ω), the associated states are such that N1(y):=ZZQ [ρ2|y|2+ρ2 5|∇y|2+ρ2 7(|yt|2+|∆y|2)+ρ2 9|∇yt|2+ρ2 11(|ytt|2+|∆yt|2)] dx dt<+∞, N2(y):= sup 0≤t≤TZΩ [ρ5|y|2+ρ2 7|∇y|2+ρ2 9(|yt|2+|∆y|2)+ρ2 11|∇yt|2]dx<+∞. (17) Proof: Let us introduce some notation: •Ly := yt−a(0)∆y, L∗ϕ:= −ϕt−a(0)∆ϕ, •P0={ϕ∈C2(Q) : ϕ= 0 on Σ}, π(ϕ, ˜ϕ) := ZZQρ−2L∗ϕL∗˜ϕ+ 1Ωρ−2 3ϕ˜ϕdx dt, •P= the completion of P0for the scalar product π(·,·). •b(ϕ) := ZΩ y0(x)ϕ(x, 0) dx +ZZQ h(x, t)ϕ dx dt. Then, bis continuous linear form on the Hilbert space P(thanks to the Carleman inequality (10)). From the results in [12], we know that, for any y0∈L2(Ω) and any hwith ZZQ ρ2 3|h|2<+∞, there exist controls v∈L2(ω×(0, T )) and associated solutions to (3) satisfying (2). The couples (y, v) can be found by minimizing ZZQ ρ2|y|2dx dt +ZZω×(0,T ) ρ2 3|v|2dx dt in the family of the state-control pairs (y, v) with v∈L2(ω×(0, T)). Furthermore, we have ZZQ ρ2|y|2dx dt +ZZω×(0,T ) ρ2 3|v|2dx dt ≤Cky0k2+ZZQ ρ2 3|h|2dx dt and ZZQ e−2sαhφ−1(|ϕt|2+|∆ϕ|2) + φ|∇ϕ|2idx dt ≤Cky0k2+ZZQ ρ2 3|h|2dx dt. Accordingly, we can write y=ρ−2L∗ϕ, v =−ρ−2 3ϕ|ω×(0,T ),(18) where ϕis the unique solution to the linear problem π(ϕ, ˜ϕ) = b( ˜ϕ)∀˜ϕ∈P, ϕ ∈P. Let us set z:= −ρ−2 3ϕ. Then ρ7z=ρm7/2z=−ρ−1m1/2e3λα0ϕand, after some computations, we see that L∗(ρ7z) = −ρ−1m1/2e3λα0L∗ϕ −∂t(ρ−1m1/2e3λα0)ϕ+ 2∇(ρ−1m1/2e3λα0)· ∇ϕ+ ∆(ρ−1m1/2e3λα0)ϕ =J1+J2+ 2J3+J4. (19) 8 In view of (18), one has: |J1| ≤ ρ−1m1/2e3λα0ρ2|y| ≤ Cρ|y|.(20a) Also, in view of the facts that m−a≤Cm−bif a<band m0is bounded and the definition of ρ, we have that |J2|+|J4| ≤ Cρ−1m−3/2|ϕ|(20b) and |J3| ≤ Cρ−1m−1/2|∇ϕ|.(20c) From the Carleman estimate (10) written for ϕand the fact that y=ρ−2L∗ϕ, we get: ZZQ ρ−2λ2(sφ)|∇ϕ|2dx dt +ZZQ ρ−2λ4(sφ)3|ϕ|2dx dt ≤ZZQ ρ2|y|2dx dt +ZZω×(0,T ) ρ2 3|v|2dx dt < +∞, whence ZZQ ρ−2m−1|∇ϕ|2dx dt +ZZQ ρ−2m−3|ϕ|2dx dt < +∞.(20d) From (20a)-(20d), we deduce that J1+J2+ 2J3+J4∈L2(Q). Consequently, taking into account the PDE satisfied by ρ7zand the fact that (ρ7z)(·, T) = 0, we see that ρ7z∈L2(0, T;H2(Ω)) ∩C0([0, T]; H1 0(Ω)),(ρ7z)t∈L2(Q).(21) In particular, (16) holds. Notice that, up to now, we have only used that ρ3h∈L2(Q) and y0∈L2(Ω). In order to get (17), we will establish several estimates: Estimates I: Multiplying (3) by ρ2 5yand integrating in Ω, we have that 1 2 d dt ZΩ ρ2 5|y|2dx +a(0) ZΩ ρ2 5|∇y|2dx ≤CZω ρ2 3|v|2dx +ZΩ ρ2 3|h|2dx +ZΩ ρ2|y|2dx+ZΩ ρ2 5|y|2dx, due to |ρ1,t| ≤ Cρ2,|∇ · (ρ5∇ρ5)| ≤ Cρ2and ρ2 5≤Cρ2 3. In view of Gronwall Lemma, we find that sup 0≤t≤TZΩ ρ2 5|y|2dx +ZZQ ρ2 5|∇y|2dx dt ≤C ZZω×(0,T ) ρ2 3|v|2dx dt +ZZQ ρ2 3|h|2dx dt +ky0k2!. (22) Now, multiplying (3) by ρ2 7ytand integrating in Ω, we see that ZΩ ρ2 7|yt|2dx +a(0) d dt ZΩ ρ2 7|∇y|2dx ≤CZω ρ2 3|v|2dx +ZΩ ρ2 3|h|2dx +ZΩ ρ2 5|∇y|2dx 9 The argument is the following. For any (h, y0)∈Za solution to (3) in Ycan be obtained by solving the following extremal problem:      Minimize ZZQ ρ2|y|2dx dt +ZZω×(0,T ) ρ2 3|v|2dx dt Subject to v∈L2(ω×(0, T)),(y, v) satisfies (3). (41) It is known that (41) possesses exactly one solution, given by y=ρ−2L∗p, v =−ρ2 0p|ω×(0,T ),(42) where pis the unique solution to the Lax-Milgram problem    π(p, p0) = ZZQ hp0dx dt +ZΩ y0(x)p0(x, 0) dx ∀p0∈P, p ∈P (43) (recall the notation introduced in the proof of Theorem 2.1). Accordingly, we set H0(0,0)−1(h, y0) = (y, v), with yand vrespectively given by (42) and (43). Obviously, (43) is the weak formulation of the following boundary-value problem, that is secondorder in time and fourth-order in space:        L(ρ−2L∗p) + ρ−2 3p1ω=h(x, t) (x, t)∈Q, p= 0, ρ−2L∗p= 0 (x, t)∈Σ (ρ−2L∗p)|t=0 =y0(x),(ρ−2L∗p)|t=T= 0 x∈Ω. Let us recall the algorithm proposed in Section 1: ALG 1: 1. Choose (y0, v0)∈Y. 2. Then, for given n≥0 and (yn, vn)∈Y, compute (yn+1, vn+1)=(yn, vn)− H0(0,0)−1(H(yn, vn)−(0, y0)).(44) For each n, the task reduces to the solution of a problem of the kind (43) with h=−∇ · ((a(0) −a(yn))∇yn)).(45) The convergence of ALG 1 is established in Theorem 1.2. Let us give the proof. Proof of Theorem 1.2: First, note that for any initial (y0, v0)∈Y, the iterates in ALG 1 are well defined. We will use a standard argument, appropriate for methods of this kind, that rely on the C1regularity of H; see for instance [2]. Thus, let us assume that ky0kH3≤ε(εis furnished by Theorem 1.1) and (y, v)∈Ysatisfies H(y, v) = (0, y0); let us set CH:= kH0(0,0)−1kL(Z;Y)and let us assume that 0 < θ < 1/(2CH). Since His continuously differentiable, there exists δ > 0 such that (˜y, ˜v)∈Y, k(˜y, ˜v)kY≤δ⇒ kH0(˜y, ˜v)− H0(y, v)kL(Y;Z)≤θ. 16 We will assume that k(y, v)kY≤δand we will prove that there exists κ > 0 such that, if (y0, v0)∈Y and k(y0, v0)−(y, v)kY≤κ, (46) then the (yn, vn) satisfy (7). Let κbe such that, if (˜y, ˜v)∈Y, k(˜y, ˜v)−(y, v)kY≤κ, then kH(˜y, ˜v)− H(y, v)− H0(y, v)((˜y, ˜v)−(y, v))kZ≤θk(˜y, ˜v)−(y, v)kY and let us choose (y0, v0)∈Ysatisfying (46). Then, if we introduce en:= (yn, vn)−(y, v), the following holds: en+1 =en− H0(0,0)−1(H(yn, vn)− H(y, v)) − H0(0,0)−1(H(yn, vn)− H(y, v)− H0(y, v)en) − H0(0,0)−1(H0(y, v)− H0(0,0))en. Therefore, ken+1kY≤CHkH(yn, vn)− H(y, v)− H0(y, v)enkZ +CHkH0(y, v)− H0(0,0)kL(Y;Z)kenkY. Since ke0kY≤κ, this inequality for n= 0 yields ke1kY≤2CHθke0kY and, in particular, we also have ke1kY≤κ. By induction, we then see that kenkY≤2CHθken−1kY≤ · · · ≤ (2CHθ)nke0kY for all n≥1. This proves that en→0 in Yand (7) holds with θ= 2CHθ.2 In order to solve numerically the problems (43), it suffices in principle to construct explicity finite dimensional spaces Ph⊂P. Note however that this is possible but needs a considerable work. This is because the functions in Pmust satisfy L∗p∈L2 loc(Q). Thus, an approximation based on a standard triangulation of Qrequires spaces Phof functions that must be globally C0 in all the variables and globally C1in space. This construction can be complex and too expensive and it is convenient to use instead a mixed formulation, as in [10, 11]. Let us briefly indicate how this can be done. Let us introduce the new variables z=ρ−1L∗p and m=ρ−1 3p, the spaces Z:= {(z, m)∈L2(Q)×L2(Q):(ρ3m)t∈L2(Q),∇(ρ3m)∈L2(Q)N} and Λ := {λ:ρλ ∈L2(Q),∇(ρλ)∈L2(Q)}, the bilinear forms α((z, m),(z0, m0)) := ZZQ zz0dx dt +ZZω×(0,T ) mm0dx dt, β((z, m), λ) := ZZQhλz+ρ−1((ρ3m)t)− ∇(ρ−1λ)· ∇(ρ3m)idx dt 17 and the linear form h`, (z, m)i:= ZZQ ρhm dx dt +ZΩ ρ3(x, 0)y0(x)m(x, 0) dx. It is not difficult to check that α(·,·), β(·,·) and `are well defined and continuous, respectively, on Z×Z, Z ×Λ and Z. Then, an appropriate mixed formulation of (43) is the following: Find (z, m)∈Zand λ∈Λsuch that (α((z, m),(z0, m0)) + β((z0, m0), λ) = h`, (z0, m0)i ∀(z0, m0)∈Z, β((z, m), λ0) = 0 ∀λ0∈Λ.(47) What we have to do is, consequently, to solve numerically (47) and then take y=ρ−1z, v =−ρ−1 3m|ω×(0,T ). Contrarily to P, it is not difficult to construct finite dimensional subspaces of Zand Λ. This leads to “natural” mixed approximations and allows efficient and computationally reasonable computations. Thus, let Thbe a triangulation of the cylinder Q= Ω ×(0, T), and let us set Zh:= {(zh, mh)∈C0(Q)×C0(Q) : zh|K, mh|K∈P1(K)∀K∈ Th, zh|Σ= 0, mh|Σ= 0, zh(·, t) = mh(·, t) = 0 for t∈[T−h, T ]} and Λh:= {λh∈C0(Q) : λh|K∈P1(K)∀K∈ Th, λh(·, t) = 0 for t∈[T−h, T]}. Then, we can approximate (47) as follows: Find (zh, mh)∈Zhand λh∈Λhsuch that (α((zh, mh),(z0 h, m0 h)) + β((z0 h, m0 h), λh) = h`, (z0 h, m0 h)i ∀(z0 h, m0 h)∈Zh, β((zh, mh), λ0 h) = 0 ∀λ0 h∈Λh.(48) In the following sections, we present the results of some experiments. 4.1. A first test (Test 1) The quasi-Newton method has been applied to the solution to the null controllability problem for (1) with the following data: •N= 2, Ω = (0,1) ×(0,1), ω = (0.2,0.8) ×(0.2,0.8), T = 0.5. •y0(x1, x2) = sin(πx1) sin(πx2). •a(s) = exp(−2 exp(−0.3s)). 18 Figure 1: The mesh. Number of vertices: 7425. Number of tetrahedrons: 38976. Note that, this choice of ahas a sense. The function s7→ exp(−2 exp(−0.3s)) can be viewed as the density profile of an isolated population under some particular circumstances. Thus, we can motivate our tests by the control of the spread of a desease in a related habitat, with diffusion depending on the number of individuals. At each step of ALG 1, the problem has been rewritten in the form (47). Then, the finite element approximation (48) has been introduced. The computations have been performed with the FreeFem++ package; for a detailed description, see http://www.freefem.org//ff++ . In this and the other tests, the stopping criterion has been kyn+1 −ynkL2 kyn+1kL2 ≤κ where yn=ρ−1zn,znis (together with some mnand λn) the solution to (48) and κ= 10−5. The mesh is displayed in Fig. 1. Starting from (y0, v0) = H0(0,0)−1(0, y0),convergence was reached after 11 iterates, see Fig. 2. The initial state can be viewed in Fig. 3. The computed control and the associated state are depicted in Figs. 4–7. 19 Figure 2: Evolution of the error at logarithmic scale for Test 1. Figure 3: The initial state y0= sin(πx1)·sin(πx2). 20 (a) Computed control (b) Associated state Figure 4: The computed control and the associated state of Test 1 at t= 0.1. Newton method. (a) Computed control (b) Associated state Figure 5: The computed control and the associated state of Test 1 at t= 0.25. 21 (a) Computed control (b) Associated state Figure 6: The computed control and the associated state of Test 1 at t= 0.4. (a) Computed control (b) Associated state Figure 7: The computed control and the associated state of Test 1 at x1= 0.68. The L2norms of the computed control and the state as functions of tare displayed in Fig. 8. 22 (a) Computed control (b) Associated state Figure 8: Evolution in time of the L2norms of the control and the state for Test 1. Figure 9: Evolution of the error at logarithmic scale for Test 2. 4.2. A second test (Test 2) In a second experiment, we have taken the same Ω, ωand Tand we have fixed •y0(x, y) = sin(πx1) sin(2πx2). •a(s) = a0(1 + 5 ·sin(50s)),with a0=e−2. This time, starting from (y0, v0) = H0(0,0)−1(0, y0),convergence was reached after 33 iterates (see Fig. 9). The initial state and the computed control and associated state are depicted in Figs. 10–13. 23 (a) Computed control (b) Associated state Figure 10: The computed control and the associated state of Test 2 at t= 0.1. (a) Computed control (b) Associated state Figure 11: The computed control and the associated state of Test 2 at t= 0.25. 24 (a) Computed control (b) Associated state Figure 12: The computed control and the associated state of Test 2 at t= 0.4. (a) Computed control (b) Associated state Figure 13: The computed control and the associated state of Test 2 at x1= 0.68. Finally, the evolution in time of the L2norms of the control and the state is depicted in Fig. 14. 25