Carmen Camacho

Agustín Pérez-Barahona

Land use dynamics and the environment

Keywords: Land use, Spatial dynamics, Pollution. Journal of Economic Literature: Q5, C6, R1, R14

This paper builds a benchmark framework to study optimal land use, encompassing land use activities and environmental degradation. We focus on the spatial externalities of land use as drivers of spatial patterns: land is immobile by nature, but local actions a↵ect the whole space since pollution flows across locations resulting in both local and global damages. We prove that the decision maker problem has a solution, and characterize the corresponding social optimum trajectories by means of the Pontryagin conditions. We also show that the existence and uniqueness of time-invariant solutions are not in general guaranteed. Finally, a global dynamic algorithm is proposed in order to illustrate the spatial-dynamic richness of the model. We find that our simple set-up already reproduces a great variety of spatial patterns related to the interaction between land use activities and the environment. In particular, abatement technology turns out to play a central role as pollution stabilizer, allowing the economy to reach a time-invariant equilibrium that can be spatially heterogeneous.

Introduction

Land use activities are usually defined as the transformation of natural landscapes for human use or the change of management practices on human-dominated lands [START_REF] Foley | Global Consequences of Land Use[END_REF]. It is widely accepted that these activities have greatly transformed the planet's surface, encompassing the existence and evolution of spatial patterns (for instance, [START_REF] Plantinga | The E↵ect of Agricultural Policies on Land Use and Environmental Quality[END_REF][START_REF] Kalnay | Impact of Urbanization and land-use change on climate[END_REF][START_REF] Chakir | Predicting land use allocation in france: A spatial panel data analysis[END_REF]. In this regard, Spatial Economics analyses the allocation of resources over space as well as the location of economic activity and, thus, the formation of spatial patterns. In particular, great e↵ort has been devoted to understanding firms' location, transport costs, trade, and regional and urban development [START_REF] Duranton | Urban evolutions: The fast, the slow, and the still[END_REF]. However, the mechanisms behind the interaction between land use and the environment that can induce spatial patterns, designated in our paper as spatial drivers, are still far for being understood.

In this paper we contribute to the theoretical foundations of land use change and the environment by considering the interaction between land use activities and pollution.

To this end we will develop a theoretical model that focuses on the spatial externalities of land use as drivers of spatial patterns.

There is an abundant literature on the interaction between land use and pollution.

Agricultural research in particular has devoted great attention to the e↵ects of pollution on agricultural land use (for instance, [START_REF] Hamilton | The Benefits of Pollution control: The Case of Ozone and U.S. Agriculture[END_REF][START_REF] Deschênes | The Economic Impacts of Climate Change: Evidence from Agricultural Output and Random Fluctuations in Weather[END_REF]. About the environmental influence of land use, many papers have identified significant environmental impacts of land use (among others, [START_REF] Matson | Agricultural intensification and Ecosystem Properties[END_REF]and Kalnay and Cail, 2003). Moreover, [START_REF] Foley | Global Consequences of Land Use[END_REF] point out that the e↵ects of environmental degradation due to land use are global but also regional/local.

Although this literature has been very fruitful, the dominant approach has been empirical. There is indeed a general agreement about the lack of explicit modelling of the spatial drivers behind the interaction between land use and pollution. Closely related to the integrated assessment approach, bottom-up models of agricultural economics (for instance, [START_REF] De Cara | Methane and nitrous oxide emissions from agriculture in the EU: A spatial assessment of sources and abatement costs[END_REF] have contributed to the understanding of the spatial drivers of land use. However, these models focus on partial equilibrium (mainly the supply side) and do not completely consider the intertemporal dimension of the problem. In this paper we use an alternative approach based on the Dynamic Spatial Theory (see [START_REF] Desmet | On spatial dynamics[END_REF], for a survey).

Within this theory, considering the forward-looking dimension of agents' decisions, the natural spatial generalization of the Ramsey model is presented in [START_REF] Brito | The dynamics of growth and distribution in a spatially heterogenous world[END_REF] and [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF]2013a). Both include a policy maker who decides the trajectory for consumption at each location. The main feature of these models is the spatial dynamics of capital, which flows in space to meet optimal decisions according to a partial di↵erential equation (PDE). Although these sophisticated models are very promising, several technical problems have been identified (Boucekkine et al., 2013b).

In particular, the application of parabolic PDEs in this new field has opened a set of questions still not solved by the mathematical literature. To date, there have been few pragmatic approaches that provide alternative set-ups. For instance, [START_REF] Costello | Optimal harvesting of stochastic spatial resources[END_REF] provide a dynamic framework to study the optimal harvesting of renewable resources in a stochastic spatial (partial equilibrium) model. Taking advantage of the special structure of the problem, they are able to analytically characterize the equilibrium. [START_REF] Desmet | Spatial Growth and Industry Age[END_REF], 2010and 2013), more in line with the spatial Ramsey model, follow the idea of imposing enough structure to the spatial problem (through factors' mobility, di↵usion of technology, and land and firm ownership) as well. Agents are assumed to be myopic. While each location solves a static problem, their model is dynamic in time. Indeed, each location decides the optimal amount to consume, how much to invest in R&D, and how much to save, taking land revenues, prices and salaries as given. Finally, all savings are coordinated by a cooperative that invests along the space. Even if this approach allows us to understand some important geographic features, the structure of their framework makes the planner's problem intractable (see also [START_REF] Desmet | On the Spatial Economic Impact of Global Warming[END_REF]. Another interesting alternative is the one followed by Brock and[START_REF] Brock | Di↵usion-induced instability and pattern formation in infinite horizon recursive optimal control[END_REF]2010). Considering [START_REF] Derzko | Necessary and su cient conditions for optimal control of quasilinear partial di↵erential systems[END_REF], they approximate (linear quadratic) the original nonlinear optimal control problem, around a time-invariant equilibrium. However, as we will show later, neither existence nor uniqueness of time-invariant solutions are ensured in an environmental spatial Ramsey framework.

We use in this paper the spatial generalization of the Ramsey model in order to understand the spatial drivers behind land use and the environment. To the best of our knowledge, our paper provides a first analytically tractable general equilibrium framework of land use that, without approximating the original optimal control problem, encompasses (i) spatial and time dimensions which are presented in a continuous manner, (ii) spatial externalities due to pollution and abatement activities, and (iii) the social optimum. Our starting point is the Spatial Ramsey model in [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF]2013a). We propose a benchmark framework in continuous time and space to study optimal land use. Each location is endowed with a fixed amount of land, which is allocated among production, pollution abatement, and housing. Although the unique production input (land) is spatially immobile by nature, this is a model of spatial growth where local actions a↵ect the entire space through pollution. Indeed, we assume that the production generates local pollution, which flows across locations. In this regard, we illustrate the di↵usion mechanism by means of the well-known Gaussian Plume equation (see Sutton, 1947a,b). Finally, we consider that local pollution damages production due to its negative e↵ect on land productivity. Moreover, we assume that pollution as a whole (global pollution) may also reduce production. This indirect consequence of pollution can be linked, for instance, to the negative e↵ect of anthropogenic GHGs on climate change. 1We prove the existence of a social optimum when the planning horizon is finite. The policy maker decisions are characterized by the Pontryagin conditions. We additionally extend our analytical results to the time-invariant equilibrium. As observed above, this particular equilibrium is crucial to apply solution methods based on approximations of the original problem around a time-invariant equilibrium. We show in this respect that the existence and uniqueness of time-invariant solution are not guaranteed in general.

Finally, to illustrate the richness of our model, we undertake numerical simulations.

To this end we adapt the methodology first developed in [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF] to the current problem. Our algorithm is an alternative framework to other numerical tools that focus on the local dynamics around a time-invariant solution. This numerical analysis is actually global, where we simulate the entire trajectory of the states, controls, and co-states from their initial distributions until they eventually reach, or not, a timeinvariant equilibrium. With the numerical tool in hand, we study the di↵erent drivers of spatial heterogeneity. We find, among other things, that the abatement technology stands out as a fundamental element to achieve time-invariant solutions, which are compatible with the emergence of long-run spatial patterns. Moreover, even if our paper focuses on land use dynamics, many simulated scenarios are consistent with the predictions of spatial models of natural resources such as the harvesting stochastic spatial approach of [START_REF] Costello | Optimal harvesting of stochastic spatial resources[END_REF].

The paper is organized as follows. Section 2 describes the economic model. Section 3 is devoted to the analytical results of our paper. In Section 4 we introduce the algorithm that is applied in the numerical exercises of Section 5. Finally, Section 6 concludes.

The model

We assume that there exists a continuum of locations along a unidimensional region R ✓ R. We also consider that R is an open and connected real set. 2 Each location has a unit of land, which can be devoted to three di↵erent activities: production, housing and pollution abatement. 3 For simplicity, we shall assume that the space required for housing at each location is equal to its population density f (x). We also consider no population growth in this paper. There exists a unique consumption good the production of which only requires land and that we denote by F (l). The remainder of the land is used to abate pollution G(1 l f (x)).

Pollution has two dimensions in our model. The local dimension (local pollution) comes directly from the production of the consumption good. For the sake of simplicity, we assume that each unit of production generates one unit of pollution. It damages production due to the negative e↵ect on land productivity. Moreover, even if land is spatially immobile, local decisions a↵ect the whole space since the pollutant travels across space. We describe the spatial dynamics of pollution by means of a well-known model in physics called the Gaussian plume. It is a standard mathematical description of the dispersion of airborne contaminants (for instance, [START_REF] Arya | Air Pollution Meteorology and Dispersion[END_REF][START_REF] Stockie | The Mathematics of Atmospheric Dispersion Modeling[END_REF].

But it is also used to model the spread of pollutants in aquifers and porous soils and rocks, as well as for nuclear contaminants. According to this model, the dynamics of the pollution at location x in time t, p(x, t), is given by the following second-order partial di↵erential equation (PDE) of parabolic type:

p t (x, t) p xx (x, t) = E(x, t), (1)
where p t and p xx denote, respectively, @p(x, t)/@t and @ 2 p(x, t)/@x 2 , and E(x, t) are the emissions in time t 0 of a single source located at x. The interpretation of equation ( 1) is the following (see [START_REF] Smith | The economics of spatial-dynamic processes: Applications to renewable resources[END_REF], for a detailed description from an environmental economic perspective). The Gaussian plume model comprises two common dispersal mechanisms of pollutants: di↵usion and /or advection. Di↵usion is the spread of a pollutant through regions where its concentration is high to regions

2 Results can be easily extended to the case R ✓ R n , n > 1, and to the case in which R is not connected but the union of connected subsets in R.

3 In this simplified set-up, the land devoted to abatement may be interpreted as being pollution removal due to the presence of, for instance, prairies and forests (see [START_REF] Nowak | Air pollution removal by urban trees and shrubs in the United States[END_REF][START_REF] Ragot | The optimal carbon sequestration in agricultural soils: do the dynamics of the physical process matter?[END_REF]. In general, one can also consider that abatement activities require physical space (land in our model).

of lower concentration (Fick's law), while advection is the flux of contaminants due to wind, ocean currents, etc. As in [START_REF] Brock | Di↵usion-induced instability and pattern formation in infinite horizon recursive optimal control[END_REF]and 2010[START_REF] Smith | The economics of spatial-dynamic processes: Applications to renewable resources[END_REF], we focus on di↵usion. The term p xx (x, t) in (1) reflects indeed the spread due to concentration di↵erential. We pay attention to this dispersal mechanism because our approach is about growth and the long-term response of the economy: the elements behind advection (e.g., wind velocity and direction) are extremely variable, in particular in the short-run, and the time horizon usually considered in this type of problems minimizes this e↵ect. For advection, the other polar case, see for instance [START_REF] Costello | Optimal harvesting of stochastic spatial resources[END_REF] in the context of the spatial economics of natural resources (fish). 4 Additionally, pollution may also harm production as a global pollutant (e.g., anthropogenic GHGs). We then allow for the distinction between local and global pollution, where global pollution is naturally defined as: 5

P (t) = Z R p(x, t)dx.
(2)

We introduce pollution damages in production using a damage function ⌦(p, P, x), where 1 ⌦ represents the share of foregone production due to local and global pollution. 6

If we denote by A(x, t) the total factor productivity at location x at time t, we have that this location produces ⌦(p, P, x)A(x, t)F (l) units of final good when it devotes an amount l of land to production. For simplicity reasons we shall assume that the abatement technology is not a↵ected by pollution. In the remaining of the paper we make the following standard assumptions regarding the production functions:

(H1) Functions F and G are positive, increasing, concave, and their first and second derivatives exist and are non negative, that is:

F (•) 2 C 2 , F(0) = 0, F 0 (•) > 0, F 00 (•)  0, lim s!0 F 0 (s) = 1, lim s!1 F 0 (s) = 0, G(•) 2 C 2 , G(0) = 0, G 0 (•) > 0, G 00 (•)  0, lim s!0 G 0 (s) = 1, lim s!1 G 0 (s) = 0.
4 The Gaussian plume model can be also used in natural resource management. [START_REF] Smith | The economics of spatial-dynamic processes: Applications to renewable resources[END_REF] observe that advection can be eventually modeled "through di↵erences in rates of dispersal", i.e., D(x, t)p xx (x, t), where D(x, t) is the di↵usion coe cient. However, this would require further physical assumptions that are beyond the scope of our paper.

5 Well-known pollutants with mostly global e↵ects are CO 2 and CFCs (see, among others, [START_REF] Nordhaus | Economic Growth and Climate: The Carbon Dioxide Problem[END_REF][START_REF] Akimoto | Global Air Quality and Pollution[END_REF]. However, air contaminants in general (including tropospheric ozone and

NO

x ) are examples of local pollutants that flow among locations.

6 Notice that the productivity loss may also encompass the negative e↵ect of pollution on individuals' health (among others, [START_REF] Pope | Epidemiology of fine particulate air pollution and human health: biologic mechanisms and who's at risk?[END_REF][START_REF] Evans | Do new health conditions support mortality-air pollution e↵ects?[END_REF]. However, we do not explicitly consider this e↵ect in our paper.

(H2) ⌦(p, P, x) is twice di↵erentiable with respect to p and P ; and decreasing in each factor: ⌦ 1 (p, P, x) = @⌦(p,P,x) @p < 0, ⌦ 2 (p, P, x) = @⌦(p,P,x) @P < 0. Function ⌦(p, P, x)

is defined on R + ⇥ R + ⇥ R and takes values in [0, 1].
Assumption (H1) is the usual hypothesis of positive and non-increasing marginal products, together with the Inada conditions. (H2) assumes that both local and global pollution a↵ect negatively production. Moreover, it is also considered that this damage is a smooth function. [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF]2013a) assume that each location produces its own consumption in the social optimum. Social welfare, however, may still increase under the possibility of spatial reallocation of production. We therefore enlarge the set of feasible abatement and production decisions by allowing for consumption "imports".

Indeed, we assume that the policy maker collects all production and re-allocates it across locations at no cost:

Z R c(x, t)f (x)dx = Z R ⌦(p, P, x)A(x, t)F (l)dx, (3)
where c(x, t) denotes consumption per capita at location x and time t. 7

The policy maker chooses consumption per capita and the use of land at each location, which maximize the discounted welfare of the entire population. Following [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF], we introduce two discount functions. The spatial discount represents the weight that the policy maker gives to each location. Alongside their paper, we identify this function as the population density f (x) in order to avoid any subjective spatial preferences. Moreover, as in the standard Ramsey model, we consider the usual temporal discount e ⇢t , with ⇢ > 0. The policy maker maximizes the lifetime discounted utility max {c,l}

Z T 0 Z R u(c(x, t))f (x)e ⇢t dxdt + Z R (p, P )(x, T )e ⇢T dx (4)
subject to

7 In this paper we do not consider transportation costs. However, it is possible to introduce them if proportional to the shipped amount of final good, or under other further assumptions. These could be for instance the compulsory gathering of output in a specific location.

P 8 > > > > > > > > > > > < > > > > > > > > > > > : p t (x, t) p xx (x, t) = ⌦(p, P, x)A(x, t)F (l(x, t)) G(1 l f (x)), R R c(x, t)f (x)dx = R R ⌦(p, P, x)A(x, t)F (l)dx, P (t) = R R p(x, t)dx, p(x, 0) = p 0 (x) 0, lim x! R p x (x, t) = 0, (5)
where (x, t) 2 R ⇥ [0, T ] and denotes R's boundaries. In particular, if

R = R then R = { 1, 1}. Moreover, R = {a, b} if R is an open interval (a, b).
Following the standard approach, we consider an instantaneous utility function u : R + ! R that is increasing and concave. As in [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF], the function is measurable and everywhere finite. It accounts for the planner's concern about the state of pollution at the end of the planning period. In the standard Ramsey model, if the policy maker does not state any function , then the optimal solution is such that savings are zero at T and it is the end of the economy (for instance, [START_REF] Acemoglu | Introduction to Modern Economic Growth[END_REF]. Similarly, if we show no concern about the pollution at the end of the planning period, then pollution will be infinite at T and its shadow price will be zero. Finally, as in [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF], the last expression in (5) is the usual boundary condition: there is no pollution flow in the boundaries of the space. 8

3 The social optimum

In this section we present the theoretical contribution of the paper. We first show that there exist a solution to our problem. Moreover, the optimal trajectories are characterized by the Pontryagin conditions, involving a system of PDEs. Section 3.2 finally focuses on the time-invariant solution, which is defined as the situation when all variables remain constant in time. 9 We prove that both existence and uniqueness of this solution (that can be spatially heterogenous) are not in general guaranteed. In this 8 No pollution flow means that lim x! R p x (x, t) is equal to a constant (this is called the Neumann problem). Without loss of generality, we can assume this constant equal to zero. Notice that, as in Boucekkine et al. (2013a), if space was a circle no boundary condition would be required since the space does not actually have boundaries.

regard, we provide a su cient condition that will be used in the numerical part of the paper.

Optimal trajectories

Let us start by showing, in Proposition 1, that there exist at least a solution to the social optimum problem. In this regard, we prove that P has a unique solution for every choice of the couple (c, l). Notice that this outcome is not a direct application of existing results [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF] because of some special features of P. In particular our model includes a global variable P , defined as the spatial integral of p. Moreover, in contrast to the previous literature, we consider that the policy maker gathers all production to distribute it later, adding the aforementioned supplementary integral constraint on consumption. Consequently, we first have to transform these two integral constraints into partial di↵erential equations in the proof of the proposition.

Afterwards, by imposing the following Assumption (H3), we can apply Theorem 12.1 in Chapter 8 in [START_REF] Pao | Nonlinear Parabolic and Elliptic Equations[END_REF] to close the proof: (H3) For all (x, t) 2 R ⇥ (0, T ], there exist some real constants p 1 > 0, ! > 0, ! 1 > 0, ! 2 > 0 and b < 1/4T , such that, as x ! R,

0 < p(x, t)  p 1 e b|x 2 | , 0 < ⌦(x, t)  !e b|x 2 | , 0 < |⌦ 1 (x, t)|  ! 1 e b|x 2 | , 0 < |⌦ 2 (x, t)|  ! 2 e b|x 2 | .
As in [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF], and [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF], this is a technical assumption that allows us to avoid explosive solutions in the frontiers of the space. Moreover, we should also observe that the exponential terms in (H3) make this hypothesis not very restrictive. For ease of exposition, we report all proofs details of the paper in the Appendices.

Proposition 1. Under assumptions (H1)-(H3), the problem (4)-( 5) has a solution in (x, t) 2 R ⇥ (0, T ], for every T < 1.

Once we know that there exists at least a solution to the social optimum, let us characterize the optimal trajectories. In this regard, we use the Ekeland method of variations in Raymond and[START_REF] Raymond | Pontryagin's principle for state-constrained control problems governed by parabolic equations with unbounded controls[END_REF][START_REF] Pope | Epidemiology of fine particulate air pollution and human health: biologic mechanisms and who's at risk?[END_REF] to obtain the Pontryagin conditions of problem (4)-(5). 10 Following this procedure, we write the associated value function V as a function of c, l, p and P as follows:

V (c, l, p, P ) = R T 0 R R u (c(x, t)) f (x)e t dxdt + R R (p, P )(x, T )e T dx R T 0 R R q(x, t) [p t (x, t) p xx (x, t) ⌦(p, P, x)A(x, t)F (l(x, t)) + G(1 l f (x))] dxdt R T 0 m(t) ⇥ P (t) R R p(x, t)dx ⇤ dt R T 0 n(t) ⇥R R c(x, t)f (x)dx R R ⌦(p, P, x)A(x, t)F (l(x, t))dx ⇤ dt, (6)
where q, m and n are auxiliary functions. We present in Proposition 2 the corresponding Pontryagin conditions, which include the dynamics of the shadow price of pollution q, together with a static equation associated with the optimal land allocation at each (x, t).

Finally, the set of first order conditions also contains a spatial boundary condition on pq

x and a terminal condition on q: Proposition 2. The Pontryagin conditions of problem ( 4)-( 5) are:

8 > > > > > > > > > > > > > > > > > > > > > > > < > > > > > > > > > > > > > > > > > > > > > > > : p t (x, t) p xx (x, t) = ⌦(p, P, x)A(x, t)F (l(x, t)) G(1 l f (x)), q t (x, t) + q xx (x, t) + ⇣ ⌦ 1 (p, P, x) + 1 f (x) ⌦ 2 (p, P, x) ⌘ A(x, t)F (l) [u 0 (c(x, t)) + q(x, t)] + ⇢q = 0, [u 0 (c(x, t)) + q(x, t)] ⌦(p, P, x)A(x, t)F 0 (l) + q(x, t)G 0 (1 l f (x)) = 0, R R c(x, t)f (x)dx = R R ⌦(p, P, x)A(x, t)F (l)dx, P (t) = R R p(x, t)dx, p(x, 0) = p 0 (x) 0, lim x! R p x (x, t) = 0, lim x! R p(x, t)q x (x, t) = 0, q(x, T ) = p (x, T ), R R P (x, T )dx = 0, (7)
for (x, t) 2 R ⇥ [0, T ] and T < 1.

The first condition in ( 7) is the equation of the Gaussian plume, which describes the dynamics of pollution in our set-up. With respect to the standard Ramsey framework, this condition corresponds to the law of motion of the state variable of our problem (pollution in this model), including the additional term p xx (x, t) that represents the di↵usion mechanism of pollution. As in the spatial Ramsey model of [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF]2013b), the second expression is the adjoint equation corresponding to as a deviation from the optimal solution (c ⇤ , l ⇤ , p ⇤ , P ⇤ ) as (c, l, p, P ) = (c ⇤ , l ⇤ , p ⇤ , P ⇤ ) + ✏(, L, ⇡, ⇧), where , L, ⇡, ⇧ are real functions of (x, t) 2 R⇥R + . Then, the optimal solution results from minimizing V with respect to the deviation ✏.

the shadow price (co-state variable) of pollution. Parallel to the Gaussian plume, it is a PDE as well, reflecting that the shadow price varies in time and space because pollution moves in time and space. Moreover, we clearly see in this equation the double dimension of pollution (local and global), which is indeed captured by the marginal e↵ects ⌦ 1 and ⌦ 2 . The third condition represents the trade-o↵ between consumption and pollution. All things being equal, an increase of land devoted to production will raise consumption. However, a greater l will also imply higher marginal social cost due to the increase of pollution (⌦(p, P, x)A(x, t)F 0 (l)) and lower availability of land to abatement (G 0 (1 l f )). We can prove, furthermore, that l(x, t) is uniquely determined by p(x, t) and q(x, t): Proposition 3. l(x, t) is a unique function of p(x, t) and q(x, t).

As in the previous section, the fourth equation represents the re-allocation of production across locations, the expression for P (t) is our definition of global pollution, and p(x, 0) is the given spatial distribution of pollution in time t = 0.

Notice that Proposition 2 also includes two boundary conditions as in the spatial Ramsey model (for further interpretation of the boundary conditions in a spatial dynamic framework, see [START_REF] Smith | The economics of spatial-dynamic processes: Applications to renewable resources[END_REF]and Boucekkine et al., 2013b). The first one is the boundary condition in (5), and the second corresponds to the shadow price of pollution. As in [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF], if we focus on interior solutions, it becomes the standard boundary condition lim

x! R q x (x, t) = 0 implied by the asymptotic constraint on pollution flow in (5).

Moreover, the last two expressions are the terminal conditions of the problem. As in [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF], the first one states that, at the end of the planning period, the shadow price of pollution is equal to the policy maker's marginal concern about the pollution left behind. The second condition says that the spatial aggregate of the marginal concern with respect to P (T ) is zero. In particular, if P does not change sign in R, then this condition amounts to P (x, T ) = 0, for all x 2 R. Note that if our original problem had a dynamic law describing the evolution in time and space of global pollution, then this condition would be similar to the terminal condition linking the final state of local pollution and the marginal concern about it. We provide next a simple example of a function and derive the associated terminal conditions. If (p, P )(x, T ) = p(x, T ) with > 0, then q(x, T ) = and

P (p, P )(x, T ) = 0. Furthermore, Z R (p, P )(x, T )e ⇢T dx = P (T )e ⇢T .
Hence, the policy maker would care about aggregated welfare plus the negative e↵ect of the discounted level of global pollution left after T .

Finally, as a corollary of the Pontryagin conditions (7), we can show that consumption per capita is identical across locations. This is a direct consequence of two main features of our "spatial" policy maker: she re-allocates production across locations, at no cost, and does not have any subjective spatial preferences. Therefore, the instantaneous consumption per capita c(x, t) is spatially homogeneous:

Corollary 1. Consumption per capita is spatially homogeneous, i.e., c(x, t) = c(t).

Notice that in this paper the spatial re-allocation of production (3) does not involve any transportation cost (see also Footnote 7). In fact, the aim of this assumption is to highlight the possible emergence of "specialized" areas. These are defined, in the context of our model, as locations where the majority of their available land is devoted to production or abatement. We study this type of spatial heterogeneity in Section 5. In this regard, the assumption above, together with the homogeneity of residents' preferences, allows us to provide a simple illustration of spatial re-allocation of production, where consumption "imports" are implicitly considered. Even if the number of residents is uniformly distributed in the space, we will see later that the possibility of production re-allocation gives to the social planner the option of specializing some areas for specific activities (abatement or production in this paper), depending on their relative technological advantage.

Let us observe as well that we consider that the residents are homogenous across space. As in [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF], the spatial discount function f (x) "stands for the location's x population density" (p. 24). However, one can also interpret f (x)

as the spatial distribution of individuals' tastes. Assuming one resident per location, this would allow us to consider (spatially) heterogeneous agents, where the individual preferences of a resident of location x are given by U (x, t) ⌘ u(c(x, t))f (x). Following the previous corollary, c(x, t) = c(t) in all locations. This outcome is a direct consequence of the preferences' separability between consumption and the individual taste for it. Nevertheless, we should also notice that residents with greater preference for consumption (i.e., a large f (x)) enjoy a higher level of utility than the individuals of other locations. An interesting line to explore, outside the objectives of the current paper, is to consider heterogeneous agents with non-separable preferences. One could study in this respect how a spatial-dynamic environment would induce and modify an eventual spatially heterogeneous consumption per capita.

The time-invariant solution

We define time-invariant solution as an equilibrium where all variables do not change over time. Therefore, considering the Pontryagin conditions ( 7), let us study the twodimensional system S defined below. We shall actually focus on the solution of the system as a couple (p, q) because, as it is clear from Proposition 3, the third variable at stake l is a unique function of p and q.

If a time-invariant solution (p, q) exists, then it verifies the following system:

S 8 > < > : p xx (x) = ⌦(p, P, x)A(x)F (l(x)) G(1 l f (x)), q xx (x) = ⇣ ⌦ 1 (p, P, x) + 1 f (x) ⌦ 2 (p, P, x) ⌘ A(x)F (l) [u 0 (c) + q(x)] + ⇢q(x),
where P = R R p(x)dx, and l(x) is the unique solution to

[u 0 (c) + q(x)] ⌦(p, P, x)A(x)F 0 (l) + q(x)G 0 (1 l f (x)) = 0, with c = R R ⌦(p, P, x)A(x)F (l)dx/ R R f (x)dx.
Note that abusing of notation, we have kept the same notation for all variables, removing their dependence of time.

We can then prove that the solution to system S is unique in a certain set. In this regard, we provide and apply a less constraining version of Theorem 3.4 in [START_REF] Pao | Nonlinear Parabolic and Elliptic Equations[END_REF].

This result allows us to establish su cient conditions for existence and uniqueness of time-invariant solution:

Proposition 4. Assume space is a bounded interval in R. Given a spatial population distribution f (x), we define a set Z of time-invariant functions

Z = {(p, q) : ⌦ 11 (p, P ), ⌦ 21 (p, P ) > 0 and AF ( l)[⌦(p, P ) p⌦ 1 (p, P )] > G(1 f l)}.
Under the assumptions (H1)-(H3) there exists a unique solution (p, q) to system S in Z.

Together with (H1)-(H3), Proposition 4 establishes further conditions to the damage function ⌦. On the one hand, we have diminishing marginal damages, so that more pollution, holding everything else constant, decreases output by less and less (i.e., ⌦ 11 (p, P ), ⌦ 21 (p, P ) > 0). On the other hand, the remaining production must be large enough, taking into account the abatement technology and the marginal damage of pollution. Let us observe that they are su cient conditions. The proof of the proposition actually allows us to establish alternative conditions for other particular specifications.

The conditions of Proposition 4 are provided for the sake of illustration, bearing in mind the functional forms that we will use in the numerical exercises.

The main message of this result is that the existence and uniqueness of time-invariant equilibrium is not guaranteed in an environmental spatial Ramsey framework. We can identify sets of functions (for instance, Z in Proposition 4) that include the unique time-invariant solution. However, one can not ensure in general that these sets are non-empty. Proposition 4 does not allow either to fully characterize the time-invariant equilibrium. This analytical characterization is very challenging because of the lack of mathematical results involving non-linear PDE systems such as S. But we can make use of the numerical analysis in this respect. Moreover, this analysis also allows us to study the corresponding transition dynamics. This is indeed what we do in Section 5 (together with situations without time-invariant equilibrium), applying the computational method that we present in the next section. From this perspective, we should observe that Proposition 4 is quite helpful: it allows us to conclude, for some cases, if an eventual (simulated) time-invariant equilibrium is the unique time-invariant one in a specific scenario.

Computational setting

Due to the complexity of the Pontryagin conditions (7), we illustrate the richness of our model by means of simulations. Before presenting the details of our method, let us point out that this numerical approach is global. Consequently, the results provided in the subsequent sections are not constrained to economies starting in the neighborhood of any particular equilibrium point. Our simulations, moreover, will also allow us to enrich Section 3.2 by means of studying the convergence to time-invariant solutions. As it is clear from that section, the existence and uniqueness of time-invariant solutions is a demanding mathematical problem. But the convergence of the trajectories to this equilibrium is even more challenging, and still an open question. In this regard, our paper provides a numerical inspection of the convergence. We describe below the computational setting, together with our algorithm to solve (7). This numerical method will be applied in Section 5, where we investigate the emergence and dynamics of spatial patterns in our environmental context.

Let us first rewrite the Pontryagin conditions, reversing time in the equation that describes the dynamic behaviour of q in time and space. Notice that we are allowed to do this operation because the planning horizon is finite. Even if this preliminary action is not necessary, it is convenient for the ease of presentation of the discretization of the Pontryagin conditions and the algorithm. Calling h(x, t) ⌘ q(x, T t), we obtain the following system of parabolic di↵erential equations where we have removed the independent variables (x, t) for simplicity reasons, writing (x, T t) when necessary:

8 > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > < > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > : p t p xx = ⌦(p, P, x)AF (l) G(1 l f ), h t h xx = = h ⌦ 1 (p(x, T t), P (x, T t), x) + 1 f (x) ⌦ 2 (p(x, T t), P (x, T t), x) i ⇥ ⇥AF (l) [u 0 (c(T t)) + h] + ⇢h, [u 0 (c) + h(x, T t)] ⌦(p, P, x)AF 0 (l) + h(x, T t)G 0 (1 l f ) = 0, c(t) = R R ⌦(p,P,x)AF (l)dx R R f (x)dx , P (t) = R R pdx, p(x, 0) = p 0 (x) 0, lim x! R p x (x, t) = 0, lim x! R h x (x, t) = 0, h(x, 0) = p (p(x, T )) , P (p, P )(x, T ) = 0, (8)
for x 2 R = (0, r) and t 2 [0, T ].

The finite di↵erence approximation

The main di culty to simulate the system above is to discretize the two PDEs of the Pontryagin conditions. In this respect, the idea is to implement a finite di↵erence approximation, where we replace the second derivative with respect to space with a central di↵erence quotient in x, and substitute the derivative with respect to time with a forward di↵erence in time. In order to implement this discretization we need to set up a grid in our space (0, r) ⇥ [0, T ]. The points in this grid are couples (j x, n t) for j = 0, 1, ..., J and n = 1, 2, ..., N, where J x = r and N t = T . Then, if v is a function defined on the grid, we write v(j x, n t) = v n j . Let us provide an example. If we want to use a finite di↵erence approximation for the parabolic di↵erential equation @v @t = @ 2 v @x 2 , we write:

11 v n+1 j v n j t = 1 x 2 v n+1 j+1 2v n+1 j + v n+1 j 1 . (9)
We can write (8) as

p n+1 j p n j t 1 x 2 p n+1 j+1 2p n+1 j + p n+1 j 1 = ⌦(p n j , P n j , j)AF (l n j ) G(1 l n j f n j ), ( 10
)
h n+1 j h n j t 1 x 2 h n+1 j+1 2h n+1 j + h n+1 j 1 = = ✓ ⌦ 1 (p T n j , P T n j , j) + 1 f j ⌦ 2 (p T n j , P T n j , j) ◆ AF (l T n j ) ⇥ u 0 (c T n ) + h n j ⇤ + ⇢h n j , (11) ⇥ u 0 (c n ) + h T n j ⇤ ⌦(p n j , P n j , j)AF 0 (l n j ) + h T n j G 0 (1 l n j f n j ) = 0, (12)
with

c n = R J j=0 ⌦(p n j , P n j , j)AF (l n j ) dj R J j=0 f (j)dj . ( 13
)
Abusing of the use of the integral sign, we compute in (13) the integral of a discrete quantity. We treat ⌦(p n j , P n j , j)AF (l n j ) as the J available observations of the continuous variable ⌦(p, P )AF (l). To these equations, we add the border conditions p n J 1 = p n J and h n J 1 = h n J , 8n = 1, 2, ..., N, and the definition of P : P n = P J j=0 p n j .

The algorithm

Our algorithm looks for the solution of the model as the fixed point of an iterative process. We start from an initial guess for the reversed-time shadow price of pollution, {h 0n j } n=1...N j=1...J . Based on this guess and using the discrete time version of the Pontryagin conditions (10) and ( 12), we compute the associated distributions of pollution and land {p n j } n=1...N j=1...J and {l n j } n=1...N j=1...J . Using these resulting distributions in equation ( 10) we can compute the induced distribution for the reversed-time shadow price of pollution {h n j } n=1...N j=1...J . Next, we compute the distance between two iterations of the reverse-time shadow price, that is between {h 0n j } n=1...N j=1...J and {h n j } n=1...N j=1...J .

The optimal solution to (10)-( 12) coincides with the fixed point of this iterative process. Hence, if the solution to two consecutive iterations is close enough, we say that we have reached the fixed point, i.e., the optimal solution, and the algorithm stops. If it is not, then we update the initial guess for h, {h 0n j } n=1...N j=1...J , with the last distribution obtained, {h n j } n=1...N j=1...J , and iterate again until two consecutive iterations become close enough.

To reach our goal we adapt the algorithm developed in [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF] to problem (8). There are still some important di↵erences. First, problem (8) includes a control variable whose dynamics are not described by a PDE. We need to provide a guess for the value of the matrices {h n j } n=1...N j=1...J and {l n j } n=1...N j=1...J , which could compromise the convergence of the algorithm. Indeed, depending on these guesses, we obtain a first approximation to consumption. To increase the convergence speed we run an intermediate loop that improves the initial guess for c and l.

Second, the current problem includes an integral constraint. We opt here for a simple solution. Rather than computing the integral at the current time, n t, we compute it at (n 1) t, that is

P n = R J j=0 p n 1 j dj.
Although this is just an approximation, let us underline that the distance between P (n) and P (n 1) is infinitesimal since P is a continuous function and the distance between points in the grid is su ciently small. In the same manner, using preceeding values for pollution, we compute c n using (13).

As afore mentioned, the convergence of our algorithm crucially depends on the initial guesses {h 0n j } n=1...N j=1...J and {l 0n j } n=1...N j=1...J . To increase the convergence speed, we add an intermediate search step to improve the estimation of c and the initial guess for l. In [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF], the algorithm only required the initial guess for one variable, consumption. Using the initial guess for consumption and the problem's Pontryagin conditions, the algorithm obtained a first estimate for physical capital. In turn, the guess for consumption was actualised using this time the Pontryagin conditions for optimality. The process continued until the distance between two iterations was small enough. In the present problem, we initiate the algorithm with 2 guesses.

Being aware of the severe dependence of the algorithm's convergence on the initial guess, we introduce an intermediate step. At every n, with the available information for global and local pollution, and with the guess for land allocation for moment n, the algorithm finds the associated consumption vector at time n. Then, always taking as given the spatial distribution for local and global pollution at n, the guess for land allocation is recomputed using the new values for consumption. The algorithm iterates until it finds a fixed point, that is, until the distance between two iterations is close enough. In the terms of the algorithm, this step can be described as: initialize the algorithm with a couple {l 0 j n , h 0 j } n=1...N j=1...J . At every n, using {l 0 j n } and the initial values for local and global pollution, we compute a first approximation for c n using (13). Taking {p n 1 j } and {P n 1 } as given, we iterate between equations ( 12) and ( 13) to obtain the fixed point of these equations. This step improves the guess for {l n j }, accelerating the convergence of the next step.

We provide below a synthetic view of the algorithm in its entirety:

Step 1: Initialization

We choose an initial distribution for air pollution p 0 = {p 0,j } and three stopping parameters ✏ i for i = 1, 2, 3. We compute P 0 = P J j=0 p 0 j . We assume an initial guess for the costate variable {h 0n j } n=1...N j=1...J and for land allocation {l 0n j } n=1...N j=1...J .

Step 2: Improvement of the first guess For every n = 1, . . . , N and given p n 1 j , l n 1 j , P n 1 , we compute

c n = R J j=0 ⌦(p n 1 j , P n , j)AF (l n 1 j ) dj R J j=0 f (j)dj .
We repeat the following scheme until the euclidean distance between two consecutive matrices h is smaller than ✏ 1 or until the number of iterations equals a fixed number K.

With c n and the guess {h 0n j } j=1,...,J , we obtain a new guess for {l n j } using (12). We recompute c n with {l n j } instead of {l n 1 j }. We iterate the process until the euclidean distance between two consecutive outcomes for c n is smaller than ✏ 2 .

Then with the resulting c n and {l n j }, we compute p n j for j = 1, ..., J using the upwind algorithm applied to equation (10).

Step 3: Using the values of {l n j }, {p n j } and {c n } in the previous step, we compute a new guess for {h n j } n=1...N j=1...J according to (11). Compute its distance to {h 0n j } n=1...N j=1...J . If the distance is smaller than ✏ 3 , then STOP. If not, we repeat step 2 using as initial guesses {h n j } n=1...N j=1...J and {l n j } n=1...N j=1...J just computed.

Numerical exercises

In this section we apply the numerical method introduced before, focusing on the emergence of spatial patterns and its corresponding drivers. As observed in the Introduction, many papers have empirically identified that the interaction between land-use activities and the environment can explain the dynamics of spatial patterns. In this regard, a central concern of this literature is the optimal spatial allocation of land-use activities. We will use our framework to better understand the mechanisms behind this problem and, in particular, the dynamics of the corresponding spatial heterogeneity. Even if the main objective of the paper is to provide a benchmark set-up, we will show that our simplified model already reproduces an ample variety of spatial heterogeneity scenarios related to the interaction between land-use and the environment. In particular, we will analyse the persistence in time of spatial heterogeneity and the subsequent emergence of specialized areas, which are defined in our model as locations where the majority of their available land is devoted to production or abatement. This spatial dynamic outcome is similar to the creation of (temporal or permanent) biological reserves in [START_REF] Costello | Optimal harvesting of stochastic spatial resources[END_REF]. From a di↵erent perspective (harvesting, partial equilibrium and dispersal of natural resources due to advection), they show that preventing areas from harvesting is optimally justified. This result is equivalent in our context to the specialization of some locations (areas) in abatement. As in their paper, the spatial connectivity (due to the dispersal process) and the particular characteristics of each location will play a fundamental role in this respect. Finally, our model will also point out that abatement technology stands out as an important ingredient to reach time-invariant solutions. The underlying mechanism is the idea of "flux equilibrium" in [START_REF] Smith | The economics of spatial-dynamic processes: Applications to renewable resources[END_REF]: variables are constant, but which maintains the equilibrium is a di↵usion flux. From an economic point of view, this particular equilibrium is interesting because pollution is stabilized (i.e., both local and global pollution become constant) in an economy that can eventually sustain a constant consumption per capita. Moreover, this type of equilibrium can be compatible with the formation of long-run spatial patterns as well.

The numerical exercises are divided in two parts. Sections 5.1-5.3 consider that population is uniformly distributed, while Section 5.4 assumes a Gaussian distribution in order to study the e↵ect of population agglomeration. The parameter values are provided in Table 1. For illustration purposes we consider that the land endowment of each location, L(x), is equal to 300, and that the total population of our economy is equal to 110. 12 We would like to underline that the values provided in this table aim at illustrating our model, and they do not correspond to any specific situation since we shall focus on the qualitative properties of our set-up. We assume that the space is a line of length 5 divided into 500 locations. The time horizon varies from 10 to 40 depending on the convergence speed of the variables. Agents preferences are given by a logarithmic utility function. We consider a Cobb-Douglas production function, where the net productivity is B + A⌦(p, P, x) with ⌦(p, P, x) = e 2 p 1 P s(x) . Following Weitzman (2009), ⌦ is an exponential damage function, taking values in the interval [0, 1]. We consider that both local and global pollution harm productivity, where 1 and 2 are constants: for given a (p, P ), the fraction 1 ⌦(p, P, x) represents the foregone productivity at location x. For the sake of simplicity we assume that A and B are both constant in space and time. Moreover, s(x) stands for the sensitivity of location x to global pollution. Assuming a linear abatement technology, we have G(l) = l. Finally, consistently with the example of function provided in Section 3, (p, P )(x, T ) = p(x, T ), we consider = 5000 in order to emphasize the policy maker's concern about the pollution left at the end of the planning period.

Notice that we consider in all scenarios that initial pollution is uniformly distributed.

We believe of no particular interest the case when the only spatial feature is the initial distribution of pollution. Indeed, any di↵erence in the initial endowment of pollution vanishes with time if all other variables are spatially homogeneous.

The benchmark scenario

We begin our analysis with the benchmark scenario in which population is evenly distributed on space. It is the objective of this benchmark illustration to underline the trade-o↵ between production and abatement. Accordingly, we have reduced the amount of land devoted to housing by means of considering a uniform distribution of population that results in 0.22 individuals per location. This implies that each location needs 0.22 units of land for housing, which is clearly not critical when the total land endowment is 300. 13 We further assume that the spatial sensitivity to pollution is constant in space,

i.e., s(x) = 1 for all x. Figure 1 shows the results. Given that there are no spatial disparities, it is not surprising that the optimal trajectories are uniform in space. The allocation of land to production starts at its highest possible level and remains unchanged until the environmental damage is large enough.

At this point, land to production is optimally reduced and, consequently, the economy devotes part of the land endowment to abatement activities. Consumption shows a decreasing trajectory due to the pollution damage of production and the replacement of land to production by abatement. Notice moreover that it eventually reaches a constant level, while local and global pollution continuously increases.

The optimal land trajectory attains a homogenous constant value too. Despite of using 2/3 of land to production, the economy cannot keep its initial level of consumption in the long-run due to the damage caused at the beginning. Both types of pollution cause indeed everlasting and increasing damage that the current abatement is not able to completely eliminate. However, if the e ciency of abatement is large enough our model shows that pollution can be stabilized. This outcome is illustrated in Figure 2,

where the abatement e ciency parameter ( ) is equal to 0.9. As we can see in this figure, the economy reaches a time-invariant solution. In contrast to Figure 1, (local and global) pollution becomes constant after some periods, together with consumption per capita and land to production. Consumption per capita, moreover, always increases during the transition, while a decreasing trajectory arises in the benchmark scenario due to the growing contamination. Notice also that in Figure 2 we have chosen a particularly high abatement e ciency (nine times the e ciency considered in Figure 1) in order to point out e↵ect pollution abatement. A direct consequence of this assumption is that pollution eventually disappears. Finally observe that, even if the focus of our simulations is the quality behavior of the economy, the substantially lower consumption per capita in Figure 2 can be explain be means of taking into account the planner's concern about the pollution at the end of the planning period. In the benchmark scenario, the welfare reduction due to the large amount of contamination left is compensated with a high levels of consumption. This compensation is not important though in scenarios where the pollution is very reduced such as in the scenario considered in Figure 2.

Let us study in the next sections the emergency of spatial patterns and the implications of the di↵erent elements of our model in this regard.

Role of abatement technology

We consider here a simple case of heterogeneous abatement technology in which abatement e ciency continuously deteriorates as we get afar from x = 0:14 (x) = 0.1 + 0.19 1 + e x 2.5 . This logistic form can be interpreted as a continuous representation of a step function, where some locations are better suited for abatement activities than others. In our particular example the abatement e ciency parameter (x) monotonically decreases, ranging from about 0.3 to 0.1. The results are displayed in Figure 3.

We can observe that the heterogeneity in abatement induces heterogeneity in land allocation from the beginning. Indeed at time t = 0, the less advanced locations in abatement specialise in production, whereas locations close to x = 0 focus on abatement.

Due to this specialization, the areas devoted to production face greater levels of (local) pollution than the locations where abatement activities are intensified. Notice also that we have improved the abatement technology in all locations, with an e ciency that more than doubles for the most suited areas. As a result and in contrast to the benchmark case, locations compensate for emissions and the economy reaches a timeinvariant solution. This outcome actually points out the role of abatement as pollution stabilizer. Moreover, in the same direction of Figure 2, pollution is reduced in areas where the abatement is su ciently e cient.

Long-term consumption takes exactly the same value as in the benchmark, although consumption monotonically increases from the start as a direct consequence of the accrued abatement. Locations specialised in abatement produce little. This greater abatement e↵ort allows them to compensate for their relatively unimportant emissions and the incoming pollution from other locations that are better qualified for production activities. However, regardless of locations' production and/or abatement, the possibility of having consumption "imports", as described in (3), enables homogeneous consumption and specialisation. Finally, since the time-invariant equilibrium is spatially heterogenous, we can also conclude that permanent di↵erences in abatement technology allow for lasting heterogeneity in land allocation and local pollution.

• Local and global damage

We have considered in the previous scenarios that both local and global pollution causes the same damage per unit, i.e., 1 = 2 . Consistently however with the examples provided in sections 1 and 2, our model also allows us to study the case of contaminants with only local or only global e↵ects. 15When the damage is only local 1 is equal to zero in ⌦. Since in this case the damage does not depend on global pollution, which is the largest pollutant by definition, the total damage of pollution is lower than in the previous scenario. As a consequence, one can see in Figure 4 that at first no location abates. Nevertheless, specialization emerges when the level of local pollution is high enough. The economy eventually reaches a time-invariant equilibrium, which is qualitatively identical to the previous case. However, the levels of local and global pollution are higher because of a lower damage of pollution. This result also points out that if the contamination damage is high enough (Figure 3) pollution will be optimally lower than when its damage is small (notice that local pollution is even reduced in some areas of the space in Figure 3). Parallel to the transition dynamics in figures 1 and 2, consumption per capita decreases until its time-invariant level, while this trajectory is increasing in Figure 3, where pollution is lower.

One should also observe that the rise of spatial heterogeneity is postponed until the economy accumulates enough contamination. This can also explain why consumption is initially higher than in the previous case: land devoted to production is higher and pollution damage is lower. We therefore conclude that, due to a lower pollution harm, the absence of global damage can delay the emergence of spatial patterns. This is an interesting dynamic property of our framework. Among the di↵erent scenarios that our simple set-up can reproduce, we can also include cases of delayed spatial heterogeneity.

From an economic perspective, this outcome points out the spatial-dynamic dimension of the problem. Even if there is spatial connectivity, the accumulation e↵ect (of pollution in our model) should be taken into account in order to fully understand how a particular element (abatement e ciency in this scenario) can induce spatial heterogeneity.

Let us consider the situation where the damage is only global ( 2 = 0). For our parameters values, and among the two previous exercises, this situation corresponds to an intermediate case of pollution damage. On the one hand, Figure 3 represents the case where pollution has both local and global e↵ects, so the resulting damage is the highest. On the other hand, in Figure 4 the damage of pollution is the lowest because its e↵ect is only local. We can thus expect that the response of the economy when the damage is only local will be a combination of these two cases. As it is clear from Figure 5, consumption and global pollution initially behave like in case of the lowest pollution damage, that is when pollution has only local e↵ects. The reason of this similarity is that pollution takes time to accumulate so that productivity losses are postponed, and the pollution damage is lower than in the case where pollution has both local and global e↵ects. Still, in contrast to Figure 4, land specialisation emerges from the beginning because (by definition) global pollution is always greater than the local one. Notice that, because the damage of pollution is lower, the econ-omy produces and consumes more than in the case considered in Figure 3 despite the abatement specialization of part of the space.

Figure 5 also shows that consumption reduces in the starting periods. But later on, when pollution is large enough, the economy evolves like in the case where pollutions has both local and global e↵ects (see Figure 3). Since the locations specialised in abatement produce very little, pollution reduces in these areas, while it rises in the locations that focus on production. Nevertheless, even if its e↵ect is only global, pollution moves gradually in space until it reaches the areas devoted to abatement. Consequently, after some periods, local and global pollution stabilize, and consumption increases until its time-invariant level as in the case corresponding to Figure 3. Just to conclude this section, let us point out two interesting features that the last numerical exercise reveals. First, one could believe that the global nature of pollution tends to homogenize space. Our example shows, however, that spatial heterogeneity can emerge even when pollution only has global e↵ects, due to pollution di↵usion and the spatial specificity of abatement activities. Second, in contrast to the previous two cases, pollution di↵usion can generate transitional dynamics that are non-monotonic (see Figure 5). This property underlines that, despite the simplicity of our model, we can provide scenarios with spatial heterogeneity and complex dynamics in time.

Spatially heterogeneous sensitivity to pollution

We consider the situation where some areas of the space are more sensitive to pollution than others. This case illustrates, for instance, the impact of pollution on global warming and the subsequent rise of the sea level. Many authors have recognized the importance of this negative e↵ect of pollution and, in particular, the associated degradation of soil quality (for instance, [START_REF] Nicholls | Sea-Level Rise and Its Impact on Coastal Zones[END_REF]. Similarly, the desertification of drylands gives us another example of spatial heterogeneity related to di↵erences on pollution sensibility (among others, [START_REF] Smith | Global Desertification: Building a Science for Dryland Development[END_REF]. In both cases, global warming is usually associated with the increase of global pollution such as the anthropogenic GHGs.

In our simplified set-up, we can study this problem by means of assuming that the sensitivity to global pollution s(x) in the damage function ⌦ is spatially heterogenous. We specifically consider the following sensitivity function:

s(x) = 10 1 + e 0.025 x 4 ⇣ 1 x 5 ⌘ .
As before, we consider a logistic function. But in this scenario locations are more sensitive to global pollution as they get afar from x = 0, with a sensitivity parameter ranging from about 1 to 10. Moreover, we assume greater concavity in order to emphasize the environmental sensitivity e↵ect (i.e., a relative large number of "fragile" locations). The numerical results are presented in Figure 6. We can observe that production is initially larger in less sensitive locations. Land devoted to production decreases indeed as one gets afar from x = 0. However, di↵erences in land allocation eventually vanish and the amount of land assigned to production reaches a constant spatially homogeneous level. This result goes against the a priori belief that the most sensitive regions would produce less than the others (and, consequently, "import" most of their consumption) in order to preserve their environmental quality. The explanation of this homogeneity outcome is the following. Since pollution flows across locations, even the regions with non-existent or little production will experience positive levels of local pollution.

Moreover, the pollution as a whole (global pollution) damages production too. Due to these two sources of pollution damage, the less sensitive locations optimally reduce their production, devoting as well land to abatement. If the most sensitive locations were endowed with better abatement technology they would then dedicate more land to abatement relatively to the less sensitive locations.

Let us finally point out that the numerical exercise considered in this section also illustrates an interesting spatial dynamic feature: an initial spatial heterogeneity (in land assigned to production, due to the spatial di↵erences in pollution sensitivity) can vanish in the long-run. In contrast to the previous scenarios, we are showing here that the di↵usion forces can also drive spatial homogeneity. Our numerical illustration actually presents an economy that eventually converges to the benchmark scenario, which is mainly characterized by its spatial homogeneity (see Figure 1).

The e↵ect of population agglomeration

Up to now our numerical exercises focused on the trade-o↵ between production and abatement. Accordingly, we have minimized in the previous sections the constraint of retaining some land for housing. As a last experiment, let us then analyse the e↵ect of population agglomeration and the resulting housing requirement. We consider in this regard that population is distributed following a Gaussian function over the interval [0, 5], that is, population agglomerates around the center of the space, x = 2.5. In order to underline the e↵ect of population agglomeration, we set total population to 10500 so that population in x = 2.5 is almost 130. Consequently, although the land endowment of each location is still equal to 300, in the central area of the space the proportion of L devoted to housing is much larger than in previous scenarios due to accrued population. 16 This contrast with the locations far away from the center, where the weight of population is similar to that in the benchmark scenario.

Let us first compare the optimal trajectories under population agglomeration with the benchmark scenario. Figure 7 shows that, due to population concentration, locations in the central area cannot devote as much land to production as the locations at the far ends. This means that agglomerations optimally "import" most of their consumption from the neighbouring areas, which are more specialised in production. One could arguably think that agglomerations would be less locally polluted because most of the production comes from the periphery and housing does not involve emissions in our simplified framework. However, by the same token, agglomerations cannot devote as much land to abatement as the rest of locations. We thus observe a heterogeneous 16 Notice that, in the previous scenarios, this increase in total population is sizable. However, a homogenous distribution of 10500 people over 500 locations would imply 21 individuals per location.

In our simplified setup, 21 individuals would need 21 units of land for housing, which still is a small figure with respect to the land endowment of each location. distribution of land to production, with almost no abatement in the center of the space.

Consequently, local pollution in the central area is not lower than in other locations.

We actually find that slight spatial disparities persist since agglomerations cannot abate pollution coming from neighbouring regions.17 This point is reinforced in the experiment considered in Figure 8.

In this second exercise we have doubled abatement e ciency in all locations (i.e., (x) = 0.2 for all x). In e↵ect, due to this technological improvement, all locations devote some land to abatement from the beginning. Both local and global pollution levels decrease then, allowing for a greater consumption per capita in the long-term.

However, spatial disparities are amplified since the abatement capacity of the central area is constrained by its housing requirement.

We should finally observe that in this last scenario all variables reach a time-invariant solution, which is characterized by lasting spatial heterogeneity in both land allocation and local pollution. As in Section 5.2, this result points out again the role of abatement as pollution stabilizer. Abatement e ciency indeed enhances consumption and enables the economy to reach a time-invariant equilibrium, which can be spatially heterogenous.

Further comments

Let us conclude the simulations of our paper with an additional discussion about the time-invariant solution and, in particular, its uniqueness as stated in Proposition 4.

We will also complete the numerical analysis with a robustness check of the algorithm regarding the optimal trajectories.

In Section 3.2 we have studied the time-invariant equilibrium. As it is clear from the previous simulations, we have found several cases where the economy ends up in a time-invariant solution. Our algorithm actually provides a numerical tool to analyse the convergence to this kind of equilibrium. However, the multiplicity of this kind of long-term solution cannot be a priori ruled out. Still, in this regard, Proposition 4 in Section 3.2 turns out to be very useful since it allows us to identify a su cient condition for uniqueness of time-invariant solutions. Since this condition was originally stated for the case when production is described by A⌦F , we must adapt it to the simulated case where production is given instead by ( B + Ã ⌦)F . Taking A⌦ = ( B + Ã ⌦), the condition is rewritten as:

⌦11 (p, P ), ⌦21 (p, P ) > 0 and BF ( l) + ÃF ( l)[ ⌦(p, P ) p ⌦1 (p, P )] > G(1 f l) at every x. We can then apply it to our numerical illustrations. The condition is in fact verified in our simulations for every case where convergence towards a time-invariant interior solution is observed. Proposition 4 therefore ensures that this equilibrium is the unique time-invariant solution.

About the optimal paths, regardless the convergence to a time-invariant equilibrium, we should point out that the uniqueness property of the trajectories is still a mathematical open question. Therefore, since our problem may have more than one optimal solution, we may wonder to which extent the solutions presented in this section depend on the set of initial guesses. We have then performed several robustness checks in this respect. In these exercises we modify the initial guesses for the shadow price of pollution, land devoted to production, and aggregated consumption, in configurations with homogeneous or heterogenous distribution of population, abatement technology, and sensitivity to pollution. The results confirm that our algorithm is robust and always generates the same optimal trajectories.

More specifically, recall that in our numerical exercises the reversed-time shadow price of pollution {h 0n j } n=1...N j=1...J was set to -5000. We have run simulations where {h 0n j } n=1...N j=1...J ranges from 4750 to 100, leaving all else equal. We have similarly varied the initial guess for land to production from 25 to 250. Finally, we have considered in the simulations that the initial guess for aggregated consumption was about 135. We so vary this value from 25 to 200. In all cases, the solution trajectories for local and global pollution, as well as for land distribution, coincide with those presented in the previous subsections. We thus conclude from these results that our algorithm is robust with respect to the initial guesses.18

Concluding remarks

The main objective of this paper is to present a benchmark framework to study optimal land use, encompassing land use activities and pollution. In our model, although land is immobile by nature, local actions a↵ect the whole space through pollution, which flows across locations resulting in both local and global damages. We find that our benchmark model reproduces a great variety of spatial patterns related to the interaction between land use activities and the environment. In particular, we identify the central role of abatement technology as pollution stabilizer, allowing the economy to achieve invariant solutions with time, which can be spatially heterogeneous.

Several extentions can be made to our basic set-up, ranging from considering endogenously distributed population [START_REF] Papageorgiou | Agglomeration as local instability of spatially uniform steady-states[END_REF][START_REF] Marchiori | When Nature Rebels: International Migration, Climate Change and Inequality[END_REF] to the empirically usage of this type of models to estimate structural spatial-dynamic parameters [START_REF] Smith | The economics of spatial-dynamic processes: Applications to renewable resources[END_REF]. But let us particularly mention that the decentralisation of the social optimum, in the spatial Ramsey model, has not been explored yet in the literature. In this regard, a challenging extension could study the possibility of optimal tax/subsidy schemes that will evolve with time but also across the space. The spatial dependence is indeed consistent with numerous papers suggesting that the optimal policies should take spatial information into account (for instance, [START_REF] Hochman | The E↵ects of Pollution Taxation on the Pattern of Resource Allocation: The Downstream Di↵usion Case[END_REF][START_REF] Weitzman | Landing fees vs harvest quotas with uncertain fish stocks[END_REF][START_REF] Costello | Optimal harvesting of stochastic spatial resources[END_REF]. This property would raise another important question about the implementation of optimal policies in a spatial dynamic context: in terms of social welfare, how far away a second best solution, without spatially di↵erentiated taxes/subsidies, would be from a spatially heterogeneous first best (see Smith et al., 2009, among others).

Appendices A Proposition 1 proof

We shall show that the system of PDEs constraining the policy maker's objective function has a unique solution for every choice of feasible functions (c, l). This proves indeed the existence of at least a solution to the social optimum problem. On this matter, we begin by converting the set of constraints into a system of parabolic di↵erential equations.

First, notice that we can take the derivative of P with respect to t and we then use the law of motion for p in P to obtain:

P t (t) = Z R p t (x, t)dx = Z R [p xx (x, t) + ⌦(p, P, x)A(x, t)F (l(x, t)) G(1 l f (x))] dx, which implies that P t (t) = lim x! R p x (x, t) lim x! R p x (x, t)+p x (x, t)| R + Z R [⌦(p, P, x)A(x, t)F (l(x, t)) G(1 l f (x))] dx, where p x (x, t)| R = lim x! R p x (x, t) lim x! R p
x (x, t), and R and R denote, respectively, the upper and lower boundary of our unidimensional region R. Since lim

x! R p x (x, t) = 0, we have that

P t (t) = Z R [⌦(p, P, x)A(x, t)F (l(x, t)) G(1 l f (x))] dx. (A.1)
Now, our initial set of constraints can be written as a system of parabolic equations. We can indeed interpret (A.1) as a partial di↵erential equation, with the second order operator equal to zero. We would need to artificially transform P into a two dimensional function, P (x, t) = P (t), 8x 2 R. Then:

(P 0 ) 8 > > > > > > > > > > > > > > > < > > > > > > > > > > > > > > > : p t (x, t) p xx (x, t) = ⌦(p, P, x)A(x, t)F (l(x, t)) G(1 l f (x)), P t (x, t) = R R [⌦(p, P, x)A(x, t)F (l(x, t)) G(1 l f (x))] dx, p(x, 0) = p 0 (x) 0, lim x! R p x (x, t) = 0, P (x, 0) = R R p 0 (x)dx, lim x! R P x (x, t) = 0, (A.2)
for all (x, t) 2 R ⇥ R + . As in [START_REF] Camacho | On the dynamics of capital accumulation across space[END_REF] and [START_REF] Boucekkine | Bridging the gap between growth theory and the new economic geography: The spatial Ramsey model[END_REF], after transforming the integral term in each dynamic equation, we make use of [START_REF] Pao | Nonlinear Parabolic and Elliptic Equations[END_REF] to prove the existence of a solution to this kind of equations for any (x, t) 2 R ⇥ (0, T ], with T < 1.

In this regard, we proceed with the following change of variable ⇥(x, t) = e t P (x, t) for any > 0, and we obtain:

⇥ t (x, t) + ⇥(x, t) = e t Z R ⇥ ⌦(p, e t ⇧, x)A(x, t)F (l(x, t)) G(1 l f (x)) ⇤ dx.
Then, we define a function ✓(t) as

✓(t) = e t Z R ⇥ ⌦(p, e t ⇥, x)A(x, t)F (l(x, t)) G(1 l f (x)) ⇤ dx.
Notice that since the integrand is globally Lipschitz continuous, so it is function ✓. We have to study now the existence of solution of the following system of equations:

8 > > > > > > > > > > > > > > > < > > > > > > > > > > > > > > > : p t (x, t) p xx (x, t) = ⌦(p, e t ⇧, x)A(x, t)F (l(x, t)) G(1 l f (x)), ⇥ t (x, t) + ⇥(x, t) = ✓(t), p(x, 0) = p 0 (x) 0, lim x! R p x (x, t) = 0, ⇥(x, 0) = R R p 0 (x)dx, lim x! R ⇥ x (x, t) = 0. (A.3)
In this respect, we apply Theorem 12.1 in Chapter 8 in [START_REF] Pao | Nonlinear Parabolic and Elliptic Equations[END_REF] in order to ensure the existence of a unique solution to the system of parabolic equations for every choice of the couple (c, l).

B Proposition 2 proof

Function ( 6) is the value function associated to problem P. V is a function of c, l, p and P .

If there exists an optimal solution (c ⇤ , l ⇤ , p ⇤ , P ⇤ ), then any other solution to problem (4)-( 5) can be written as a deviation from the optimal solution as

c(x, t) = c ⇤ (x, t) + ✏(x, t), l(x, t) = l ⇤ (x, t) + ✏L(x, t), p(x, t) = p ⇤ (x, t) + ✏⇡(x, t), P (t) = P ⇤ (t) + ✏⇧(t). (B.1)
We can take the first order derivative of the value function V with respect to ✏ in order to minimize the deviation of the trajectory from the optimal. Beforehand and using integration by parts, we re-arrange some integral terms in V :

R T 0 R R q(x, t)p xx (x, t)dxdt = R T 0 q(x, t)p x (x, t)| R dt R T 0 q x (x, t)p(x, t)| R dt + R T 0 R R q xx (x, t)p(x, t)dxdt, (B.2)
and as usual:

R T 0 R R q(x, t)p t (x, t)dxdt = R R p(x, t)q(x, t)| T 0 dx R T 0 R R p(x, t)q t (x, t)dxdt = R R p(x, T )q(x, T )dx R R p(x, 0)q(x, 0)dx R T 0 R R p(x, t)q t (x, t)dxdt. (B.3)
We then obtain:

@V (c,l,p,P ) @✏ = R T 0 R R u 0 (c(x, t))f (x)e ⇢t (x, t)dxdt + R R p (p(x, T ), P (T ))⇡(x, T )e ⇢T dx + R R P (p(x, T ), P (T ))⇧(T )e ⇢T dx + R T 0 R R ⇡(x, t) [q t (x, t) + q xx (x, t)] dxdt R R q(x, T )⇡(x, T )dx R T 0 ⇡(x, t)q x (x, t)| R dt + R T 0 R R q(x, t)⌦ 1 (p, P, x)A(x, t)F (l(x, t))⇡(x, t)dxdt + R T 0 R R q(x, t) [⌦ 2 (p, P, x)A(x, t)F (l(x, t))⇧(t) + ⌦(p, P, x)A(x, t)F 0 (l(x, t))L(x, t)] dxdt + R T 0 R R q(x, t)G 0 (1 l f (x))L(x, t)dxdt R T 0 m(t) ⇥ ⇧(t) R R ⇡(x, t)dx ⇤ dt R T 0 n(t) ⇥R R (x, t)f (x)dx ⇤ dt + R T 0 n(t){ R R [⌦ 1 (p, P, x)AF (l)⇡(x, t) + ⌦ 2 (p, P, x)AF (l)⇧(t) + ⌦(p, P, x)AF 0 (l)L(x, t)] dx}dt.
To get the necessary conditions, we can group the elements multiplying , ⇡, L and ⇧, and equate them to zero. If all factors multiplying deviations from optimal values for c, p, P and l are equal to zero, then @V @✏ = 0. We would need:

8 > > > > > < > > > > > :  : u 0 (c)e ⇢t = n(t), ⇡ : q t + q xx + (q + n) ⌦ 1 AF (l) + m = 0, ⇧ : m(t) = 1 f (x) ⌦ 2 AF (l) (q + n) , L : q (⌦AF 0 + G 0 ) + n(t) (⌦AF 0 ) = 0. (B.4)
To these conditions, we need to add the following spatial boundary and transversality condi-

tions: 19 8 > > < > > : lim x! R q x p = 0, q(x, T ) = p (x, T ), R R P (x, T )dx = 0.
After detrending the co-state variables and substituting m(t) by 1 f R⌦ 2 AF (l) (q + n) into the dynamic equation for q, we obtain the set of necessary conditions presented in the proposition. As usual, abusing of notation, we denote in the statement of Proposition 2 the detrended co-state variables as the original ones.

19 Notice that, if we had assumed a di↵erent boundary condition for p in (5), the necessary condition on the border would have been: lim

x! R [q(x, t)p x (x, t) p(x, t)q x (x, t)] = 0.

if and only if u ⇤ xx f (x, u ⇤ ). Similarly, u ⇤ is a lower solution of (E.1) if and only if

u ⇤xx  f (x, u ⇤ ).
Indeed, it is evident from the system S that the time-invariant solution of p and q is described by two parabolic equations as in (E.1). We can then state and prove a two-dimensional version of Theorem 3.4 in [START_REF] Pao | Nonlinear Parabolic and Elliptic Equations[END_REF]:

Theorem 1. Let (p ⇤ , q ⇤ ), (p ⇤ , q ⇤ ) be ordered upper and lower solutions of p xx = f 1 (p, q) and q xx = f 2 (p, q), such that p ⇤ p ⇤ 0 and q ⇤ q ⇤ 0.

Let (p, q), (p, q) be positive maximum and minimum solutions with p, p 2< p ⇤ , p ⇤ > and q, q 2< q ⇤ , q ⇤ >. Assume that f = (f 1 , f 2 ) 2 C 1 in < p ⇤ , p ⇤ > ⇥ < q ⇤ , q ⇤ >. If f 1 satisfies either 8 < : @f 1 @q 0 @(f 1 /p) @p 0 or 8 < : @f 1 @q  0 @(f 1 /p) @p  0 and f 2 verifies either 8 < : @f 2 @p 0 @(f 2 /q) @q 0 or 8 < : @f 2 @p  0 @(f 2 /q) @q  0 then p = p, q = q and is the unique positive solution in < p ⇤ , p ⇤ > ⇥ < q ⇤ , q ⇤ >.

Proof. Let L denote the parabolic operator, that is Ly = y xx . We know that since (p, q), (p, q) are two solutions to S: Lp = f 1 (p, q), Lq = f 2 (p, q), and Lp = f 1 (p, q), Lq = f 2 (p, q). Now we treat the time-invariant problem for p and q separately. Let us multiply the first equation above for p by p, and the second one for p by p: pLp = pf 1 (p, q), pLp = pf 1 (p, q).

We subtract the second to the first: pLp pLp = pf 1 (p, q) pf 1 (p, q), and taking integrals in R:

Z R pLp pLp dx = Z R pf 1 (p, q) pf 1 (p, q) dx = Z R pp ✓ f 1 (p, q) p f 1 (p, q) p ◆ dx.
By Green's Theorem R R pLp pLp dx = 0. Let us prove the unicity of the time-invariant equilibrium for p in the first case, when @f 1 @q 0 and @(f

1 /p) @p 0. Since 0 = R R pp ⇣ f 1 (p,q) p f 1 (p,q) p ⌘ dx R R pp ⇣ f 1 (p,q) p f 1 (p,q) p ⌘ dx 0 if @(f 1 /p) @p
0, then we have that p = p. By a similar procedure, one can prove that the equality holds under the second set of conditions for f 1 , namely, @f 1 @q  0 and @(f 1 /p) @p  0. Moreover, following the same reasoning, we can prove that q = q under either set of conditions on f 2 .

Let us apply Theorem 1 to our problem. In this respect, we shall perform a change of variable q = q since all results apply to positive solutions and we do know that q is negative.

We can then study the existence and uniqueness of time-invariant solutions to S 0 : Then, in problem S 0 , f 1 (x, p, q) = ⌦(p, P, x)A(x)F (l(x)) G(1 l f (x)) and f 2 (x, p, q) = ⇣ ⌦ 1 (p, P, x) + 1 f (x) ⌦ 2 (p, P, x) ⌘ A(x)F (l) [u 0 (c) q(x)] + ⇢q(x).

In order to compute the derivatives of f 1 and f 2 with respect to p and q, we need to compute first @l/@p and @l/@ q. Let us rewrite the condition that defines l in (D.1) (see Appendix D):

J(l, p, q) = [u 0 (c) q]⌦AF 0 (l) qG 0 .

Applying the Implicit Function Theorem: dl dp = @J/@p @J/@l = (u 0 q)⌦ 1 AF 0 (u 0 q)⌦AF 00 +qG 00 < 0, dl dq = @J/@ q @J/@l = ⌦AF 0 +G 0 (u 0 q)⌦AF 00 +qG 00 < 0.

Since @f 1 @ q = (⌦AF 0 + G 0 ) dl dq  0, we need @(f 1 /p) @p  0 in order to satisfy Theorem's 1 hypothesis. Let us compute this derivative:

@(f 1 /p) @p = h AF ⌦ 1 + (⌦AF 0 + G 0 ) dl dp i p (⌦AF G) p 2 .
Thus, we need to check whether

 AF ⌦ 1 + (⌦AF 0 + G 0 ) dl dp p < ⌦AF G.
It su ces to impose AF (⌦ p⌦ 1 ) > G to ensure the negativeness of @(f 1 /p) @p . In particular, if ⌦AF G 0, then the above condition is trivially verified.

Next, let us compute the partial derivative of @(f 2 /q) @ q , which under the problem assumptions is negative: @(f 2 /q) @ q = A(⌦ 1 + ⌦ 2 /f ) q  (u 0 q)F 0 dl dq + F u 0 q  0.

Thus, since @(f 2 /q) @ q  0 we need to find conditions under which the derivative of f 2 with respect to p is negative in order to comply with the theorem assumptions.

@f 2 @p = ✓ ⌦ 11 + 1 f ⌦ 21 ◆ AF (u 0 q) ✓ ⌦ 1 + 1 f ⌦ 2 ◆ AF 0 (u 0 q) dl dp
is negative if we assume, for instance, that ⌦ 11 , ⌦ 21 > 0.

Summing up, we have shown that @f 1 @ q  0, @(f 1 /p) @p  0, @f 2 @p  0 and @(f 2 /q) @ q  0 if (p, P , f ) is a three dimensional vector in R 3 , which gathers the values of p(x), P (x) and f (x) that maximize ⌦AF . We then define % as the evaluation of where ⇣ = sup l(x)2[0,f (x)] G(1 l(x) f (x)) = sup x G(1 f (x)), and ⌘ > 0 is a real number. One can easily check that the proposed couples do define upper and lower solutions to S 0 . It remains however to check that the lower solution is smaller than the upper solution. First, q⇤ (x) < q⇤ (x) since ⌘ can take su ciently large real values. However, in order to ensure that p ⇤ (x) < p ⇤ (x) we need to impose the following constraint on the constant ◆: ⇣ 2 < (◆ + 1).

h ⇢ ⇣ ⌦ 1 + 1 f ⌦ 2 ⌘ AF i u 0 (c)
Since we can prove the existence of at least an upper and a lower solution to the stationary problem S 0 , and that f 1 and f 2 comply with the requirements of Theorem 1 under the condition (E.2), then we can ensure the existence of a unique solution to the stationary problem S 0 and thus to S.

Figure 1 :

 1 Figure 1: Benchmark scenario.

Figure 2 :

 2 Figure 2: Pollution stabilization.

Figure 3 :

 3 Figure 3: Role of abatement technology.

Figure 4 :

 4 Figure 4: Damage function only depends on local pollution ( 1 = 0).

Figure 5 :

 5 Figure 5: Damage function only depends on global pollution ( 2 = 0).

Figure 6 :

 6 Figure 6: Spatially heterogeneous sensitivity to pollution.

Figure 7 :

 7 Figure 7: The role of population agglomeration (Gaussian distribution).

Figure 8 :

 8 Figure 8: Population agglomeration with abatement e ciency doubling.

  Weitzman M. (2009), "On modeling and interpreting the economics of catastrophic climate change", The Review of Economics and Statistics, 91(1), 1-19.

  = ⌦(p, P, x)A(x)F (l(x)) G(1 l f (x)), qxx (x) = ⇣ ⌦ 1 (p, P, x) + 1 f (x) ⌦ 2 (p, P, x) ⌘ A(x)F (l) [u 0 (c) q(x)] + ⇢q(x).

  order to apply Theorem 1, it only remains to find an upper and a lower solution to the stationary problem S 0 under (E.2). We can define the upper solution as (p ⇤ , q⇤ )P , f ) , with (p, P , f ) = arg max {⌦AF } .Note that the supremum that defines exists since the function ⌦ takes values in [0, 1].

  at this supremum values for (p, P, f ).Next, let us define a lower solution (p ⇤ , q⇤ ):

Table 1 :

 1 Parameters values for the numerical exercises.

	B Minimum productivity	0.5
	A Max. productivity increase 10
		Abatement e ciency	0.1
	⇢	Time discount rate	0.05
	1	P damage	0.005
	2	p damage	0.005
	↵ Cobb-Douglas parameter	0.75
	p 0 Initial pollution at x	100

According to[START_REF] Akimoto | Global Air Quality and Pollution[END_REF], tropospheric ozone, methane and CO are well-known examples of pollutants that flow across locations. Methane and CO have both local and global e↵ects. Moreover, CO a↵ects the oxidizing capacity of the atmosphere, raising the lifetime of GHGs.

Since we consider a finite planning period, we prefer to use the term "time-invariant". The designation of "steady-state" is commonly employed in infinite horizon contexts.

The Ekeland variational principle[START_REF] Ekeland | On the variational principle[END_REF] ensures the existence of a maximum value for V in R, when V is upper semicontinous. Hence, any feasible solution (c, l, p, P ) verifying (5) can be written

This method is called the implicit finite di↵erence approximation (for instance,[START_REF] Smith | Numerical solution of partial di↵erential equations[END_REF][START_REF] Sewell | The Numerical Solution of Ordinary and Partial Di↵erential Equations[END_REF]. Other approximation schemes exist but the implicit one is unconditionally stable, meaning that it is stable without restrictions on the relative size of t and x. It also allows us to use a larger time step and to save this way computational time.

Notice that space is finite in numerical exercises. This implies that total population does not need to be equal to 1 since the convergence of the integral in the objective function is ensured. Therefore, taking advantage of this property, we increase both total population and land endowment in order to enlighten our numerical results.

We will consider the e↵ect of population agglomeration and the subsequent accrued need for housing in Section 5.4.

For empirical evidence of di↵erences in abatement technology see, for instance,[START_REF] De Cara | Methane and nitrous oxide emissions from agriculture in the EU: A spatial assessment of sources and abatement costs[END_REF] and[START_REF] Nowak | Air pollution removal by urban trees and shrubs in the United States[END_REF].

The results of these scenarios are qualitatively equivalent to the case of pollutants with mainly local ( 2 > 1 ) or mainly global ( 2 < 1 ) e↵ect. Obviously, if 1 = 2 = 0 no land will be devoted to abatement since pollution does not damage our economy. Therefore, consumption will stay at its maximum constant level (after taking housing into account, the remaining land will be completely assigned to production), where both local and global pollution increase steadily.

Pollution due to housing and/or transportation would amplify this e↵ect. These additional sources of contamination may have potential interesting implications, in particular if labour were a spatially mobile production factor.

Detailed results of the robustness checks can be provided upon request.

C Corollary 1 proof

As we can see in the first equation of the system (B.4), proof of Proposition 2, u 0 (c(x, t)) = n(t)e ⇢t , for all (x, t). Hence, neither u 0 (c(x, t)) nor c(x, t) depend on space.

D Proposition 3 proof

We show next that l is a unique function of p, P and q, or since P (t) = R R p(x, t)dx, that l is a unique function of p and q. From the last equation of (B.4), we can verify that q(x, t)  0 under the assumptions (H1) and (H2), and provided that n(t) = u 0 (c(t)) 0. Since

we can then identify a lower bound for q: q(x, t) u 0 (c(x, t)) for all (x, t). (D.1) has a unique solution for l as a function of p, P and q if its left hand side (LHS) and right hand side (RHS) cross once. On the one hand, lim l!0 LHS = +1 and, when l = 0, RHS is equal to a non-negative constant. On the other, LHS is equal to a non-negative constant, when l = 1 f (x), and lim l!(1 f (x)) RHS = +1. Since @LHS/@l  0 and @RHS/@l 0, both LHS and RHS are monotone functions. Consequently, LHS and RHS cross only once implying that l is uniquely determined by p and q in every t 0 (notice that, by definition, P is a unique function of p).

E Proposition 4 proof

To prove the existence and uniqueness of solution to S, we provide a version of Theorem 3.4 in [START_REF] Pao | Nonlinear Parabolic and Elliptic Equations[END_REF]. In this regard, we should introduce first the notion of upper and lower solution of parabolic equations since our time-invariant equilibrium is described by this type of equations.

Afterwards, we state and prove a less constraining version of the original theorem by Pao.

Finally, we make use of this results to demonstrate our proposition.

Notice that, in the sequel of the present and in the next appendix, the problem needs to be defined on a bounded domain. Therefore, we assume without loss of generality that R is a finite real interval in R, R = (0, r), where r is a real finite number.

Let us establish the definition of upper and lower solutions:

Definition. u ⇤ is an upper solution of equation u xx = f (x, u), for x 2 R and u 2 R n , (E.1)